Towards the end

Hi, I finally managed to do some amount of work for the past ten days, and I’m happy to day that all that I had wanted to do this GSoC has been pushed in and I’m waiting for comments on my final Pull Request. These were the changes that I had made.

1. Making stuff similar to the series function
While playing around with the code, I found out that, doing series(eq, n=terms), gives the series expansion upto O(x^n) rather than the nth term. This makes the code look a lot neater too. Also the default power is 6, so I changes that too. For instance, take this sample SymPy session

from sympy import dsolve, series
from sympy.abc import x, y
f = Function("f")
pprint(series(exp(x))
         2    3    4     5        
        x    x    x     x     ⎛ 6⎞
1 + x + ── + ── + ── + ─── + O⎝x ⎠
        2    6    24   120       
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       

I also added support for power series solutions of homogeneous differential equations at ordinary points, and regular singular points. A homogeneous second differential equation is of the form P(x)\frac{d^2f}{dx^2} + Q(x)\frac{df}{dx} + R(x). (Any DE mentioned below is homogeneous second order unless otherwise specified, because typing DE is easier than typing a homogeneous second order differential equation).
a] A point x0 is said to be ordinary, if \frac{Q(x)}{P(x)} and \frac{R(x)}{P(x)} are analytic at the point.
b] It is said to be regular singular, if (x - x0)\frac{Q(x)}{P(x)} and (x - x0)^2\frac{Q(x)}{P(x)} are analytic. For simplicity in the series expansions, assumptions are made such that P(x), Q(x) and R(x) are polynomials. Okay, now for a bit more detail

1. Ordinary points : A DE has a power series solution, at x0. This can be found by substituting \sum_{k=0}^n a_{n}x^n in the DE, and equating the nth coefficient in order to obtain a recurrence relation. However this is not as trivial as it sounds and was one of the more tougher things I had to do.

Take the general case of P(x)\frac{d^2f}{dx^2} + Q(x)\frac{df}{dx} + R(x). Substituting y as \sum_{n = 0}^\infty a_{n}x^n in the DE. (One has to expand each of the terms P(x), Q(x) and R(x)), make transformations such that for each term the power of x is n. The catch here is that the starting point of the summation also changes, so the necessary terms shouls be stripped off such that the starting point is same. This is better explained in one of these tutorials, http://tutorial.math.lamar.edu/Classes/DE/SeriesSolutions.aspx , so you could read that after you’ve read the rest of the blogpost.

A small sample session.

from sympy import *
from sympy.abc import x, y
eq = (1 + x**2)*(f(x).diff(x, 2)) + 2*x*(f(x).diff(x)) -2*f(x)
pprint(eq)
                          2               
    d          ⎛ 2    ⎞  d                
2⋅x⋅──(f(x)) + ⎝x  + 1⎠⋅───(f(x)) - 2⋅f(x)
    dx                    2               
                        dx
pprint(dsolve(eq))
                 ⎛   4         ⎞        
                 ⎜  x     2    ⎟    ⎛ 6⎞
f(x) = C₁⋅x + C₀⋅⎜- ── + x  + 1⎟ + O⎝x ⎠
                 ⎝  3          ⎠
eq = f(x).diff(x, 2) + x*(f(x).diff(x)) + f(x)
pprint(eq)
                      2      
  d                  d       
x⋅──(f(x)) + f(x) + ───(f(x))
  dx                  2      
                    dx
pprint(dsolve(eq))
            ⎛   2    ⎞      ⎛ 4    2    ⎞        
            ⎜  x     ⎟      ⎜x    x     ⎟    ⎛ 6⎞
f(x) = C₁⋅x⋅⎜- ── + 1⎟ + C₀⋅⎜── - ── + 1⎟ + O⎝x ⎠
            ⎝  3     ⎠      ⎝8    2     ⎠  

2. Regular singular points :
1. Try expressing (x - x0)P(x) and ((x - x0)^{2})Q(x) as power series
solutions about x0. Find p0 and q0 which are the constants of the
power series expansions.
2. Solve the indicial equation f(m) = m(m - 1) + m*p0 + q0, to obtain the
roots m1 and m2 of the indicial equation.
3. If m1 - m2 is a non integer there exists two series solutions. If
m1 = m2, there exists only one solution. If m1 - m2 is an integer,
then the existence of one solution is confirmed. The other solution may
or may not exist.

The power series solution is of the form x^{m}\sum_{n=0}^\infty a_{n}x^n. The
coefficients are determined by the following recurrence relation.
a_{n} = -\frac{\sum_{k=0}^{n-1} q_{n-k} + (m + k)p_{n-k}}{f(m + n)}. For the case
in which m1 - m2 is an integer, it can be seen from the recurrence relation
that for the lower root m, when n equals the difference of both the
roots, the denominator becomes zero. So if the numerator is not equal to zero,
a second series solution exists. (Oh and by the way I just copy – pasted the docstring for the
last part, because I was lazy to write the whole thing again :P)

    eq = x**2*(f(x).diff(x, 2)) - 3*x*(f(x).diff(x)) + (4*x + 4)*f(x)  # One solution
    pprint(dsolve(eq))
                ⎛      3                 ⎞        
              2 ⎜  16⋅x       2          ⎟    ⎛ 6⎞
   f(x) = C₀⋅x ⋅⎜- ───── + 4⋅x  - 4⋅x + 1⎟ + O⎝x ⎠
                ⎝    9                   ⎠  
    eq = x**2*(f(x).diff(x, 2)) - x**2*(f(x).diff(x)) + (
        x**2 - 2)*f(x)  # Two solutions
             ⎛    6      5    4    2        ⎞                                    
             ⎜   x    3⋅x    x    x    x    ⎟                                    
          C₁⋅⎜- ─── - ──── - ── + ── + ─ + 1⎟         ⎛   3    2        ⎞        
             ⎝  720    80    8    2    2    ⎠       2 ⎜  x    x    x    ⎟    ⎛ 6⎞
   f(x) = ─────────────────────────────────── + C₀⋅x ⋅⎜- ── + ── + ─ + 1⎟ + O⎝x ⎠
                           x                          ⎝  60   20   2    ⎠  
    

Anyway, that sums up my official GSoC work. Now that I just need to address comments made on my PR which is at, https://github.com/sympy/sympy/pull/2395 . And I just realised that I have my exams coming up the next week, and its been a while since I touched my books (and gone to classes), and thats not a good thing. So have to catch up there. Seems like my whole life is about catching up.

Advertisements

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: