2
$\begingroup$

I'm having some trouble plotting my results from NDSolve in an appealing way. Below is how I tried to solve the heat equation for axisymmetric, transient diffusion in the radial direction only as a function of time. There is an internal heat generation term and constant temperature at the radial boundry (ro).

ro = .3;
l = 1;
s = 25;
\[Rho] = 1050; 
c = 4097;
k = .55;
q = 19305.6*Exp[.000139*t];
bc1 = T[ro, t] == s;
bc4 = (D[T[r, t], r] /. r -> 0.0001) == 0.0;
ic = T[r, 0] == s;

heat = D[T[r, t], r, r] + (1/r)*D[T[r, t], r] + q/k == (([Rho]*c)/k)*D[T[r, t], t];

sol = NDSolve[{heat, bc1, bc4, ic}, T[r, t], {r, .0001, ro}, {t, 0, 10000}]

When I plot my results, I get cartesian coordinates with radius on one axis, time on another axis, and temperature on the third axis.

Plot3D[T[r, t] /. sol, {r, .0001, .3}, {t, 0, 10000}]

enter image description here

This result is correct and interpretable; however, I'm looking to get a circular plot that shows me temperature at each radial position as time evolves. Any thoughts on how to achieve this?

$\endgroup$

2 Answers 2

2
$\begingroup$
sol = NDSolve[{heat, bc1, bc4, ic}, T, {r, .0001, ro}, {t, 0, 10000}]
f[a_, b_] := T[a, b] /. sol[[1]]
tab = Table[
   ParametricPlot3D[
    Evaluate[{u Cos[v], u Sin[v], f[u, j]}], {u, 0.0001, 0.3}, {v, 0, 
     2 Pi}, PlotRange -> {{-0.3, 0.3}, {-0.3, 0.3}, {20, 150}}, 
    BoxRatios -> {1, 1, 1/2}, Mesh -> None], {j, 0, 10000, 100}];

tab exported as animated gif:

enter image description here

or color-coding temperature (rescaling 20 to 120):

tab = Table[
   ParametricPlot3D[
    Evaluate[{u Cos[v], u Sin[v], f[u, j]}], {u, 0.0001, 0.3}, {v, 0, 
     2 Pi}, PlotRange -> {{-0.3, 0.3}, {-0.3, 0.3}, {20, 150}}, 
    BoxRatios -> {1, 1, 1/2}, Mesh -> None, 
    ColorFunction -> 
     Function[{x, y, z}, ColorData["TemperatureMap"][(z - 20)/100]], 
    ColorFunctionScaling -> False], {j, 0, 10000, 500}];

enter image description here

$\endgroup$
3
  • $\begingroup$ This is perfect, thank you so much! One final question: is it possible to place a bar legend that tells me the temperature represented by each color? $\endgroup$ Mar 16, 2017 at 21:35
  • $\begingroup$ @ChaseKayrouz you can add a legend using Legended. Wrap the plot in Legended, e.g. Legended[plot,BarLegend[{"TemperatureMap", {20, 120}}]]. See the documentation for legend titles etc. Good luck :) $\endgroup$
    – ubpdqn
    Mar 17, 2017 at 2:48
  • 1
    $\begingroup$ You can also use 'RevolutionPlot3D[f[r, 10000], {r, 0, ro}]' instead of 'ParametricPlot3D', so you don't have to write out the parameterisation explicitly. $\endgroup$ Nov 17, 2020 at 12:24
2
$\begingroup$

Here is sort of a quick and dirty way. We replace r with Sqrt[x^2 + y^2] and define a piece-wise function as follows:

f[rr_, tt_] = If[rr >= 0.3 || rr <= 0.0001, 25, First@T[r, t] /. sol /. {r -> rr, t -> tt}]

Using If to give the values of the function on and outside the boundaries is necessary in order that we don't evaluate the InterpolatingFunction outside of its domain. Then, we Plot3D and Manipulate over time:

Manipulate[
  Plot3D[f[Sqrt[x^2 + y^2], t],
    {x, -0.3, 0.3}, {y, -0.3, 0.3},
    PlotRange -> {20, 140},
    PerformanceGoal -> "Quality"],
  {t, 0, 10000}]

enter image description here

$\endgroup$
2
  • $\begingroup$ Thanks for the help, this is what I'm looking for! Unfortunately, when I try this I get "$Aborted" in the plot window. $\endgroup$ Mar 16, 2017 at 6:07
  • $\begingroup$ @ChaseKayrouz. I see why. See the updated code. $\endgroup$
    – march
    Mar 16, 2017 at 16:04

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.