Pensieve header: Finding a power series solution of an arbitrary ODE of the form $y’=f(x,y)$.

PowerSeriesSolve[f_,x0_, y0_, n_] := (
ϕ[0] = y0;
Do[
ϕ[k] = y0+(∫_x0^x)(Normal[Series[f /. y->ϕ[k-1], {x,x0,k-1}]]  /. x->t)\[DifferentialD]t,
{k, 1, n}
];
ϕ[n]
);

PowerSeriesSolve[Sqrt[1+y^2], 0, 0,5]

PowerSeriesSolve2[f_,x0_, y0_, n_] := Module[{ϕ=y0},
Do[
ϕ += ((D[f /. y->ϕ, {x, k-1}] /. x->x0)/(k!))(x^k),
{k, 1, n}
];
ϕ
];

PowerSeriesSolve2[Sqrt[1+y^2], 0, 0,10]

Series[Sinh[x], {x,0,10}]