Solving Ordinary Differential Equations with Maple...


Getting started...

Type maple to begin.
oscar.Princeton.EDU% maple
    |\^/|     Maple V Release 5 (WMI Campus Wide License)
._|\|   |/|_. Copyright (c) 1981-1997 by Waterloo Maple Inc. All rights
 \  MAPLE  /  reserved. Maple and Maple V are registered trademarks of
 <____ ____>  Waterloo Maple Inc.
      |       Type ? for help.

First Order Linear ODEs


# y'(x) + (1/x) * y(x) = 1 

> sol1 := dsolve(diff(y(x), x) + y(x) / x = 1, y(x));

                                                           _C1
                                    sol1 := y(x) = 1/2 x + ---
                                                            x

#This is a general solution

# Let's apply an initial condition y(1) = -1 and find the constant _C1

> dsolve({diff(y(x), x) + y(x) / x =1 , y(1) = -1} , y(x));

                                      y(x) = 1/2 x - 3/2 1/x
# Thus _C1 = -3/2 


# Another example 
# y'(x) = 8 * x^3 * y^2
 
> dsolve(diff(y(x), x) = 8 * x^3 * y(x)^2, y(x));
                                                     1
                                        y(x) = - ----------
                                                    4
                                                 2 x  - _C1


Second Order ODEs

# y''(x) + y(x) = 0
# Notice that $2 indicates a second derivative

> eq1 := diff(y(x), x$2) + y(x) = 0;
                                          / 2      \
                                          |d       |
                                   eq1 := |--- y(x)| + y(x) = 0
                                          |  2     |
                                          \dx      /

#Let's find a general solution 
> sol1 := dsolve(eq1, y(x));
                              sol1 := y(x) = _C1 sin(x) + _C2 cos(x)

# Note two constants: _C1 and _C2 

# Let's try another example...
# y''(x) + y(x) = x

> eq2 := diff(y(x), x$2) + y(x) = x;
                                          / 2      \
                                          |d       |
                                   eq2 := |--- y(x)| + y(x) = x
                                          |  2     |
                                          \dx      /

> sol2 := dsolve(eq2, y(x));
                            sol2 := y(x) = x + _C1 sin(x) + _C2 cos(x)

> 

Using Initial Conditions...


# Consider the following equation with initial conditions:
# y'' + y = sin(t)
# y(0) = 0   and  y'(0) = 1

> eq5 := dsolve({diff(y(t), t$2) + y(t) = sin(t), y(0) = 0, D(y)(0) = 1}, y(t));         

                                      3
              eq5 := y(t) = 1/2 sin(t)  + (1/2 cos(t) sin(t) - 1/2 t) cos(t) + sin(t)

# Notice that there are no arbitrary constants in this solution
# Function rhs() is used to obtain the right hand side of eq5 in the example below.
# Type this to plot the graph


> plot(rhs(eq5), t=0..4*Pi); 

  |                                                                                               
  +                                                                  AAAA                         
  |                                                                 AA   A                        
  +                                                                A      A                       
 4+                                                                A       A                      
  |                                                               A         A                     
  +                                                              A          A                     
  +                                                             A            A                    
  |                                                            AA            A                    
  +                                                            A              A                   
  +                                                            A              A                   
  |                                                           A               AA                  
 2+               AAAAAA                                     A                 A                  
  +            AAAA    AAA                                   A                 A                  
  |          AAA          AA                                 A                  A                 
  +        AA               A                               A                   A                 
  +      AA                  AA                            AA                    A                
  |    AA                     A                            A                     A                
  + AAA                        A                           A                     A                
  |AA            2              4A             6          A   8             10    A       12      
  *--+--+--+--+--+--+--+--+--+--+*-+--+--+--+--+--+--+--+*-+--+--+--+--+--+--+--+-*--+--+--+--+-- 
 0+                               A                      A                        AA              
  |                                A                     A                         A              
  +                                A                    A                          A              
  +                                 AA                 A                           A              
  |                                  A                 A                            A             
  +                                   A               A                             A             
-2+                                    A             AA                             AA            
  |                                     A            A                               A            
  +                                     AA          A                                A            
  +                                      AA        A                                 AA           
  |                                        AA    AA                                   A           
  +                                         AAAAAA                                    A           
  +                                                                                    A          
  |                                                                                    A          
-4+                                                                                     A         
  +                                                                                     A         
  |                                                                                     A         
  +                                                                                      A        
  |                                                                                      A        
  +                                                                                       A       
  +                                                                                        A      
  |                                                                                         A     
-6+                                                                                         AA    
  +                                                                                          AAAA 
  |                                                                                               
                                                                                                  

Solving Systems of ODEs


# Consider the system of ODEs
 x' = 4x + 2y
 y' = 3x + 3y

> ex1 := {diff(x1(t), t) = 4 * x1(t) + 2 * x2(t), 
>         diff(x2(t), t) = 3 * x1(t) + 3 * x2(t)}; 
                        d                             d
                ex1 := {-- x2(t) = 3 x1(t) + 3 x2(t), -- x1(t) = 4 x1(t) + 2 x2(t)}
                        dt                            dt

> 
> dsolve(ex1, {x1(t), x2(t)} );

   {x1(t) = 2/5 _C1 exp(t) + 3/5 _C1 exp(6 t) + 2/5 _C2 exp(6 t) - 2/5 _C2 exp(t),

    x2(t) = 3/5 _C1 exp(6 t) - 3/5 _C1 exp(t) + 3/5 _C2 exp(t) + 2/5 _C2 exp(6 t) }

 

Numerical Solutions of ODEs

# Some problems do not have solutions in closed-form.
# Consider the motion of a simple pendulum with unit natural frequency:
# y'' + sin(y) = 0

> pend := diff(y(x), x$2) + sin(y(x)) = 0;
                                        / 2      \
                                        |d       |
                                pend := |--- y(x)| + sin(y(x)) = 0
                                        |  2     |
                                        \dx      /

> 
# Let's try to solve it...

> dsolve(pend, y(x));


     y(x)                                             y(x)
    /                                                /
   |               1                                |               1
   |      -------------------- d_a - x - _C2 = 0, - |      -------------------- d_a - x - _C2 = 0
   |                       1/2                      |                       1/2
  /       (2 cos(_a) + _C1)                        /       (2 cos(_a) + _C1)



# Maple found two implicit solutions, expressing 'x' as a function of 'y' in terms of
# elliptic integrals.

# Let's find the numerical solution to the pendulum equations.
# Suppose that  y(0) = 0  and y'(0) = 1.

> sol := dsolve( {pend, y(0) = 0, D(y)(0) = 1}, y(x), type=numeric); 

                                  sol := proc(rkf45_x)  ...  end

# Note that the solution is returned as a procedure rkf45_x, displayed in abbreviated form.

# We will use function odeplot() to get a sketch of the result.
# This function is given here with tree arguments:
#   - the name of the procedure returned by dsolve(),
#   - a list with the names of the independent and dependent variables, and
#   - a range for the independent variable.


> plots[odeplot](sol, [x, y(x)], 0..3*Pi);  

                                                                                                  
    |                                                                                             
    |             AAAAA                                                            AAAA           
   1+            AA    A                                                         AA    AA         
    |           A       A                                                       A       A         
    |          AA        A                                                     AA        A        
    +          A          A                                                    A         AA       
    |         A            A                                                  A           A       
    |        A             A                                                 AA            A      
    +        A              A                                                A              A     
    |       AA              A                                               AA              A     
    |       A               AA                                              A               AA    
    +       A                A                                             A                 A    
    |      A                 A                                             A                  A   
    |     A                   A                                            A                      
    +     A                   A                                           A                       
    |     A                    A                                          A                       
 0.5+    A                     A                                         A                        
    |    A                      A                                        A                        
    |    A                      A                                        A                        
    +   A                       A                                       A                         
    |   A                        A                                      A                         
    |   A                        A                                      A                         
    +  A                         A                                     A                          
    |  A                          A                                    A                          
    | A                           A                                    A                          
    + A                           AA                                  A                           
    | A                            A                                  A                           
    |AA                            A                                 A                            
    +A                             A                                 A                            
    |A                              A                                A                            
    |A                              A                               A                             
  --*---+---+---+--+---+---+---+---+-*+---+---+---+---+--+---+---+--*+---+--+---+---+---+---+--+- 
   0|                  2             A    4                  6      A           8                 
    +                                A                              A                             
    |                                A                             A                              
    |                                 A                            A                              
    +                                 A                           A                               
    |                                 A                           A                               
    |                                  A                          A                               
    +                                  A                         A                                
    |                                  A                         A                                
    |                                   A                        A                                
    +                                   A                       A                                 
    |                                    A                      A                                 
    |                                    A                      A                                 
-0.5+                                    A                     A                                  
    |                                     A                    A                                  
    |                                     A                   A                                   
    +                                      A                  A                                   
    |                                      A                 A                                    
    +                                      A                 A                                    
    |                                       A               AA                                    
    |                                       A               A                                     
    +                                        A             A                                      
    |                                        A             A                                      
    |                                         A            A                                      
    +                                          A          A                                       
    |                                          AA        A                                        
    |                                           AA      AA                                        
  -1+                                            AA   AA                                          
    |                                              AAAA                                           
    |                                                                                             
                                                                                                  
                                                                                                  
> 
> 

Numerical Solutions of Systems of ODEs

# y' = x
# x' = -y - 0.1 * x * (y^2-1)
# suppose y(0) = 1/2   and x(0) = 0

> eq9 := {diff(y(t), t) = x(t),
>         diff(x(t), t) = -y(t) - 0.1 * (y(t)^2 - 1) * x(t),
>         y(0) = 1/2,  x(0) = 0};
               d               d                         2
       eq9 := {-- y(t) = x(t), -- x(t) = -y(t) - .1 (y(t)  - 1) x(t), y(0) = 1/2, x(0) = 0}
               dt              dt

>
> sol5 := dsolve(eq9, {x(t), y(t)}, type = numeric);

                                 sol5 := proc(rkf45_x)  ...  end

> 
> 
> plots[odeplot](sol5, [[t, x(t)], [t, y(t)]], 0..20);

    |                                                                                            
    +                                                                                  BB        
    |                                                                           AAA    B B       
   1+                                                                           A A   B  B       
    |                                                                          A   A B    B      
    +                                                                          A   A B    B      
    |                                                     BBBB                A     AB     B     
    +                                               AAA   B  B                A     *      B     
    |                                              AA  A BB   B              AA     *      B     
    +                                              A   A B    BB             A      *A     B     
    |                          BBB                A     AB     B             A      BA      B    
    +                   AAA   B   B               A     BA     B            A      B A      B    
    |                  A   AA B   BB             A      BA      B           A      B A      B    
 0.5*B                A     AB     B             A     BBA      B           A      B  A      B   
    | BB              A     *A      B            A     B A      B           A      B  A      B   
    +  B             AA     BA      B           A      B  A      B         A      B   A          
    |   B            A     B  A     BB          A      B  A      B         A      B   A          
    +   BB          A      B  A      B          A     B   A      B         A      B   AA         
    |    B          A     BB  A      B         A      B    A     B         A      B    A         
    +    B         AA     B    A      B        A      B    A      B       AA     B     A         
    |     B        A      B    A      B        A     B     A      B       A      B     A         
    +     B        A      B    AA     B       AA     B     A      B       A      B     A         
    |      B      A      B      A      B      A      B      A      B      A      B      A        
  --**---+-*-+---+*---+--*+----+*--+---*+---+-*-+----*---+--*-+---+*---+--*+---+*---+---*----+-- 
   0|A     B     AA     B 5      A     B      A10   B       A      B  15 A      B       A   20   
    +A      B    A      B        A     B     A      B       A      B     A      B       A        
    | A      B   A      B        A      B    A      B        A     B     A      B       A        
    + A      B  A      B         A      B   AA     B         A      B    A      B        A       
    |  A     B  A      B          A      B  A      B         A      B   A      B         A       
    +  A      BA      B           A      B  A      B          A     B   A      B         A       
    |   A     BA      B            A     B  A      B          A      B  A      B         A       
    +   A     AB     B             A      BA      B           A      B  A      B         A       
-0.5+    AA  A BB   BB             AA     BA      B           A      B A      BB          A      
    |      AAA  B  BB               A     *B     B             A     B A      B           A      
    +            BBB                A     AB     B             A      BA      B           A      
    |                                A   AA B   BB             AA     *       B           A      
    +                                 AAA   B   B               A     *      B             A     
    |                                  A     BBBB               A     AB     B             A     
    +                                         B                  A    AB    BB             A     
    |                                                            A   A BB   B               A    
    +                                                             AAA   B  BB               A    
    |                                                              A    B  B                A    
  -1+                                                                    BBB                A    
    |                                                                                        A   
    |                                                                                            
                                                                                                 
> 

Displaying Graphs in Separate Windows


> plotsetup(x11,xterm);
> plot(x^2, x=-2..2);

These examples were taken from [Robertson 1996].
Yefim Shuf | CS Department | Princeton University
Last modified: Wed Oct 21 13:20:24 EDT 1998