4
$\begingroup$
sol = NDSolve[{y''[x] + (y[x])^2 == 1, y[0] == 2, y'[0] == 3}, y[x], {x, 0, 30}]
g = Plot[Evaluate[{y[x], y'[x], y''[x]} /. sol], {x, 0, 30}, PlotStyle -> Automatic]

I'm trying to coding this equation $y''+ y^2 =1$ with initial value $y(0)= 2$, and $y'(0)= 3$ using NDSolve. Anyone could tell me what is the problem above in the code? Thank you.

$\endgroup$
7
  • $\begingroup$ It does solve the system, but I get a message which says that there is a singularity at x = 3.95. $\endgroup$
    – Lotus
    Oct 12, 2018 at 5:13
  • $\begingroup$ Also, by changing the initial condition to y'[0] == 0 gives a nice periodic solution. Is this what you want ? $\endgroup$
    – Lotus
    Oct 12, 2018 at 5:17
  • $\begingroup$ I need it to give me a graph, but I tried old method I learned before, and it didn't work. $\endgroup$
    – Mousa
    Oct 12, 2018 at 5:22
  • $\begingroup$ How to change the initial value for y ' [0]=3 to y ' [0]=0 when the given initial in the problem is this one "y ' [0]=3" $\endgroup$
    – Mousa
    Oct 12, 2018 at 5:25
  • $\begingroup$ I am not saying that we should change it. I am just pointing out that with y'[0] = 3 we get stiffness issues, which goes away with y'[0] = 0. Understanding of the physics of the problem may explain this behaviour. $\endgroup$
    – Lotus
    Oct 12, 2018 at 5:33

3 Answers 3

5
$\begingroup$

The desired solution can be obtained semi-analytically.

s = DSolveValue[y''[x] + (y[x])^2 == 1, y[x], x] /. {C[1] -> c1, C[2] -> c2}
(* -(-1)^(2/3) 6^(1/3) 
   WeierstrassP[(-(1/6))^(1/3) (c1 + x), {2 (-1)^(2/3) 6^(1/3), c2}] *)

Next, determine the constants by applying the boundary conditions.

bc1 = s /. x -> 0;
bc2 = D[s, x] /. x -> 0;
sc = FindRoot[{bc1 == 2, bc2 == 3}, {c1, I}, {c2, 1}] // Chop
(* {c1 -> -0.523953 + 1.79242 I, c2 -> -10.3333} *)

and plot.

Plot[Chop[s /. sc], {x, 0, 10}, ImageSize -> Large]

enter image description here

The solution is periodic in x. In light of the singularities in the solution, it is not surprising that NDSolve failed for x near 4.

Addendum

The period can be increased by decreasing y'[0] and becomes infinite for y'[0] == 0. A WorkingPrecision of at least 45 is needed in that case.

$\endgroup$
2
$\begingroup$

May not be exactly what OP wants but still:

sol = NDSolve[{y''[x] + (y[x])^2 == 1, y[0] == 2, y'[0] == 0}, 
  y, {x, 0, 30}]

g = Plot[Evaluate[{y[x], y'[x], y''[x]} /. sol], {x, 0, 30}, 
  PlotStyle -> Automatic]

Gives the graphs as required. With y'[0] == 3 as given by OP there seems to be a singularity/stiffness at x = 3.95

$\endgroup$
2
$\begingroup$

just a supplement to @Lotus answer:

Using ParametricNDSolveValue

Y = ParametricNDSolveValue[{y''[x] + (y[x])^2 == 1, y[0] == 2,y'[0] == ys}, y, {x, 0, 30}, ys]

you could do some parameter studies varying ys=y'[0]

Show[{
Plot[Table[Y[p][x], {p, .1, 3. , .1}] , {x, 0, 25},PlotRange -> {-5, 5}] ,
Plot[ Y[0][x] , {x, 0, 25}, PlotStyle -> {Thickness[.01], Red}], 
Plot[ Y[3][x] , {x, 0, 25}, PlotStyle -> { Green },PlotRange -> {-5, 5}]
}]

enter image description here

which show, that only y'[0]==0 gives convergent solution! (The green solution is y'[0]==3 OP asked for )

$\endgroup$
4
  • $\begingroup$ Increasing WorkingPrecision will give a different plot for y'[0]==0. $\endgroup$
    – bbgodfrey
    Oct 12, 2018 at 15:34
  • $\begingroup$ @bbgodfrey Thanks, trying ParametricNDSolveValue[... , WorkingPrecision -> 50] Plot[Y[3][x], {x, 0, 10}, PlotStyle -> {Green}, PlotRange -> {-5, 5} ] doesn't finish ... $\endgroup$ Oct 12, 2018 at 15:54
  • $\begingroup$ Works fine for me with NDSolveValue and Mathematica 11.3. Solution is flat at about -1 for large x, which agrees with the symbolic solution. $\endgroup$
    – bbgodfrey
    Oct 12, 2018 at 16:05
  • $\begingroup$ @bbgodfrey Thanks, MMA 11.0.1. NDSolve gives the solution in the range 0<x<3.95 $\endgroup$ Oct 12, 2018 at 19:25

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.