The Nusselt Problem
Author's Data: Housam BINOUS
Department of Chemical Engineering
National Institute of Applied Sciences and Technology
Tunis, TUNISIA
Email: binoushousam@yahoo.com

Acknowledgement :

Governing equations and analytical solution of the nusselt problem can be found in: Rice and Do, Applied Mathematics and Modeling for Chemical Engineers, Wiley, 1995.

Off[General::"spell1"]

Solution using NDSolve

sol = 
      NDSolve[{D[θ[τ, ξ], {τ, 
                1}] ==  (1/ξ) 
              D[ξ D[θ[τ, ξ], {ξ, 1}], {ξ, 
                  1}], 
θ[0, ξ] == 
            Exp[-1000 (1 - ξ)],  θ[τ, 1] ==  
            1, (D[θ[τ, ξ], {ξ, 1}] /. ξ -> 
                  10^-8) ==  
            0}, 
{θ[τ, ξ]}, {τ, 0, 
          10}, {ξ, 10^-8, 1}]

τ = 1;

plt1 = 
      Plot[First[θ[τ, ξ] /. sol], {ξ, 10^-8, 1}, 
        PlotRange ->  All, PlotStyle ->  RGBColor[1, 0, 0], 
        DisplayFunction ->  Identity]

τ = 0.1;

plt2 = 
      Plot[First[θ[τ, ξ] /. sol], {ξ, 10^-8, 1}, 
        PlotRange ->  All, PlotStyle ->  RGBColor[0, 0, 1], 
        DisplayFunction ->  Identity]

τ = 0.01;

plt3 = 
      Plot[First[θ[τ, ξ] /. sol], {ξ, 10^-8, 1}, 
        PlotRange ->  All, PlotStyle ->  RGBColor[0, 1, 0], 
        DisplayFunction ->  Identity]

Show[plt1, plt2, plt3, 
      DisplayFunction ->  $DisplayFunction]

Analytical Solution gives same result as the one obtained using NDSolve

n = 1;

While[
      n < 10, {(Subscript[λ, n]) = FindRoot[
              BesselJ[0, x] ==  0, {x, n Pi}][[1, 2]], 
        Print[(Subscript[λ, n])], n++}]

θanal[τ_, ξ_] := 
      1 - 2 (∑_(n = 1)^9) BesselJ[
                    0, (Subscript[λ, n]) ξ]/((Subscript[λ, n]) 
                    BesselJ[1, (Subscript[λ, n])]) 
              Exp[-(Subscript[λ, n])^2 τ]

τ = 1;

plt4 = 
      Plot[First[θanal[τ, ξ] /. sol], {ξ, 10^-8, 1}, 
        PlotRange ->  All, PlotStyle ->  RGBColor[1, 0, 0], 
        DisplayFunction ->  Identity]

τ = 0.1;

plt5 = 
      Plot[First[θanal[τ, ξ] /. sol], {ξ, 10^-8, 1}, 
        PlotRange ->  All, PlotStyle ->  RGBColor[0, 0, 1], 
        DisplayFunction ->  Identity]

τ = 0.01;

plt6 = 
      Plot[First[θanal[τ, ξ] /. sol], {ξ, 10^-8, 1}, 
        PlotRange ->  All, PlotStyle ->  RGBColor[0, 1, 0], 
        DisplayFunction ->  Identity]

Show[plt1, plt2, plt3, plt4, plt5, plt6, 
      DisplayFunction ->  $DisplayFunction]