3
$\begingroup$

In this question, the equation $uu'=\nu u''$, $x\in (-1,1)$, $u(-1)=1+\delta$, $u(1)=-1$, was considered. The question asked about how to solve this differential equation problem in Mathematica. There were two answers. The first one suggested using the option Method -> FiniteElement in NDSolve for version 12.0. The second solution proposed the trapezoidal finite difference method, which was programed from scratch.

Let us analyze the approximations at $x=0.5$, for $\nu=0.05$ and $\delta=0.1$:

test = ParametricNDSolve[{u''[x] - (1/nu)*u[x]*u'[x] == 0, 
    u[-1] == 1 + delta, u[1] == -1}, u, {x, -1, 1}, {nu, delta}, 
   Method -> FiniteElement];

u[0.05, 0.1][0.5] /. test
(* 1.09905 *)
fdmODE[nu_, delta_, n_] := 
  Module[{h, mesh, f, u, v, eqns, sv, froot, sol},
   h = 2/n;
   mesh = -1 + h*Range[0, n];
   f[{u_, v_}] = {v, (1/nu)*u*v};
   eqns = 
    Flatten[Join[{u[0] == 1 + delta, u[n] == -1}, 
      Table[Thread[{u[i], v[i]} == {u[i - 1], v[i - 1]} + 
          0.5*h*(f[{u[i - 1], v[i - 1]}] + f[{u[i], v[i]}])], {i, 1, 
        n}]]];
   sv = Flatten[Table[{{u[i], 0}, {v[i], 0}}, {i, 0, n}], 1];(* 
   starting value, initial guess *)
   froot = FindRoot[eqns, sv];
   sol = Table[u[i], {i, 0, n}] /. froot;
   Return@Thread[{mesh, sol}];
   ];

fdmODEinterp[nu_, delta_, n_, x_] := 
 Interpolation[fdmODE[nu, delta, n], InterpolationOrder -> 1][x]

fdmODEinterp[0.05, 0.1, 10000, 0.5]
(* 1.09923 *)

fdmODEinterp[0.05, 0.1, 11000, 0.5]
(* 1.09923 *)

fdmODEinterp[0.05, 0.1, 12000, 0.5]
(* 1.09923 *)

It seems that Method -> FiniteElement is not sufficiently accurate. Is there any better option to tackle this problem in Mathematica?

$\endgroup$
  • 2
    $\begingroup$ Have you tried refining the mesh? $\endgroup$ – user21 Dec 11 '19 at 14:56
5
$\begingroup$

You could refine the mesh:

test = ParametricNDSolve[{u''[x] - (1/nu)*u[x]*u'[x] == 0, 
    u[-1] == 1 + delta, u[1] == -1}, u, {x, -1, 1}, {nu, delta}, 
   Method -> {"FiniteElement", 
     "MeshOptions" -> {"MaxCellMeasure" -> 0.01}}];

u[0.05, 0.1][0.5] /. test
1.09923
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.