Dear wolfram community,

I hope my problem is clear and easy to solve.

I have already solved the following heat equation over a domain:


pde = 1/r D[r*D[T[t, r, z], {r, 1}], {r, 1}] + 
   D[D[T[t, r, z], {z, 1}], {z, 1}] - D[T[t, r, z], {t}];

\[CapitalDelta]z = 50*10^(-4);(*[m]*)
R = 5*10^(-2);(*[m]*)
\[Tau] = 3.5*10^(-3)(*[s]*);
Tw = 200;(*[K]*)
\[CapitalOmega] = 
 ImplicitRegion[0 <= r <= R && -2*\[CapitalDelta]z <= z <= 0, {r, z}];

bc = {DirichletCondition[T[t, r, z] == Tw, z == -2*\[CapitalDelta]z], 
   DirichletCondition[T[t, r, z] == Tw, r == R]};
(*Note that the unspecified boundaries are set to Neuman zero by \

ic = {T[0, r, z] == 300 (*[K]*)};

sol = NDSolve[{pde == 0, bc, ic}, 
  T, {t, 0, \[Tau]}, {r, z} \[Element] \[CapitalOmega]]

I optained an Interpolationfunction.

Manipulate[(*Plot mit relativen Farbspectrum*)
 Plot3D[T[t, r, z] /. sol, {r, 0, R}, {z, -2*\[CapitalDelta]z, 0},
  PlotRange -> {180, Tliq*1.1}, Mesh -> None,
  Axes -> True, AxesLabel -> {"r", "z", "T"} , 
  AxesStyle -> Thickness[0.001],
  ColorFunction -> "TemperatureMap",
  ColorFunctionScaling -> True,
  ImageSize -> Large,
  PlotLegends -> Automatic],
 {t, 0, \[Tau], \[Tau]*0.01}]

Now I want to use the temperature distribution T[[Tau],r,z] as part of the initial value (T[0,r,z]) for a new domain (omeganew ) and then solve the heat equation over omeganew. The new domain is essentially the old domain but enlarged by an amount of 4*Δz:

 omeganew = ImplicitRegion[0 <= r <= R && -4*Δz <= z <= 0, {r, z}];

Here is the tricky part: For the new solution I want the Temperature distribution from the previous solution to be "shifted" to the lower part of omeganew and a different function fc[r,z] (lets assume for the purpose of simplicity that fc[r,z]= constant) to be valid within the region that is newly added. So I need a shifted version of the previous solution with respect to the z-coordinate.

In a nutshel, I would like to extract the information (Temperaturevalues at an abitrary time t) obtained by NDsolve (e.g. sol) to construct a Unitstep function (is there a better alternative?) that I would then use as the initial condition for the new domain (omeganew).

Can anyone please help me?

  • $\begingroup$ What are the values of tau, R, deltaz, etc? Without them it's not possible to simulate. Please provide them. $\endgroup$ – user21 Jan 14 '19 at 11:22
  • $\begingroup$ [CapitalDelta]z = 50*10^(-4);([m]) ;[Tau] = 3.5*10^(-2); R = 5*10^(-2); ([m]) $\endgroup$ – Gustavo Meyagan Jan 14 '19 at 11:39
  • $\begingroup$ Please add them to your post. $\endgroup$ – user21 Jan 14 '19 at 11:43
  • $\begingroup$ ok I have added the additional information. Thanks in advance! $\endgroup$ – Gustavo Meyagan Jan 14 '19 at 11:48

Here is a way to do it. Let' set up the model:

Δz = 50*10^(-4);(*[m]*)R = 5*10^(-2);(*[m]*)τ = 
pde = 1/r D[r*D[T[t, r, z], {r, 1}], {r, 1}] + D[D[T[t, r, z], {z, 1}], {z, 1}] - D[T[t, r, z], {t}];
Ω = ImplicitRegion[0 <= r <= R && -2*Δz <= z <= 0, {r, z}];
bc = {DirichletCondition[T[t, r, z] == Tw, z == -2*Δz], DirichletCondition[T[t, r, z] == Tw, r == R]};

If you now call NDSolveValue you will get a solution (looks like it's zero but that is a different issue)

sol = NDSolveValue[{pde == 0, bc, T[0, r, z] == 0}, T, {t, 0, τ}, {r, z} ∈ Ω];

If you evaluate the inteprolating function out side of the region you will get a warning and an Indetermiante as an answer.

sol[0, -1, 3]

enter image description here


To change that you can use:

sol = NDSolveValue[{pde == 0, bc, T[0, r, z] == 0}, 
   {t, 0, τ}, 
   {r, z} ∈ Ω, 
   {"ExtrapolationHandler" -> {5 &, "WarningMessage" -> False}}

Now you will get the extrapolation value specified (5) and no warning:

sol[0, -1, 3]

With this you can then call NDSolveValue on a different domain with a different initial value like so:

sol2 = NDSolveValue[{pde == 0, bc, T[0, r, z] == sol[τ, r, z]}, 
   {t, 0, τ}, 
   {r, z} ∈ ImplicitRegion[ 0 <= r <= R && -10*Δz <= z <= 0, {r, z}]
| improve this answer | |
  • $\begingroup$ OHH I also forgot to provide the values for Tw and ic (The "wall temperature" and the initial temperature distribution respectively . They are now included in the post! $\endgroup$ – Gustavo Meyagan Jan 14 '19 at 12:10
  • $\begingroup$ I want sol[tau, r, z] to be the intial value for the domain1: 0 <= r <= R && -4*Δz <= z <= -Δz (e.g.sol[tau, 0, 0]=T[0,0,-Δz ],sol[tau, 0, -0.5*Δz ]=T[0,0,-0.5*Δz -Δz ] ,and so forth) I want to essentially shift the function downwards. And the initial values of the domain2 0 <= r <= R && -Δz <= z <= 0 to be specified by a completely different function ( lets say a constant) Of course sol2 is to be solved for the union of domain1 and domain2 $\endgroup$ – Gustavo Meyagan Jan 14 '19 at 12:50
  • 1
    $\begingroup$ @GustavoMeyagan, you need to specify everything in the question. Try to do this in a clean, clear way such that people can follow your line of thought. $\endgroup$ – user21 Jan 14 '19 at 15:14
  • $\begingroup$ @User21I think I've figured out the problem. The interpolation function obtained from NDSolve does not enable you to obtain the value at arbitrary points, for exaple sol[0, -1, 3] will not give the value 5. One actually needs to solve the pde with NDSolveValue instead. From there you can use sol[t,r,z] like any other function including embedding it in a unitstep function.Thank you very much for your help! You guys rock! $\endgroup$ – Gustavo Meyagan Jan 15 '19 at 10:25
  • $\begingroup$ @GustavoMeyagan, good that you figured it out. NDSolve and NDSolveValue essentially return the same result. NDSolve will give you a Rule where the second argument is the interpolating function. To use the result from NDSolve you first need to replace the dependent variable to get to the interpolating function. With NDSolveValue there is not such need as it directly returns the Interpolating function. $\endgroup$ – user21 Jan 15 '19 at 11:19

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.