2
$\begingroup$

I have found a transient solution to the 1D heat equation where the initial condition is discontinuous. The results are accurate except for when time is small. The initial condition looks like this:

Initial Condition

$\qquad $ Initial Condition

After solving with NDSolveValue using "PDEDiscretization" the solution at time 0 looks like this:

Solution at Time t=0

$\qquad $ Solution at Time t = 0

There is clearly some instability near the discontinuity that is resulting in inaccurate evaluation of the temperature near the simulated boundary. The solution stabalizes with very little time such that the solution at time 0.01 looks like this:

Solution at Time t=0.01

$\qquad $ Solution at Time t= 0.01

I am happy with the long term results but for my purposes i need an accurate solution within an extremely small time-frame (micro seconds, t ~ .000001). Can anyone please suggest ways to eliminate these issues that are being caused by the discontinuity?

<< "NDSolve`FEM`"

ts = .000250; tsl = .000250; (* m *)
ρs = 3980;  ρsl = 958; (* kg/m3 *)
ks = .035;  ksl = .00067; (* kW/m/K *)
cs = .75; csl = 4.22; (* kJ/kg/K *)

ρ = If[0 <= z < ts, ρs, ρsl];
k = If[0 <= z < ts, ks, ksl];
c =  If[0 <= z < ts, cs, csl];

eqn1 = k*D[T1[z], z] + 28; 
Tbl = 100;
Subscript[Γ1, D] = DirichletCondition[T1[z] == Tbl, z == ts + tsl];

BCr = NDSolveValue[{eqn1 == 0, Subscript[Γ1, D]}, T1, {z, 0, ts + tsl}];

Plot[BCr[z], {z, 0, ts + tsl}, GridLines -> {{ts}, {0}}]

Ti[z_] :=  Piecewise[{{BCr[z], 0 <= z < ts}, {100, True}}];
Plot[Ti[z], {z, 0, ts + tsl}, GridLines -> {{ts}, {0}}]

eqn2 = ρ*c*D[T[t, z], t] - k*D[T[t, z], z, z];
Subscript[Γ, D] = DirichletCondition[T[t, z] == Tbl, z == ts + tsl];
Subscript[Γ, N] = NeumannValue[-28, z == 0];

soln = 
  NDSolveValue[
    {eqn2 == Subscript[Γ, N], Subscript[Γ, D], T[0, z] == Ti[z]}, 
    T, {t, 0, 1}, {z, 0, ts + tsl}, 
    Method -> 
      {"PDEDiscretization" -> {"MethodOfLines", 
       "SpatialDiscretization" -> {"FiniteElement"}}}]

Animate[
  Plot[soln[t, z], {z, 0, ts + tsl}, 
    GridLines -> {{ts, ts + .00005}, {0}}, 
    PlotRange -> {99, 111}], 
  {t, 0, 1}]
Plot[soln[t, ts + .00005], {t, 0, 1}, 
  GridLines -> {{0}, {107}}, 
  PlotRange -> Full]
Plot[soln[0.01, z], {z, 0, ts + tsl}, GridLines -> {{ts}, {0}}]

EDIT: I added my code.

Improve this question
$\endgroup$
5
  • 1
    $\begingroup$ You can always try increasing the number of points in the spatial discretization with the "MinPoints" option or try a lower "DifferenceOrder" to decrease overshooting. Method Of Lines has some details. $\endgroup$
    – Thies Heidecke
    Commented Nov 26, 2015 at 1:06
  • 2
    $\begingroup$ You can experiment with "SpatialDiscretization" -> {"FiniteElement", "MeshOptions" -> {"MeshOrder" -> 1}} or things like "SpatialDiscretization" -> {"FiniteElement", "MeshOptions" -> {"MaxCellMeasure" -> (ts + tsl)/1011}} $\endgroup$
    – user21
    Commented Nov 26, 2015 at 1:26
  • 2
    $\begingroup$ @ThiesHeidecke, OP uses the FEM and equivalent of MinPoints is changing the MaxCellMeasure and the difference order reduction can be done with specifying a "MeshOrder". $\endgroup$
    – user21
    Commented Nov 26, 2015 at 1:28
  • 1
    $\begingroup$ @user21 Thanks for clarifying! $\endgroup$
    – Thies Heidecke
    Commented Nov 26, 2015 at 12:55
  • $\begingroup$ Thanks for the suggestions! the "MaxCellMeasure" option seems to have helped significantly. I am getting much better solutions. $\endgroup$
    – dowlguest
    Commented Nov 29, 2015 at 3:05

0

Reset to default

Your Answer

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