It's because NDSolve
is using different method for these 2 problems. For the first problem, Shooting
method, which is for nonlinear boundary value problem of ordinaray differential equation, is used; while for the second problem, FiniteElement
method, which supports nonlinear problem as of v12, is used. Different solvers of NDSolve
are not equally powerful. For the second problem, we can obtain a result by adjusting the option a little. (The result is not perfect though):
sol = First@
NDSolve[{eqns, ic}, {r[R, t], TRR[R, t], γRR}, {R, 0.1, 1}, {t, 0, 1},
InitialSeeding -> {r[R,t] == 3, TRR[R,t] == 3, γRR[R,t] == 2},
Method -> {FiniteElement, MeshOptions -> MaxCellMeasure -> 0.05}]
PlotSol@sol
"Then, is it possible to avoid using pure FiniteElement
in NDSolve
?" Sadly the answer seems to be no, at least now, because adding Method -> {MethodOfLines, TemporalVariable -> t}
causes tvic
warning, which seems to indicate NDSolve
is having difficulty in noticing the system can actually be discretized to a DAE system of t
. Still, if we discretize the system outside of NDSolve
(I'll use pdetoode
for the task):
γθθ[R_, t_] := 4
bleft = 1/10; bright = 1;
With[{r = r[R, t], TRR = TRR[R, t], γRR = γRR[R, t], γθθ = γθθ[R, t]},
eqns = {r D[r, R] == γRR γθθ R,
r^4 R γθθ D[TRR, R] == 2 γRR (r^4 - R^4 γθθ^4),
D[γRR, t] == 0};
bc = {r == 0.1 /. R -> bleft, TRR == 0 /. R -> bright};
ic = γRR == 3 /. t -> 0;];
Block[{γRR}, γRR[R_, t_] := 3;
solinitlst = NDSolveValue[{eqns // Most, bc}, {r[R, t], TRR[R, t]}, {R, 0.1, 1}]];
icadd = {r[R, 0], TRR[R, 0]} == solinitlst // Thread;
tend = 1; points = 100;
grid = Array[# &, points, {bleft, bright}]; difforder = 2;
(* Definition of pdetoode isn't included in this post,
please find it in the link above. *)
ptoofunc = pdetoode[{r, TRR, γRR}[R, t], t, grid, difforder];
ode = {D[{Rest@ptoofunc@eqns[[1]], Most@ptoofunc@eqns[[2]]}, t], ptoofunc@eqns[[3]]};
odebc = bc // ptoofunc;
odeic = {ic // ptoofunc, icadd[[1]] // ptoofunc // Rest, icadd[[2]] // ptoofunc // Most};
var = Outer[#@#2 &, {r, TRR, γRR}, grid];
sollst = NDSolveValue[{ode, odebc, odeic}, var, {t, 0, tend}];
{solr, solTRR, solγRR} = rebuild[#, grid, 2] & /@ sollst;
Plot[{solr, solTRR, solγRR}[R, tend] // Through // Evaluate, {R, bleft, bright}]
NDSolve
will manage to find an accurate enough result:
The approach is excessive for this toy example, but can be useful for more complicated problem.
#Remark
Without the
D[…, t]
in definition ofode
, the program will still find a solution, but it's not that accurate, maybe becauseNDSolve
can handle differential equation involving derivative oft
better.solinitlst
is necessary, probably because the DAE solver can't find accurate enough initial conditions automatically in this case.