3
$\begingroup$

I am trying to find the positive solution to this nonlinear ODE. Theoretically, it has been proven that this problem has a unique positive solution for lam large. But, Mathematica can detect only the zero solution (zero is a solution for all lam's). Any techniques to detect the positive solution?

lam = 30;
sol = NDSolve[{-D[f[x], x, x] == lam*f[x]*(1 - f[x]), 
(D[f[x], x] /. x -> 1) == -Sqrt[lam]*0.1*f[1], (*first boundary condition*) 
(D[f[x], x] /. x -> 0) == Sqrt[lam]*0.1*f[0]}, (*second boundary condition*) 
f, {x, 0, 1}];
Plot[f[x] /. sol, {x, 0, 1}, AxesLabel -> {x, f}, LabelStyle -> Directive[Black, Bold, 12]]
$\endgroup$

1 Answer 1

3
$\begingroup$

Let us try a shooting method. The following code produces the solution for a given value of f[0]:

lam = 30;
sol[f0_?NumericQ]:=NDSolveValue[{-f''[x]==lam*f[x]*(1-f[x]),
                                 f[0]==f0,f'[0]==Sqrt[lam]*0.1*f0},f,{x,0,1}];

We want to understand for which f[0] the boundary condition at x==1 is satisfied. This condition is that the following must be zero:

cond1[f_]:=f'[1]-(-Sqrt[lam]*0.1*f[1])

Let us plot this:

Plot[cond1[sol[f0]],{f0,-0.53,0.91},AxesLabel->{"f[0]","cond1"}]

enter image description here

I picked the plot domain after some playing around. We see four values of f[0] where the condition at x==1 is also satisfied. I think the one OP asks for is the one on the very right, so let us see where it is:

f0x=f0/.FindRoot[cond1[sol[f0]],{f0,0.9}]
(* 0.905557 *)

Plot the corresponding solution:

With[{s=sol[f0x]}, Plot[s[x],{x,0,1},AxesLabel->{"x","f[x]"}]]

enter image description here

It is everywhere positive. It is also not difficult to guess that this solution is symmetric, $f(x) = f(1-x)$, and a numerical check confirms this. No error or warning message was produced, so I am assuming we can trust the numerical evaluation. But of course, one should check.

Note. The above is a manual implementation of the shooting method. Actually NDSolve has a shooting method built in. To see how this works, see this tutorial. Here:

lam=30;
sol2[f0guess_?NumericQ]:=NDSolveValue[{
                     -f''[x]==lam*f[x]*(1-f[x]),
                     f'[0]==Sqrt[lam]*0.1*f[0],
                     f'[1]==-Sqrt[lam]*0.1*f[1]},
                     f,x,Method->{"Shooting","StartingInitialConditions"->
                          {f[0]==f0guess,f'[0]==Sqrt[lam]*0.1*f0guess}}];

Let us use this to plot the four solutions mentioned above:

With[{sols=Map[sol2[#][x]&,{-0.5,0,0.85,0.9}]},
  Plot[Evaluate[sols],{x,0,1},PlotStyle->{Black,Blue,Green,Magenta},AxesLabel->{"x","f[x]"}]]

enter image description here

Note 2. The quantity 1/2*f'[x]^2+lam*(1/2*f[x]^2-1/3*f[x]^3) is independent of x for every solution.

$\endgroup$
0

Your Answer

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

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