3
$\begingroup$

Consider the one-dimensional heat equation $$ \frac{\partial u(x,t)}{\partial t}=\alpha\frac{\partial^2 u(x,t)}{\partial x^2}, \quad (x,t)\in (0,L)\times (0,\infty), $$ subject to the Dirichlet boundary conditions $$ u(0,t)=u(L,t)=0, \quad t\in (0,\infty), $$ and initial condition $$ u(x,0)=f(x), \quad x\in (0,L). $$ A solution of the above problem is given by $$ u(x,t)=\sum_{n=1}^\infty A_n\sin\left(\frac{n\pi}{L}x\right)e^{-\alpha n^2\pi^2t/L^2} $$ I am hoping to simulate the above solution to this problem on a rod and export this as a gif, that is, we can visually see the heat propogating along the rod. I have sketched below an illustration of the simulation. Is this possible to do using Mathematica? enter image description here

Edit

As per my conversation with Nasser, the animation is of a thin 2D rod where the heat only flows along the bar in the horizontal direction (i.e. along the length of the bar).

$\endgroup$
9
  • 2
    $\begingroup$ What exactly is the problem? Did you have a look in the documentation? Try searching for heat transfer. $\endgroup$
    – user21
    Aug 6, 2023 at 5:15
  • $\begingroup$ @user21 I have searched extensively but am still having trouble constructing this simulation $\endgroup$
    – Bell
    Aug 6, 2023 at 5:17
  • 1
    $\begingroup$ Then share the code of what you have tried, and explain what did not go as expected. $\endgroup$
    – user21
    Aug 6, 2023 at 5:18
  • 1
    $\begingroup$ A useful tutorial from Brown University $\endgroup$
    – Syed
    Aug 6, 2023 at 5:20
  • 1
    $\begingroup$ Maybe this or this are useful to get started. $\endgroup$
    – user21
    Aug 6, 2023 at 5:39

1 Answer 1

4
$\begingroup$

This is basic animation using Manipulate. This can be solved analytically by Mathematica. Once you obtain the analytical solution then Manipulate is used to replace $f(x),L,k$ by any other values without having to solve the PDE again.

Manipulate is used for just changing the time and length and thermal conductivity. In this $f(x)=\sin(x)$ is used for initial conditions. You can change that as needed. Number of terms used in $20$ in the series. The more the more accurate it will be, but $20$ is more than enough. To save to animation gif file, there are many questions and answers on this site how to save Manipulate to gif file.

enter image description here

Code

ClearAll["Global`*"];
pde=D[u[x,t],t]==k*D[u[x,t],{x,2}];
bc={u[0,t]==0,u[L,t]==0};
ic=u[x,0]==f[x];
sol=DSolveValue[{pde,ic,bc},u[x,t],{x,t}];
sol=sol/. {K[1]->n,Infinity->20};

Mathematica graphics

Manipulate[
 Module[{sol2, g},
  sol2 = Activate[sol /. {f[x] -> Sin[x], L -> L0, k -> k0}];
  g = Grid[{{Row[{"time = ", t0}]},
     {Quiet@
       Plot[ sol2 /. t -> t0, {x, 0, L0}, 
        PlotRange -> {Automatic, {-1.2, 1.2}}, ImageSize -> 300, 
        GridLines -> Automatic, GridLinesStyle -> LightGray, 
        PlotStyle -> Red]}}];
  g
  ],
 {{L0, 3, "Length?"}, 1, 10, .1, Appearance -> "Labeled"},
 {{k0, .01, "Thermal conductivity?"}, 0.01, .1, .01, 
  Appearance -> "Labeled"},
 {{t0, 0.01, "time?"}, 0.01, 100, .1, Appearance -> "Labeled"},
 TrackedSymbols :> {L0, k0, t0}
 ]

Update

Here is a possibility, using DensityPlot. I am sure this can be improved more.

enter image description here

Code

ClearAll["Global`*"];
pde = D[u[x, t], t] == k*D[u[x, t], {x, 2}];
bc = {u[0, t] == 0, u[L, t] == 0};
ic = u[x, 0] == f[x];
sol = DSolveValue[{pde, ic, bc}, u[x, t], {x, t}];
sol = sol /. {K[1] -> n, Infinity -> 20};

Manipulate[
 Module[{currentSol, g},
  currentSol = Activate[sol /. {f[x] -> Sin[x], L -> L0, k -> k0}];
  g = Grid[{{Row[{"time = ", t0}]},
     {Quiet[
       DensityPlot[
        Evaluate[currentSol /. t -> t0], {x, 0, L0}, {y, 0, 1},
        ColorFunction -> Function[{x}, Hue[x]],
        ColorFunctionScaling -> False, AspectRatio -> 1/10, 
        ImageSize -> 300]]}}
    ];
  g
  ],
 {{L0, 3, "Length?"}, 1, 10, .1, Appearance -> "Labeled"},
 {{k0, .1, "Thermal conductivity?"}, 0.01, .1, .01, 
  Appearance -> "Labeled"},
 {{t0, 0.01, "time?"}, 0.01, 100, .1, Appearance -> "Labeled"},
 TrackedSymbols :> {L0, k0, t0}
 ]
$\endgroup$
7
  • $\begingroup$ Thanks for your response! I understand the manipulate command and the output you have produced. For my purposes, I want to animate the analytical solution on a cylinder not a plot (to represent the heat conduction along a thin metal rod). Does this make sense? $\endgroup$
    – Bell
    Aug 6, 2023 at 6:56
  • $\begingroup$ @Bell This is a 1-D pde you have given. So it is assumed that the bar is of very small thickness. Otherwise it is becomes a 2D diffusion problem actually. So the above shows how the heat changes along the length of the bar (1-D). $\endgroup$
    – Nasser
    Aug 6, 2023 at 7:13
  • $\begingroup$ +1 Good points, thank you. Out of curiosity, if we make the assumption that the rod is very thin, could this still be done in the 1-D case? I am most interested in this. $\endgroup$
    – Bell
    Aug 6, 2023 at 7:18
  • $\begingroup$ @Bell in 1-D diffusion it is always assumed the bar is very thin. So I am not sure I understand the question. Another way to look at it, is that there is no heat change across the width of the bar, and only heat changes/flows along the length of the bar. If you want to make 2D plot of this, ofcourse this is possible. Just have to make sure there is no change across the width of the bar and heat only flows in one direction which is along its length. I did not do that, as that will require making 2D graphics. May be someone else who have more time can do that. $\endgroup$
    – Nasser
    Aug 6, 2023 at 7:28
  • $\begingroup$ This is EXACTLY what I want to simulate - apologise for not making this clear! If we assume that heat can only move in one direction, along the length of the rod, this should work. I will continue to try and animate this. $\endgroup$
    – Bell
    Aug 6, 2023 at 7:50

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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