1
$\begingroup$

I am have trouble with what should be a fairly simple NDSolve operation for the unsteady state heat transfer in a rod. The rod starts off at 20oC, and its surface temperature is fixed at 37oC. I am fairly new to Mathematica, so I'm not sure what is going wrong! I keep getting the error `"Infinite expression 1/0. encountered."'

My code is:

 a = 0.010
    Cp = 2000
    k = 0.1
    rho = 900
    alpha = k/(rho*Cp)

    Tval[r_, t_] = NDSolve[{

   (* PDE - 1D Radial Heat Equation *)
   (1/r)*(D[(r*D[T[r, t], r]), r]) == (1/alpha)*D[T[r, t], t],

   (* Boundary Conditions *)
   T[a, t] == 37,
   (D[T[r, t], r] /. r -> 0) == 0,

   (* Initial Condition *)
   T[r, 0] == 20},

(* Define Variables *)
  T[r, t], {r, 0, a}, {t, 0, 1800}]

Plot[{Tval[r,0],Tval[r,100],Tval[r,200]},{r,0,a}]

If anyone can work out why I am getting division by zero that would be much appreciated!

Improve this question
$\endgroup$
8
  • $\begingroup$ I'm not sure, but I would tend to think your equation has an actual singularity at $r=0$ and that it is wrong. By rod, do you mean disk? Maybe this is useful: mathworld.wolfram.com/HeatConductionEquationDisk.html. $\endgroup$
    – anderstood
    Apr 4, 2017 at 14:50
  • 1
    $\begingroup$ Replace 0 by some small number, say 10^-4, as the inner boundary in r, and NDSolve will work. Define Tval[r_, t_] = T[r, t] /. Flatten@NDSolve[..., and Plot will work too. $\endgroup$
    – bbgodfrey
    Apr 4, 2017 at 15:34
  • $\begingroup$ However, T[a, t] == 37 and T[r, 0] == 20 together have the effect of assigning T[a, 0] two different values, which could be a problem. $\endgroup$
    – bbgodfrey
    Apr 4, 2017 at 15:42
  • $\begingroup$ @anderstood: The ODE is correct. It does have a singular point at $r = 0$, but it is a regular singular point, which means that you can find well-behaved series solutions in the neighborhood of $r = 0$. This is a common problem that has to be dealt with when doing PDEs in curvilinear coordinates; it arises from the fact that $\vec{\nabla}r$ is not well-defined at that point. $\endgroup$
    – Michael Seifert
    Apr 4, 2017 at 17:05
  • $\begingroup$ @MichaelSeifert So there is a singularity, induced by the parametrisation and not the physics, interesting. Thank you for the informative comment & link. $\endgroup$
    – anderstood
    Apr 4, 2017 at 18:09

1 Answer 1

Reset to default
2
$\begingroup$

This gives a solution when boundary and initial condition get matched with a very steep rise of the temperature at r=a and r starts from a very small value, not 0:

    Tval[r_, t_] = 
                  T[r, t] /. 
            First@NDSolve[{(*PDE-1D Radial Heat Equation*)(1/
                 r)*(D[(r*D[T[r, t], r]), r]) == (1/alpha)*D[T[r, t],t],
    T[a, t] == 37 (1 - (37 - 20)/37 E^(-8 t)), (D[T[r, t], r] /. r -> 10^-8) ==
                0, T[r, 0] == 20}, T, {r, 10^-8, a}, {t, 0, 1800}]
Improve this answer
$\endgroup$
0

Your Answer

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

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