Previous Page Next Page Contents

ode -- the domain of ordinary differential equations

Introduction

ode(eq, y(x)) represents an ordinary differential equation (ODE) for the function y(x).

ode({eq1, eq2, ...}, {y1(x), y2(x), ...}) represents a system of ODEs for the functions y1(x), y2(x) etc.

Call(s)

ode(eq, y(x))
ode({eq <, inits>}, y(x))
ode({eq1, eq2, ... <, inits>}, {y1(x), y2(x), ...})

Parameters

eq, eq1, eq2, ... - equations or arithmetical expressions in the unknown functions and their derivatives with respect to x. An arithmetical expression is regarded as an equation with vanishing right hand side.
y, y1, y2, ... - the unknown functions: identifiers
x - the independent variable: an identifier
inits - the initial or boundary conditions: a sequence of equations

Returns

an object of type ode.

Further Documentation

Section 8.3 of the Tutorial.

Related Functions

numeric::odesolve, numeric::odesolve2, plot::ode

Details

Example 1

In the following, we show how to create and solve a scalar ODE. First, we define the ODE x^2*y'(x) + 3*x*y(x) = sin(x)/x. We use the quote token ' to represent derivatives:

>> eq := ode(x^2*y'(x) + 3*x*y(x) = sin(x)/x, y(x))
                /            sin(x)    2                     \
             ode| 3 x y(x) - ------ + x  diff(y(x), x), y(x) |
                \              x                             /

We get an element of the domain ode which we can now solve:

>> solve(eq)
                             {   C1 + cos(x) }
                             { - ----------- }
                             {        3      }
                             {       x       }
>> delete eq:

Example 2

An initial value problem is defined as a set consisting of the ODE and the initial conditions:

>> ivp := ode({f''(t) + 4*f(t) = sin(2*t), 
               f(0) = a, f'(0) = b}, f(t))
      ode({f(0) = a, D(f)(0) = b, 4 f(t) + diff(f(t), t, t) -
      
         sin(2 t)}, f(t))
>> solve(ivp)
            { sin(2 t)                b sin(2 t)   t cos(2 t) }
            { -------- + a cos(2 t) + ---------- - ---------- }
            {    8                        2            4      }
>> delete ivp:

Example 3

With some restrictions, it is also possible to solve systems of ODEs. First, we define a system:

>> sys := {x'(t) - x(t) + y(t) = 0, y'(t) - x(t) - y(t) = 0}
          {y(t) - x(t) + D(x)(t) = 0, D(y)(t) - y(t) - x(t) = 0}

A call to solve yields the general solution with arbitrary parameters:

>> solution := solve(ode(sys, {x(t), y(t)}))
      {[y(t) = C7 exp((1 + I) t) + C8 exp((1 - I) t),
      
         x(t) = I C7 exp((1 + I) t) - I C8 exp((1 - I) t)]}

To verify the result, we substitute it back into the system sys. However, for the substitution, it is necessary to rewrite the system into a notation using the diff function:

>> eval(subs(rewrite(sys, diff), op(solution)))
                                  {0 = 0}
>> delete sys, solution:

Example 4

In this example, we point out the various return formats of ode's solve facility. First, we solve an ODE with an initial condition. The solution involves a symbolic integral:

>> solve(ode({y'(x) + x*y(x) = cos(x), y(0) = 3}, y(x)))
          {      /    2 \                     2 1/2             }
          {      |   x  |   int(cos(t2) exp(t2 )   , t2 = 0..x) }
          { 3 exp| - -- | + ----------------------------------- }
          {      \   2  /                    2 1/2              }
          {                             exp(x )                 }

The following system is solved incompletely:

>> sys := {x'(t) = -3*y(t)*z(t),
           y'(t) =  3*x(t)*z(t),
           z'(t) = -x(t)*y(t)}:
   solution := solve(ode(sys, {x(t), y(t), z(t)}))
                                  2 1/2                        2 1/2
      {[x(t) = (C10 - C11 + 3 z(t) )   , y(t) = (- C10 - 3 z(t) )   ]
      
                                      2 1/2
         , [x(t) = (C10 - C11 + 3 z(t) )   ,
      
                                 2 1/2
         y(t) = - (- C10 - 3 z(t) )   ],
      
                                      2 1/2
         [x(t) = - (C10 - C11 + 3 z(t) )   ,
      
                               2 1/2
         y(t) = (- C10 - 3 z(t) )   ], [
      
                                     2 1/2
         x(t) = - (C10 - C11 + 3 z(t) )   ,
      
                                 2 1/2
         y(t) = - (- C10 - 3 z(t) )   ]}

In these four different solutions branches, no solution for the unknown function z(t) is provided. In fact, the partial solution above is subject to a further condition on z(t). We substitute the first of our four solutions back into the system sys:

>> eval(subs(rewrite(sys, diff), solution[1]))
      {                                  2 1/2
      { diff(z(t), t) = - (- C10 - 3 z(t) )
      {
      {
      
                            2 1/2    3 z(t) diff(z(t), t)
         (C10 - C11 + 3 z(t) )   , - -------------------- =
                                                    2 1/2
                                     (- C10 - 3 z(t) )
      
                                   2 1/2
         3 z(t) (C10 - C11 + 3 z(t) )   ,
      
           3 z(t) diff(z(t), t)                             2 1/2 }
         ------------------------ = - 3 z(t) (- C10 - 3 z(t) )    }
                            2 1/2                                 }
         (C10 - C11 + 3 z(t) )                                    }

For each solution branch, there remains is exactly one differential equation for z(t) to solve.

>> delete sys, solution:

Example 5

It may happen that MuPAD cannot solve a given equation. In such a case, a symbolic solve command is returned:

>> solve(ode(y'(x) + y(x)^2 = b + a*x, y(x)))
                                                      2
            solve(ode(- b - a x + diff(y(x), x) + y(x) , y(x)))

Example 6

MuPAD's ODE solver contains algebraic algorithms for computing Liouvillian solutions of linear ordinary differential equations over the rational functions as well. These algorithms are based on differential Galois theory. For second order equations, an algorithm similar to the Kovacic algorithm is implemented. However, instead of computing the so-called Riccati-polynomial of a solution, this algorithm computes the solutions directly using formulas. Only when both solutions are algebraic functions, then in three cases it is necessary to compute the minimal polynomial of a solution, again using formulas. Hence, for Kovacic's famous example

y''(x) + (3/(16*x^2) + 2/(9*(x - 1)^2) - 3/(16*x*(x - 1)))*y(x) = 0,

it is possible to compute the minimal polynomial of solution:

>> solve(ode(y''(x) + (3/(16*x^2) + 2/(9*(x - 1)^2)
                       - 3/(16*x*(x - 1)))*y(x), y(x))) 
                    24            8        8
      {C2 RootOf(_Y1   + 2985984 x  (x - 1)  -
      
                  4    8        6         2    16        3
         2799360 x  _Y1  (x - 1)  - 4320 x  _Y1   (x - 1)  -
      
                 5    4               1/2      1/2 7
         165888 x  _Y1  (x - 2) (I x 3    - I 3   )  +
      
                 3    12  1/2               1/2      1/2 4
         5760 I x  _Y1   3    (x - 2) (I x 3    - I 3   ) , _Y1)}

MuPAD may find Liouvillian solutions for higher order equations as well. However, there is no guarantee that all of them are found. The following third order equation can be solved completely:

>> solve(ode(diff(y(x), x, x, x)
           + diff(y(x), x, x)*(2*x^2 - 1)/(2*x^2 - 2*x)
           + diff(y(x), x)*(295*x - 491*x^2 + 196*x^3 - 98)
                            /(196*x^2 - 392*x^3 + 196*x^4)
           + y(x)*(3*x - x^2 - 1)
                   /(196*x^2 - 392*x^3 + 196*x^4), y(x)))
      {
      {                        1/2     1/14
      { C3 (2 x + 2 (x (x - 1))    - 1)     +
      {
      {
      {
      {
      
                        C4
         -------------------------------- +
                             1/2     1/14
         (2 x + 2 (x (x - 1))    - 1)
      
         /
         |                        1/2     1/14
         | C2 (2 x + 2 (x (x - 1))    - 1)
         |
         \
      
         /
         |                        1/2
         | int(exp(-x) (x (x - 1))
         |
         \
      
                             1/2     1/14
         (2 x + 2 (x (x - 1))    - 1)    , x) /
      
                             1/2     1/7
         (2 x + 2 (x (x - 1))    - 1)    -
      
            /                         1/2         \ \ \      }
            |      exp(-x) (x (x - 1))            | | |      }
         int| --------------------------------, x | | | / 28 }
            |                     1/2     1/14    | | |      }
            \ (2 x + 2 (x (x - 1))    - 1)        / / /      }
                                                             }
                                                             }

Background

Changes




Do you have questions or comments?


Copyright © SciFace Software GmbH & Co. KG 2000