Skip to main content
4 of 6
added 2 characters in body
xzczd
  • 68.4k
  • 9
  • 174
  • 489

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

enter image description here

"Then, is it possible to avoid using pure FiniteElementin 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:

enter image description here

The approach is excessive for this toy example, but can be useful for more complicated problem.

#Remark

  1. Without the D[…, t] in definition of ode, the program will still find a solution, but it's not that accurate, maybe because NDSolve can handle differential equation involving derivative of t better.

  2. solinitlst is necessary, probably because the DAE solver can't find accurate enough initial conditions automatically in this case.

xzczd
  • 68.4k
  • 9
  • 174
  • 489