# Finding power series solution for differential equation in Mathematica [duplicate]

I know this topic has been covered before, but I've tried all the solutions I can find from other users' questions and none of them have worked.

I need to find a power series solution to the following nonlinear differential equation:

y''=x+y^2


with initial conditions

y(0)=0, y'(0)=0


up to the x^4 term. Does anyone have a quick and easy way to do this?

• Welcome to Mathematica StackExchange. In order to learn how to use this site take the tour. When copying equations from a notebook to your question one should format using inline code by selecting the code and clicking the {} button above the edit window. It is recommended that you browse the Markdown help – Jack LaVigne Nov 24 '17 at 20:40
• One of the reasons to do this is to help the users who attempt to answer your question so they can directly copy your equations into a notebook. Otherwise typographical errors can easily occur. – Jack LaVigne Nov 24 '17 at 21:09
• seriesDSolve[y''[x] == (x + y[x]^2), y, {x, 0, 13}, {y[0] -> 0, y'[0] -> 0}] from my answer to (25363) solves your problem. – Michael E2 Nov 24 '17 at 22:03
• @MichaelE2 Mathematica should really have Series option to solving ODE's. This option has been there in Maple for ever. I do not understand why Mathematica does not support such a basic option to solving ODE's. – Nasser Nov 25 '17 at 1:50

I do not get an x^4 term.

findSeriesSolution[y_,x_,nTerms_]:=Module[
{roots,pt=0,ode,s0,s1,ic,eq,sol},ic={y[0]->0,y'[0]->0};
ode=y''[x]- (x+y[x]^2);
s0=Series[ode,{x,pt,nTerms}];
s0=s0/.ic;
roots=Solve@LogicalExpand[s0==0];
s1=Series[y[x],{x,pt,nTerms+2}];
sol=Normal[s1]/.ic/.roots[[1]]
]


And now

sol = findSeriesSolution[y,x,10]
(*  x^3/6 + x^8/2016 *)


Compare

Plot[sol,{x,0,1}]


solm=NDSolve[{y''[x]-(x+y[x]^2)==0,y[0]==0,y'[0]==0},y,{x,0,1}];
Plot[Evaluate[y[x]/.solm],{x,0,1}]


Analytical solution is (thanks to Maple 2017.3)

$$-\ln \left( {\frac {1}{{{\rm Ai}\left(-x\right)}{{\rm Bi}^{(1)}\left( -x\right)}-{{\rm Ai}^{(1)}\left(-x\right)}{{\rm Bi}\left(-x\right)}} \left( {\frac { \left( {\it \_C2}\,\sqrt {3}\pi+{3}^{2/3}\Gamma \left( 2/3 \right) \right) {{\rm Ai}\left(-x\right)}}{\pi}}-{\it \_C2}\,{{\rm Bi}\left(-x\right)} \right) } \right) -2\,i\pi\,{\it \_Z2 }$$

Screen shot

Mathematica can't at this time solve this analytically. May be in future versions.

DSolve[{y''[x] - ( x + y[x]^2) == 0, y[0] == 0, y'[0] == 0}, y[x], x]