Nonlinear ODE by an infinite series expansion

Printable View

September 5th 2009, 11:32 PM
javicg

Nonlinear ODE by an infinite series expansion

I have to solve the nonlinear DE y'=x²-y² by using an infinite series expansion y= $\sum_{n=0}^{\infty} a_n x^n$ , but I've tried in vain. Maybe a change of variables would make it easier, but I don't know which one.

Thanks
September 7th 2009, 10:54 AM
Coomast

Quote:

Originally Posted by javicg

I have to solve the nonlinear DE y'=x²-y² by using an infinite series expansion y= $\sum_{n=0}^{\infty} a_n x^n$ , but I've tried in vain. Maybe a change of variables would make it easier, but I don't know which one.

Thanks

I suppose you are not really interested in this but if you want I can give you the analytical solution. I did not try answering your question so you can ignore the following if you want to. The solution involves the modified Bessel function of the first kind. It is not an easy equation but it is solvable. If you want just ask and I will post the solution.

coomast
September 7th 2009, 10:59 AM
shawsend

Hi. I'm not sure about the applicability of power series to non-linear equations but I suspect if the equation is well-posed and the resulting expressions converge, then I would have hope the solution may be valid.

So what' wrong with considering:

$y'+y^2=x^2;\quad y(0)=y_0$

and let:

$y(x)=\sum_{n=0}^{\infty} a_nx^n$

Now inserting this power series into the DE I would get:

$\sum_{n=0}^{\infty} n a_n x^{n-1}+\sum_{n=0}^{\infty}\sum_{k=0}^{n} a_k a_{n-k} x^n-x^2=0$

Adjusting the indexing so as to have the same powers on x:

$\sum_{n=0}^{\infty} n a_n x^{n-1}+\sum_{n=1}^{\infty}\left(\sum_{k=0}^{n-1} a_k a_{n-1-k}\right) x^{n-1}-x^2=0$

which gives me a possible solution with:

$a_0:\quad\text{arbitrary}$

$n=1:\quad a_1=-a_0^2$

$n=2:\quad a_2=-a_0 a_1$

$n=3:\quad a_3=\frac{1}{3}\left(1-\sum_{k=0}^2 a_k a_{n-1-k}\right)$

$n\geq 4:\quad a_n=-\frac{1}{n}\sum_{k=0}^{n-1} a_k a_{n-1-k}$

I may need to study further the sum which represents the coefficients for the power series though, and I could get some indication of the validity of my proposed solution above by actually coding it in Mathematica, (say for 100 terms) and then comparing a plot of the solution with the actual solution which Mathematica can compute exactly.
September 7th 2009, 04:01 PM
shawsend

I would think another approach is to convert it to a second order linear DE via a Riccati substitution $y=\frac{u'}{u}$ . Doing that I get:

$u''-x^2u=0$

Now that one is more amendable to power series and since it's second order, likely will involve two power series. Then we can differentiate and make the substitution above to get it in terms of $y(x)$ .
September 8th 2009, 01:59 PM
shawsend

1 Attachment(s)

You guys don't mind if I work on this right? I still left him some work to do and sides, he's probably done with it anyway.

Starting with:

$y'+y^2-x^2=0,\quad y(0)=y_0$

I do the Riccati substitution and obtain:

$u''-x^2u=0,\quad y(x)=\frac{u'(x)}{u(x)}$

Letting $u(x)=\sum_{n=0}^{\infty} a_n x^n$

I get:

$\sum_{n=0}^{\infty} n(n-1) x^{n-2}-\sum_{n=4}^{\infty} a_{n-4} x^{n-2}=0$

Now, after some work, I can express the solution in terms of two power series:

$u(x)=a_0\left[1+\sum_{k=1}^{\infty}\frac{x^{4k}}{\prod_{j=0}^{k-1} (3+4j)(4+4j)}\right]+a_1\left[x+\sum_{k=1}^{\infty}\frac{x^{4k+1}}{\prod_{j=0}^{ k-1} (4+4j)(5+4j)}\right]$

Not sure how that power series is related to the Bessel functions suggested by others and that reported by Mathematica however. Then:

$y(x)=\frac{u'(x)}{u(x)}$

At this point, I'm a little unclear about expressing the initial conditions in u from those in y. I assume all that's needed is the expression $y(0)=\frac{u'(0)}{u(0)}$ and we can use any initial conditions in u to make that ratio true. But I'm not sure about that. So here goes: I considered the IVP:

$y'+y^2+x^2=0,\quad y(0)=2$

That means in term of u, I'd have $2=\frac{u'(0)}{u(0)}=\frac{2}{1}$ . I'm a little unclear on that however. I next coded this in mathematica: I first solved numerically:

$u''-x^2 u=0,\quad u(0)=1, u'(0)=2$

Then let $a_0=1, a_1=2$ in the series expression for u above. Plotted the first 25 terms and then superimposed the numerical results of u with the series results for u. The first plot shows excellent agreement. I next solved numerically:

$y'+y^2-x^2=0,\quad y(0)=2$

Then made the substitution $y(x)=\frac{u'(x)}{u(x)}$ using the power series representation for u. I then superimposed the numerical plots of y with the power series plot of y. The second plot below again shows excellent agreement. This gives me some confidence that the power series is correct -you guys do realize all the secrets of the Universe can be found in non-linear equations right? :)

Here's the Mathematica code if anyone is interested:

Code:

(* solve equation numerically in u *) sol = NDSolve[{Derivative[2][ua][x] - x^2*ua[x] == 0, ua[0] == 1, Derivative[1][ua][0] == 2}, ua, {x, 0, 5}] p1 = Plot[ua[x] /. sol, {x, 0, 3}, PlotStyle -> Blue] (* plot u as the power series *) kmax = 25; a0 = 1; a1 = 2; ub[x_] := a0*(1 + Sum[x^(4*k)/Product[ (3 + 4*j)*(4 + 4*j), {j, 0, k - 1}], {k, 1, kmax}]) + a1*(x + Sum[x^(4*k + 1)/Product[ (4 + 4*j)*(5 + 4*j), {j, 0, k - 1}], {k, 1, kmax}]); p2 = Plot[ub[x], {x, 0, 3}, PlotStyle -> Red] netinu = Show[{p1, p2}] (* solve equation in y numerically *) sol2 = NDSolve[{Derivative[1][yb][x] + yb[x]^2 - x^2 == 0, yb[0] == 2}, yb, {x, 0, 3}] p3 = Plot[yb[x] /. sol2, {x, 0, 3}, PlotStyle -> Blue] (* plot equation in y from power series in u *) ya[x_] = D[ub[x], x]/ub[x]; p4 = Plot[ya[x], {x, 0, 3}, PlotStyle -> Red] netiny = Show[{p4, p3}]