1
$\begingroup$

I am trying to solve the heat conduction equation in cylindrical coordinates with a defined heat generation within the volume (a thin disk, with different heat flux BC top and bottom). The problem is asimuthally symmetric so I only have derivatives in "r" and "z". I am trying to use the method discussed on this site for solving inhomogeneous PDEs. I was able to get a solution when using a thickness equal to the radius (zt=1), but when I reduce the thickness (zt=0.2) I get the following errors

NDSolve::mxst: Maximum number of 10000 steps reached at the point t == 1.5022438114283172`. >>

NDSolve::eerr: Warning: Scaled local spatial error estimate of 11.936241527250862` at t = 1.5022438114283172` in the direction of independent variable r is much greater than prescribed error tolerance. Grid spacing with 15 points may be too large to achieve the desired accuracy or precision. A singularity may have formed or you may want to specify a smaller grid spacing using the MaxStepSize or MinPoints method options. >>

How can I fix this (I am using Mathematica Ver 8)?

Code:

qdot = 10; (*heat input*)
k = 1;  (*thermal conductivity*)
ti = 20; (*iniital temp*)
hb = 1; (*convection coeff back*)
hf = .01; (*convection coeff front*)
tb = 20; (*coolant temp back*)
tf = 20; (*coolant temp front*)
r0 = 1; (*radius*)
zt = .2; (*thickness*)


Clear[y];
f1[r_ /; r > 0] := 1;
f1[r_ /; r == 0] := 2;
f2[r_ /; r > 0] := 1/r;
f2[r_ /; r == 0] := 0;
qdot/k + Derivative[0, 0, 2][y][t, r, z] + 
  f2[r]*Derivative[0, 1, 0][y][t, r, z] + 
  f1[r]*Derivative[0, 2, 0][y][t, r, z] == 
 Derivative[1, 0, 0][y][t, r, z]

eqns = {qdot/k + Derivative[0, 0, 2][y][t, r, z] + 
     f2[r]*Derivative[0, 1, 0][y][t, r, z] + 
     f1[r]*Derivative[0, 2, 0][y][t, r, z] == 
    Derivative[1, 0, 0][y][t, r, z], 
   Derivative[0, 1, 0][y][t, r0, z] == 0, 
   Derivative[0, 0, 1][y][t, r, zt] == -.5, 
   Derivative[0, 1, 0][y][t, 0, z] == 0, 
   Derivative[0, 0, 1][y][t, r, 0] == +1, 
   y[0, r, z] == z - z^2/(2*zt)*1.5};

y[t_, r_, z_] = y[t, r, z] /.
   First[NDSolve[eqns, 
     y[t, r, z], {t, 0, 1.6}, {r, 0, r0}, {z, 0, zt}, 
     Method -> {"MethodOfLines", Method -> "StiffnessSwitching", 
       "DifferentiateBoundaryConditions" -> {True, 
         "ScaleFactor" -> 0}}]];
$\endgroup$
  • $\begingroup$ You may be able to obtain a solution for zt = 1, but it does not look like a credible solution to me. Clearly, there are problems at r = z = 0, and perhaps elsewhere. $\endgroup$ – bbgodfrey Feb 19 '15 at 3:00

Your Answer

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

Browse other questions tagged or ask your own question.