# Plotting heat equation in a Manipulate expression [duplicate]

I have a 2D heat equation $u_t = \alpha (u_{xx} + u_{yy})$ with conditions:

• $u(x, y, 0) = 300$,
• $u_y(x, 0, t) = \mu_1(x)$,
• $u_y(x, 1, t) = \mu_2(x)$,
• $u(0, y, t) = \mu_3(y)$,
• $u(1, y, t) = \mu_4(y)$,

where:

• $\mu_1(x) = c_1$ if $|x - 0.5| < 0.25,$ else $0$,
• $\mu_2(x) = c_2$ if $|x - 0.5| < 0.25,$ else $0$.
• $\mu_3(y) = \mu_4(y) = T_0$

I am trying to solve it with the code, given below:

c1 = 100.0; c2 = -500.0; T0 = 500;

a = 0.002;

phi[x_, y_] = 300;

mu1[x_] := If[Abs[x - 0.5] <= 0.25, c1, 0];
mu2[x_] := If[Abs[x - 0.5] <= 0.25, c2, 0];

mu3[y_] := T0;
mu4[y_] := T0;

Off[NDSolve::ibcinc];

Manipulate[
Module[{plt, sol},
sol = NDSolve[
{
D[u[x, y, t], t] == a (D[u[x, y, t], x, x] + D[u[x, y, t], y, y]),
u[x, y, 0] == phi[x, y],
(D[u[x, y, t], y] /. y -> 0) == mu1[x],
(D[u[x, y, t], y] /. y -> 1) == mu2[x],
u[0, y, t] == mu3[y],
u[1, y, t] == mu4[y]
},
u[x, y, t], {x, 0, 1}, {y, 0, 1}, {t, 0, 100}];
plt = ContourPlot[
u[x, y, t] /. sol /. t -> time, {x, 0, 1}, {y, 0, 1},
Contours -> 15,
ColorFunction -> "TemperatureMap",
PlotLabel -> "Temperature profile",
ContourLabels -> True,
ImageSize -> Large]
],
{{time, 5, "Time (t)"}, 5, 100, 0.5}
]


The problem is that I set $c_1, c_2$ to different values (even functions of $x$) but I can't see any difference on the plot. I think I should see some effect. If my code is incorrect, how do I fix it?

• Why did you silence the error message? u[0, 0, 0] = phy[0,0] = 300 != u[0, 0, 0] = mu3[y] = 500 May 3, 2013 at 13:49
• If I set $T_0 = 300$, I will get plot full of anomalies and even singularities. But even if I do so (and get equally red plot), why doesn't heat flow affect the solution? May 3, 2013 at 14:27
• I am not sure what NDSolve is returning as a solution with those conflicting boundary conditions. Did you come up with the numbers yourself? Otherwise, where do they come from? May 3, 2013 at 14:34
• It is one of many options. Actually I am going to find analytical solution of this for any $c_1, c_2, T_0$. But also I need visualization of it for some non-trivial $c_1, c_2, T_0$. But even if I remove inconsistency by considering $T_0 = \mu_3 = \mu_4$, I still get the same results for any $c_1, c_2$. But actually the final solution depends on them. May 3, 2013 at 14:41
• Mixo123, I managed to find a complete answer for your question, see my update below. Sep 28, 2016 at 9:04

# Update

Recently I revisited the obscure tutorial for "MethodOfLines" and figured out how to fix your problem i.e. make NDSolve use (D[u[x, y, t], y] /. y -> 0) == mu1[x] and (D[u[x, y, t], y] /. y -> 1) == mu2[x] for solving the equation rather than ignore them. What you need to do is:

(1) Add Method -> {"MethodOfLines", "DifferentiateBoundaryConditions" ->{True, "ScaleFactor" -> 1}} to NDSolve.

This manual option adjustment is necessary because NDSolve will automatically set "ScaleFactor" -> 0 for non-Dirichlet boundary condition, when i.c. and b.c. are inconsistent, this will make the Neumann inconsistent b.c. that's not a function of t be completely ignored.

(2) Use UnitStep instead of Piecewise or If to define mu1 and mu2.

It's probably a bug, luckily we have SimplifyPWToUnitStep and PiecewiseExpand so don't need to reconstruct the function by hand.

The following is the fixed code:

c1 = -200.0; c2 = -200.0;
T0 = 500;
a = 0.002;
phi[x_, y_] = 300;

mu1[x_] = SimplifyPWToUnitStep@PiecewiseExpand@If[Abs[x - 0.5] <= 0.25, c1, 0]
mu2[x_] = SimplifyPWToUnitStep@PiecewiseExpand@If[Abs[x - 0.5] <= 0.25, c2, 0]

mu3[y_] := T0;
mu4[y_] := T0;

sol = u /. First@
NDSolve[{D[u[x, y, t], t] == a (D[u[x, y, t], x, x] + D[u[x, y, t], y, y]),
u[x, y, 0] == phi[x, y], (D[u[x, y, t], y] /. y -> 0) ==
mu1[x], (D[u[x, y, t], y] /. y -> 1) == mu2[x], u[0, y, t] == mu3[y],
u[1, y, t] == mu4[y]}, u, {x, 0, 1}, {y, 0, 1}, {t, 0, 100},
Method -> {"MethodOfLines",
"DifferentiateBoundaryConditions" ->{True, "ScaleFactor" -> 1}}]

Manipulate[ContourPlot[sol[x, y, time], {x, 0, 1}, {y, 0, 1}, Contours -> 15,
ColorFunction -> "TemperatureMap", PlotLabel -> "Temperature profile",
ContourLabels -> True, ImageSize -> Large]
, {{time, 5, "Time (t)"}, 5, 100, 0.5}
]


Now mu1 and mu2 apparently plays roles in the solution.

Well, this answer is quite incomplete because I can't fix the problem. I post this answer just to point out that the true reason for the changeless sol is the BCs (D[u[x, y, t], y] /. y -> 0) == mu1[x] and (D[u[x, y, t], y] /. y -> 1) == mu2[x] are largely ignored by NDSolve. If you try:

Clear["*"]
T0 = 500; a = 0.002; phi[x_, y_] = 300;

(*   I changed your definitions for mu1[x] and mu2[x] a little
since it doesn't work at least for version 8. *)
mu1[x_] := Piecewise[{{c1, Abs[x - 0.5] <= 0.25}}];
mu2[x_] := Piecewise[{{c2, Abs[x - 0.5] <= 0.25}}];
mu3[y_] := T0;
mu4[y_] := T0;

sol[c1_, c2_] := NDSolve[{D[u[x, y, t], t] == a (D[u[x, y, t], x, x] + D[u[x, y, t], y, y]),
u[x, y, 0] == phi[x, y],
(D[u[x, y, t], y] /. y -> 0) == mu1[x],
(D[u[x, y, t], y] /. y -> 1) == mu2[x],
u[0, y, t] == mu3[y],
u[1, y, t] == mu4[y]},
u[x, y, t], {x, 0, 1}, {y, 0, 1}, {t, 0, 100}];

Plot3D[Evaluate[D[u[x, y, t] /. sol[100, -500], y] /. y -> 1],
{x, 0, 1}, {t, 0, 100}, PlotRange -> All]
Plot3D[Evaluate[D[u[x, y, t] /. sol[100, -500], y] /. y -> 0],
{x, 0, 1}, {t, 0, 100}, PlotRange -> All]


you'll see:

Apparently mu1[x] and mu2[x] are almostly (maybe essentially…) not used. You can try several other c1 and c2 or even take off the BCs for $y=0$ and $y=1$ to confirm.

I guess some changes for the Method will help but haven't succeeded until now… Look forward to someone who's more experienced in solving PDEs！

I think that one of the main issues related to the perceived responsiveness of Manipulate is still missing and should be addressed in this thread. Same as SynchronousUpdating, TrackedSymbols is an option which must be taken into account when using Manipulate.
In general, if you have a toy application for Manipulate you won't need any of the two options, however, the more complicated your application the more likely you will need to use these two options.
TrackedSymbols is to tell Manipulate when to trigger an update. If you don't specify which symbols to track Manipulate may update much more often than what you would expect. In general you may think of Manipulate as a Mathematica expression that is wrapped with Dynamic and that will thus reevaluate its contents whenever any symbols in the expression changes. To avoid continuous and unexpected reevaluations use TrackedSymbols to specify which symbols shuld trigger an update and your Manipulate will most probably start behaving as you expect.

One more comment: the more complicated your application is the more probably Manipulate will start becoming a bad tool to address your problem. In this case you should start considering using DynamicModule which is basically the underlying backbone of Manipulate. I have worked extensively with Dynamic tools and widgets (this screencasts shows a good example) and I can tell that starting to use DynamicModule requires a fairly large amount of practice and study time.

• Welcome to Mathematica.SE! I edited your post and improved formatting, this can be helpful How do I format my posts using Markdown or HTML?. Dec 16, 2013 at 16:56
• Great to see you here, Ariel! Your expertise will surely be much appreciated in this community. Please don't be put off by certain differences in the cultures of M SE and e.g. MathGroup (as some of other MathGroup veterans alas were), this place is in fact much friendlier than it might look from some distance and on the first sight. Dec 23, 2013 at 11:10

There are two things about your Manipulate that stand out as possible sources of trouble.

1. The use of Module to localize variables within a Manipulate isn't a good idea. It is better to introduce an invisible control to create a local variable. Also note that the variable plt isn't needed at all.

2. For complex computations such as you are carrying out, it is a good idea to use the SynchronousUpdating -> False option. I suspect that not having this option set to False is the main source of your problem.

Making the required changes to your code to correct these two problems produces a Manipulate in which the plot the dynamically changes when the time slider is moved.

Manipulate[
sol = NDSolve[{
D[u[x, y, t], t] == a (D[u[x, y, t], x, x] + D[u[x, y, t], y, y]),
u[x, y, 0] == phi[x, y],
(D[u[x, y, t], y] /. y -> 0) == mu1[x],
(D[u[x, y, t], y] /. y -> 1) == mu2[x],
u[0, y, t] == mu3[y],
u[1, y, t] == mu4[y]},
u[x, y, t],
{x, 0, 1}, {y, 0, 1}, {t, 0, 100}];
ContourPlot[u[x, y, t] /. sol /. t -> time, {x, 0, 1}, {y, 0, 1},
Contours -> 15,
ColorFunction -> "TemperatureMap",
PlotLabel -> "Temperature profile",
ContourLabels -> True],
{{time, 5, "Time (t)"}, 5, 100, 0.5},
{sol, ControlType -> None},
SynchronousUpdating -> False]


Here is the plot a t =5.

Here is the plot a t =30.

• Thank you for optimizations! I will definitely use them. But the problem was in fixing incorrect results from NDSolve or understanding why are they correct. I don't get why results are the same for any $c_1, c_2$: heat flow should somehow have an influence on solution. May 3, 2013 at 14:18
• @m_goldberg: why do you say that localizing variables with Module within a Manipulate isn't a good idea? I don't think that in general there is anything wrong with doing that. Especially I can't see how defining a hidden dynamic variable has any advantages over a Module local variable for this case. Why do you think it's better? May 3, 2013 at 16:20
• @AlbertRetey Well, for one thing, using module inside manipulate creates a lot of symbols. So I would use Module around manipulate. But in that case it's better to use DynamicModule as that would retain values between sessions. Compare for example Manipulate[Module[{z}, {z, a}], {a, 0, 2}] and Module[{z}, Manipulate[{z, a}, {a, 0, 2}]]; Dec 16, 2013 at 22:40
• @Ajasja: while it ist true that Module does create a lot of symbols for the local variables if used inside a Manipulate these usually are temporary and removed after the Module finishes (with the usual exceptions as in your examples). A Manipulate inside a Module will, on the other hand, only generate symbols for the local variables once, but these will always leak if used in the Manipulate and remain existing -- which I think is potentially a larger problem... Dec 17, 2013 at 8:23
• @Ajasja: ... here are two examples with a more appropriate usage of Module which show what I am trying to say: Manipulate[Module[{z}, z = a^2; z], {a, 0, 2}] vs. Module[{y}, Manipulate[y = a^2; y, {a, 0, 2}]]. Interestingly by checking with Names["Global*"] I found that Manipulate itself does leave behind a mess of new symbols for its own local variables, which I honestly do consider a "dirty" job, especially considering that DynamicModule` doesn't show the problem... Dec 17, 2013 at 8:34