ODE Expansion about a Regular Point

Example:  Expand solution to  y' [x]=1+y [x]^2 ,  y[0]=0 around x=0

Analytic Solution

Clean up old definifions of variables just in case

Remove["Global`*"]

Define the differential equation

deq  =  y'[x] == 1 + y[x]^2

And solve it

sol  =  DSolve [{deq, y[0] == 0}, y[x], x]

Hence, the expansion we are looking for is

Series [y[x] /.  sol,  {x, 0, 8}]

However, only very rarely are analytical solutions possible. Therefore, we need methods that do not rely on this.

Iteration

Remove["Global`*"]

Initialize:

n = 0 ;  y = 0 ;

Now, every time the following line is repeated, we get another term in the local expansion

\(\(
n = n + 2; \ \ yp = Normal[1 + y^2 + O[x]^\((n - 1)\)]; \ 
    y = Integrate[yp, {x, 0, ξ}] /. ξ -> x\)\)

The iteration approach is attractive in that we do not need to assume the precise form of the expansion.

Expansion

Remove["Global`*"]

Define a series of 8 terms. The initial condition is incorporated by omitting a constant term in the series

n = 8;  y[x_] := Sum[a[i] x^i, {i, 1, n}] + O[x]^(n + 1)

Now, substitute the series into the differential equation,  and generate separate equations for the different powers of x:

LogicalExpand[y'[x]==1+y[x]^2]

and solve these equations:

Solve[%]