# Moving on to power series methods

Hi, Sorry for the late post this week. I couldn’t do much of work last week, but I did manage to write some code.

First of all, I managed to get the lie group hint finally in.  The Pull Request can be seen here, https://github.com/sympy/sympy/pull/2359 . Its nice to see the work that you have been doing, for around two months, finally in 😀

I just started work on the power series methods for first order differential equations. Contrary to what I learnt, which is substituting $\sum_{k=0}^n a_{k}x^k$ and finding out a recurrence relation, I found out a much better (straightforward) way to find the power series methods to a given differential equation. Before discussing the algorithm, when does any first order differential equation have a power series solution and what is a power series solution?.

1. Condition for power series solution: Let us take a general first order ODE $P(x, y) + Q(x, y)\frac{dy}{dx} = 0$. It has a power series solution, at a given point when \frac{Q(x, y)}{P(x, y)} is analytic at a given point. Analyticity of an expression at a point is confirmed when, either the expression is infinitely differentiable at a given point, or when the expression has a Power series solution, at a given point. Right now, it is impossible to find out (atleast using SymPy) whether,
any of the above conditions are true. So for now, checking is just being done to see if $\frac{Q(x, y)}{P(x, y)}$ and $\frac{d \frac{Q(x, y)}{P(x, y)}}{dx}$ , exist at the given point.

2. A power series solution is when you can express the solution in the form of a Taylor series. Let us say the solution to a differential equation $P(x, y) + Q(x, y)\frac{dy}{dx} = 0$, is $y = f(x)$, then a power series solution exists at $x0$ when, $f(x)$, can be written as $\sum_{k=0}^n (x - x0)\frac{f^{n}(x0)}{n!}$. As is evident, this can exist only when $f(x)$ is infinitely differentible at $x0$

The algorithm is as follows, lets say $\frac{dy}{dx} = h(x, y)$ and $y(x0) = y0$. Then

1] $\frac{dy}{dx} = F1$
2] $\frac{d^2 y}{dx ^2} = \frac{dF1}{dx} = \frac{\partial F1}{\partial x} + F1\frac{\partial F1}{\partial y} = F2$
3] $\frac{d^2 y}{dx ^2} = \frac{dF2}{dx} = \frac{\partial F2}{\partial x} + F1\frac{\partial F2}{\partial y} = F3$

and so on. Since we know the expressions, we can also find their value at a particular point.

A sample session.

from sympy import *
f = Function("f")
eq = f(x).diff(x) - sin(x*f(x))
pprint(dsolve(eq, hint='1st_power_series', ics={f(2):2, 'terms':1}))
f(x) = (x - 2)⋅sin(4) + 2
pprint(dsolve(eq, hint='1st_power_series', ics={f(2):2, 'terms':3}))
2
(x - 2) ⋅(2⋅cos(4) + 2⋅sin(4)⋅cos(4))    ⎛ 3⎞
f(x) = 2 + (x - 2)⋅sin(4) + ───────────────────────────────────── + O⎝x ⎠
2
pprint(dsolve(eq))
4 ⎛    3       ⎞            2
x ⋅⎝- C₀  + 3⋅C₀⎠        C₀⋅x     ⎛ 6⎞
f(x) = ───────────────── + C₀ + ───── + O⎝x ⎠
24                 2

eq = f(x).diff(x) - x*f(x)
pprint(dsolve(eq, hint='1st_power_series'))

2       4
C₀⋅x    C₀⋅x     ⎛ 6⎞
f(x) = C₀ + ───── + ───── + O⎝x ⎠
2       8


Thats all for now. Cheers!

1. Ok, this looks promising!

I don’t think we need to check the preconditions within sympy, at least not right now.

Are the lines mixed up in the example? The first two pprint belong the the first eq with sin, the second two to the second, no?

2. I compared your first example to Axiom/Fricas:

(1) -> y := operator y

(1) y
Type: BasicOperator
(2) -> f := operator f

(2) f
Type: BasicOperator
(3) -> eq := D(y(x), x) + sin(x*f(x))

,
(3) y (x) + sin(x f(x))

Type: Expression(Integer)
(4) -> seriesSolve(eq=0, y, x=2, [2])
Compiling function %D with type List(UnivariateTaylorSeries(
Expression(Integer),x,2)) -> UnivariateTaylorSeries(Expression(
Integer),x,2)

(4)
2
f(x)cos(2f(x)) 2 f(x) sin(2f(x)) 3
2 – sin(2f(x))(x – 2) – ————– (x – 2) + ————— (x – 2)
2 6
+
3 4
f(x) cos(2f(x)) 4 f(x) sin(2f(x)) 5
————— (x – 2) – ————— (x – 2)
24 120
+
5 6
f(x) cos(2f(x)) 6 f(x) sin(2f(x)) 7
– ————— (x – 2) + ————— (x – 2)
720 5040
+
7 8
f(x) cos(2f(x)) 8 f(x) sin(2f(x)) 9
————— (x – 2) – ————— (x – 2)
40320 362880
+
9
f(x) cos(2f(x)) 10 11
– ————— (x – 2) + O((x – 2) )
3628800
Type: UnivariateTaylorSeries(Expression(Integer),x,2)