6.13 airy’s equation (7.12.00)

6.13.1 Chuck Baker
6.13.2 Chris Eilbeck (8.12.00)
6.13.3 Bill Bauldry (8.12.00)
6.13.4 Dr Francis J. Wright (8.12.00)
6.13.5 Craig B. Watkins (8.12.00)
6.13.6 Willard, Daniel Dr (8.12.00)
6.13.7 Robert Israel (8.12.00)
6.13.8 Francois Boucher (8.12.00)

6.13.1 Chuck Baker

Is it possible to develop a power series solution, in powers of x, for Airy’s equations using Maple?

\[ y'' - x y = 0 \]

6.13.2 Chris Eilbeck (8.12.00)

Use the expansions given in Abramowitz and Stegun, Handbook of Mathematical Functions, Section 10.4.

6.13.3 Bill Bauldry (8.12.00)

A naive approach will work here to get one solution:

>  de := (D@@2)(y)(x) - x*y(x) =0; 
 
                            (2) 
                    de := (D   )(y)(x) - x y(x) = 0 
 
>  Order := 10: # if you want more than 5 terms 
>  dsolve(de, y(x), series); 
 
                                       3                 4 
   y(x) = y(0) + D(y)(0) x + 1/6 y(0) x  + 1/12 D(y)(0) x  + 
 
                     6                  7                 9      10 
         1/180 y(0) x  + 1/504 D(y)(0) x  + 1/12960 y(0) x  + O(x  )
 

But it may be more fun to work with the predefined functions:

>  Order := 6: 
>  taylor(AiryAi(x), x); 
>  taylor(AiryBi(x), x); 
 
          (1/3)          (1/6)                        (1/3) 
         3              3      GAMMA(2/3)            3         3 
   1/3 ---------- - 1/2 ----------------- x + 1/18 ---------- x  - 
       GAMMA(2/3)              Pi                  GAMMA(2/3) 
 
               (1/6) 
              3      GAMMA(2/3)  4      6 
         1/24 ----------------- x  + O(x ) 
                     Pi 
 
          (5/6)          (2/3)                        (5/6) 
         3              3      GAMMA(2/3)            3         3 
   1/3 ---------- + 1/2 ----------------- x + 1/18 ---------- x  + 
       GAMMA(2/3)              Pi                  GAMMA(2/3) 
 
               (2/3) 
              3      GAMMA(2/3)  4      6 
         1/24 ----------------- x  + O(x ) 
                     Pi
 

6.13.4 Dr Francis J. Wright (8.12.00)

From the online help for AiryAi in maple6:

series(AiryAi(x),x,4); 
 
         1/3           1/6                         1/3 
        3             3    GAMMA(2/3)             3        3      4 
 1/3 ---------- - 1/2 --------------- x + 1/18 ---------- x  + O(x ) 
     GAMMA(2/3)             Pi                 GAMMA(2/3)
 

There are probably also harder ways to do it!

6.13.5 Craig B. Watkins (8.12.00)

It’s quite easy. One example I have is:

Order:=20; 
f1:=dsolve({diff(y(x),x)=exp(x*y),y(0)=0},y(x),'type=series'); 
f2:=subs(f1,y(x)); 
f3:=convert(f2,polynom); 
plot(f3,x=-1..1); 
g1:=dsolve({diff(y(x),x$2)+x*y=0,y(0)=0,D(y)(0)=1},y(x),'type=series'); 
g2:=subs(g1,y(x)); 
g3:=convert(g2,polynom); 
plot(g3,x=-1..5,scaling=constrained);
 

This may be more than you wanted (and I use \(y'' + x y = 0\) for this example), but it’s what I like to show undergraduates.

An explanation of what this is used for (somewhat MIT-dependent in terms of subject matter is at

http://web.mit.edu/18.03-esg/www/cws00/maple/seriessol.html

6.13.6 Willard, Daniel Dr (8.12.00)

Try this (Release 5.1):

 restart;#?dsolve[series] 
 with(ODEtools); 
 ode:=diff(y(x),x$2)-x*y(x)=0; 
 Order:=20; 
 dsolve(ode,y(x),'type=series'); 
  # or dsolve({ode, y(0)=A,D(y)(0)=B},y(x),'type=series');
                                                                                    
                                                                                    
 

6.13.7 Robert Israel (8.12.00)

There were lots of replies (including mine) that used the "series" option for dsolve, but Chuck later confirmed to me that he was interested in a symbolic expression for the coefficients. In Maple 6, this can be obtained using the new "formal_series" option:

 
 dsolve((D@@2)(y)(x)-x*y(x)=0,y(x),'formal_series',coeffs=mhypergeom); 
 
                         4 
  {_C[1] x + 1/12 _C[1] x 
 
                                /infinity                               \ 
                                | -----              n  (3 n + 1)       | 
                                |  \            (1/3)  x                | 
           2/9 _C[1] Pi sqrt(3) |   )     ------------------------------| 
                                |  /       n                            | 
                                | -----   3  GAMMA(n + 4/3) GAMMA(n + 1)| 
                                \ n = 2                                 / 
         + --------------------------------------------------------------, 
                                     GAMMA(2/3) 
 
                           3 
        _C[1] + 1/6 _C[1] x 
 
                            /infinity                               \ 
                            | -----                n  (3 n)         | 
                            |  \              (1/3)  x              | 
         + _C[1] GAMMA(2/3) |   )     ------------------------------|} 
                            |  /                    n               | 
                            | -----   GAMMA(n + 1) 3  GAMMA(n + 2/3)| 
                            \ n = 2                                 /
 

6.13.8 Francois Boucher (8.12.00)

Two ideas (with mapleV4) :

> with(powseries): 
> edo:=diff(y(x),x,x)-x*y(x); 
 
                             / 2      \ 
                             |d       | 
                      edo := |--- y(x)| - x y(x) 
                             |  2     | 
                             \dx      / 
 
> sol:=powsolve({edo,y(0)=1,D(y)(0)=0}); 
 
                    sol := proc(powparm)  ...  end 
 
 
 sol(_k); 
 
                              a(_k - 3) 
                             ----------- 
                             _k (_k - 1) 
 
> ak:=rsolve({a(k)=a(k-1)/3/k/(3*k-1),a(0)=1},a(k)); 
 
                              (-k) 
                             9     GAMMA(2/3) 
                  ak := --------------------------- 
                        GAMMA(k + 1) GAMMA(k + 2/3) 
 
> sum(ak*x^(3 *k),k=0..infinity); 
 
                                               3 
                     hypergeom([], [2/3], 1/9 x ) 
 
> restart: 
> with(share): 
#See ?share and ?share,contents for information about the share library 
> readshare(gfun,analysis): 
> with(gfun): 
> edo:=diff(y(x),x,x)-x*y(x); 
 
                             / 2      \ 
                             |d       | 
                      edo := |--- y(x)| - x y(x) 
                             |  2     | 
                             \dx      / 
 
> rec:=diffeqtorec({edo,y(0)=1,D(y)(0)=0},y(x),u(m)); 
 
rec := 
               2 
    {-u(m) + (m  + 5 m + 6) u(m + 3), u(2) = 0, u(0) = 1, u(1) = 0} 
rsolve don't work on this but : 
 
> rec2:=op(1,rec); 
 
                                 2 
               rec2 := -u(m) + (m  + 5 m + 6) u(m + 3) 
 
> rec2:=subs([u(m)=v(n),u(m+3)=v(n+1),m=3*n],rec2); 
 
                                  2 
              rec2 := -v(n) + (9 n  + 15 n + 6) v(n + 1) 
 
> vn:=rsolve({rec2=0,v(0)=1},v(n)); 
 
                              (-n) 
                             9     GAMMA(2/3) 
                  vn := --------------------------- 
                        GAMMA(n + 1) GAMMA(n + 2/3) 
 
> sol:=sum(vn*x^(3*n),n=0..infinity); 
 
                                                  3 
                 sol := hypergeom([], [2/3], 1/9 x )