ODE

User Functions

These are functions that are imported into the global namespace with from sympy import *. They are intended for user use.

dsolve()

static ode.dsolve(eq, func, hint='default', simplify=True, **kwargs)

Solves any (supported) kind of ordinary differential equation.

Usage

dsolve(eq, f(x), hint) -> Solve ordinary differential equation eq for function f(x), using method hint.

Details

eq can be any supported ordinary differential equation (see
the ode docstring for supported methods). This can either be an Equality, or an expression, which is assumed to be equal to 0.
f(x) is a function of one variable whose derivatives in that
variable make up the ordinary differential equation eq.
hint is the solving method that you want dsolve to use. Use
classify_ode(eq, f(x)) to get all of the possible hints for an ODE. The default hint, ‘default’, will use whatever hint is returned first by classify_ode(). See Hints below for more options that you can use for hint.
simplify enables simplification by odesimp(). See its
docstring for more information. Turn this off, for example, to disable solving of solutions for func or simplification of arbitrary constants. It will still integrate with this hint. Note that the solution may contain more arbitrary constants than the order of the ODE with this option enabled.

Hints

Aside from the various solving methods, there are also some meta-hints that you can pass to dsolve():

“default”:
This uses whatever hint is returned first by classify_ode(). This is the default argument to dsolve().
“all”:

To make dsolve apply all relevant classification hints, use dsolve(ODE, func, hint=”all”). This will return a dictionary of hint:solution terms. If a hint causes dsolve to raise the NotImplementedError, value of that hint’s key will be the exception object raised. The dictionary will also include some special keys:

  • order: The order of the ODE. See also ode_order().
  • best: The simplest hint; what would be returned by “best” below.
  • best_hint: The hint that would produce the solution given by ‘best’. If more than one hint produces the best solution, the first one in the tuple returned by classify_ode() is chosen.
  • default: The solution that would be returned by default. This is the one produced by the hint that appears first in the tuple returned by classify_ode().
“all_Integral”:
This is the same as “all”, except if a hint also has a corresponding “_Integral” hint, it only returns the “_Integral” hint. This is useful if “all” causes dsolve() to hang because of a difficult or impossible integral. This meta-hint will also be much faster than “all”, because integrate() is an expensive routine.
“best”:
To have dsolve() try all methods and return the simplest one. This takes into account whether the solution is solvable in the function, whether it contains any Integral classes (i.e. unevaluatable integrals), and which one is the shortest in size.

See also the classify_ode() docstring for more info on hints, and the ode docstring for a list of all supported hints.

Tips
  • You can declare the derivative of an unknown function this way:
    >>> from sympy import Function, Derivative
    >>> from sympy.abc import x # x is the independent variable
    >>> f = Function("f")(x) # f is a function of x
    >>> # f_ will be the derivative of f with respect to x
    >>> f_ = Derivative(f, x)
    
  • See test_ode.py for many tests, which serves also as a set of examples for how to use dsolve().

  • dsolve always returns an Equality class (except for the case when the hint is “all” or “all_Integral”). If possible, it solves the solution explicitly for the function being solved for. Otherwise, it returns an implicit solution.

  • Arbitrary constants are symbols named C1, C2, and so on.

  • Because all solutions should be mathematically equivalent, some hints may return the exact same result for an ODE. Often, though, two different hints will return the same solution formatted differently. The two should be equivalent. Also note that sometimes the values of the arbitrary constants in two different solutions may not be the same, because one constant may have “absorbed” other constants into it.

  • Do help(ode.ode_hintname) to get help more information on a specific hint, where hintname is the name of a hint without “_Integral”.

Examples

>>> from sympy import Function, dsolve, Eq, Derivative, sin, cos
>>> from sympy.abc import x
>>> f = Function('f')
>>> dsolve(Derivative(f(x),x,x)+9*f(x), f(x))
f(x) == C1*cos(3*x) + C2*sin(3*x)
>>> dsolve(sin(x)*cos(f(x)) + cos(x)*sin(f(x))*f(x).diff(x), f(x),
...     hint='separable')
-log(sin(f(x))**2 - 1)/2 == C1 + log(sin(x)**2 - 1)/2
>>> dsolve(sin(x)*cos(f(x)) + cos(x)*sin(f(x))*f(x).diff(x), f(x),
...     hint='1st_exact')
f(x) == acos(C1/cos(x))
>>> dsolve(sin(x)*cos(f(x)) + cos(x)*sin(f(x))*f(x).diff(x), f(x),
... hint='best')
f(x) == acos(C1/cos(x))
>>> # Note that even though separable is the default, 1st_exact produces
>>> # a simpler result in this case.

classify_ode()

static ode.classify_ode(eq, func, dict=False)

Returns a tuple of possible dsolve() classifications for an ODE.

The tuple is ordered so that first item is the classification that dsolve() uses to solve the ODE by default. In general, classifications at the near the beginning of the list will produce better solutions faster than those near the end, thought there are always exceptions. To make dsolve use a different classification, use dsolve(ODE, func, hint=<classification>). See also the dsolve() docstring for different meta-hints you can use.

If dict is true, classify_ode() will return a dictionary of hint:match expression terms. This is intended for internal use by dsolve(). Note that because dictionaries are ordered arbitrarily, this will most likely not be in the same order as the tuple.

You can get help on different hints by doing help(ode.ode_hintname), where hintname is the name of the hint without “_Integral”.

See sympy.ode.allhints or the sympy.ode docstring for a list of all supported hints that can be returned from classify_ode.

Notes on Hint Names

“_Integral”

If a classification has “_Integral” at the end, it will return the expression with an unevaluated Integral class in it. Note that a hint may do this anyway if integrate() cannot do the integral, though just using an “_Integral” will do so much faster. Indeed, an “_Integral” hint will always be faster than its corresponding hint without “_Integral” because integrate() is an expensive routine. If dsolve() hangs, it is probably because integrate() is hanging on a tough or impossible integral. Try using an “_Integral” hint or “all_Integral” to get it return something.

Note that some hints do not have “_Integral” counterparts. This is because integrate() is not used in solving the ODE for those method. For example, nth order linear homogeneous ODEs with constant coefficients do not require integration to solve, so there is no “nth_linear_homogeneous_constant_coeff_Integrate” hint. You can easily evaluate any unevaluated Integrals in an expression by doing expr.doit().

Ordinals

Some hints contain an ordinal such as “1st_linear”. This is to help differentiate them from other hints, as well as from other methods that may not be implemented yet. If a hint has “nth” in it, such as the “nth_linear” hints, this means that the method used to applies to ODEs of any order.

“indep” and “dep”

Some hints contain the words “indep” or “dep”. These reference the independent variable and the dependent function, respectively. For example, if an ODE is in terms of f(x), then “indep” will refer to x and “dep” will refer to f.

“subs”

If a hints has the word “subs” in it, it means the the ODE is solved by substituting the expression given after the word “subs” for a single dummy variable. This is usually in terms of “indep” and “dep” as above. The substituted expression will be written only in characters allowed for names of Python objects, meaning operators will be spelled out. For example, indep/dep will be written as indep_div_dep.

“coeff”

The word “coeff” in a hint refers to the coefficients of something in the ODE, usually of the derivative terms. See the docstring for the individual methods for more info (help(ode)). This is contrast to “coefficients”, as in “undetermined_coefficients”, which refers to the common name of a method.

“_best”

Methods that have more than one fundamental way to solve will have a hint for each sub-method and a “_best” meta-classification. This will evaluate all hints and return the best, using the same considerations as the normal “best” meta-hint.
Examples
>>> from sympy import Function, classify_ode, Eq
>>> from sympy.abc import x
>>> f = Function('f')
>>> classify_ode(Eq(f(x).diff(x), 0), f(x))
('separable', '1st_linear', '1st_homogeneous_coeff_best',
'1st_homogeneous_coeff_subs_indep_div_dep',
'1st_homogeneous_coeff_subs_dep_div_indep',
'nth_linear_constant_coeff_homogeneous', 'separable_Integral',
'1st_linear_Integral',
'1st_homogeneous_coeff_subs_indep_div_dep_Integral',
'1st_homogeneous_coeff_subs_dep_div_indep_Integral')
>>> classify_ode(f(x).diff(x, 2) + 3*f(x).diff(x) + 2*f(x) - 4, f(x))
('nth_linear_constant_coeff_undetermined_coefficients',
'nth_linear_constant_coeff_variation_of_parameters',
'nth_linear_constant_coeff_variation_of_parameters_Integral')

ode_order()

static ode.ode_order(expr, func)

Returns the order of a given ODE with respect to func.

This function is implemented recursively.

Examples
>>> from sympy import Function, ode_order
>>> from sympy.abc import x
>>> f, g = map(Function, ['f', 'g'])
>>> ode_order(f(x).diff(x, 2) + f(x).diff(x)**2 +
... f(x).diff(x), f(x))
2
>>> ode_order(f(x).diff(x, 2) + g(x).diff(x, 3), f(x))
2
>>> ode_order(f(x).diff(x, 2) + g(x).diff(x, 3), g(x))
3

checkodesol()

static ode.checkodesol(ode, func, sol, order='auto', solve_for_func=True)

Substitutes sol for func in ode and checks that the result is 0.

This only works when func is one function, like f(x). sol can be a single solution or a list of solutions. Either way, each solution must be an Equality instance (e.g., Eq(f(x), C1*cos(x) + C2*sin(x))). If it is a list of solutions, it will return a list of the checkodesol() result for each solution.

It tries the following methods, in order, until it finds zero equivalence:

  1. Substitute the solution for f in the original equation. This only works if the ode is solved for f. It will attempt to solve it first unless solve_for_func == False
  2. Take n derivatives of the solution, where n is the order of ode, and check to see if that is equal to the solution. This only works on exact odes.
  3. Take the 1st, 2nd, ..., nth derivatives of the solution, each time solving for the derivative of f of that order (this will always be possible because f is a linear operator). Then back substitute each derivative into ode in reverse order.

This function returns a tuple. The first item in the tuple is True if the substitution results in 0, and False otherwise. The second item in the tuple is what the substitution results in. It should always be 0 if the first item is True. Note that sometimes this function will False, but with an expression that is identically equal to 0, instead of returning True. This is because simplify() cannot reduce the expression to 0. If an expression returned by this function vanishes identically, then sol really is a solution to ode.

If this function seems to hang, it is probably because of a hard simplification.

To use this function to test, test the first item of the tuple.

Examples
>>> from sympy import Eq, Function, checkodesol, symbols
>>> x, C1 = symbols('x,C1')
>>> f = Function('f')
>>> checkodesol(f(x).diff(x), f(x), Eq(f(x), C1))
(True, 0)
>>> assert checkodesol(f(x).diff(x), f(x), Eq(f(x), C1))[0]
>>> assert not checkodesol(f(x).diff(x), f(x), Eq(f(x), x))[0]
>>> checkodesol(f(x).diff(x, 2), f(x), Eq(f(x), x**2))
(False, 2)

homogeneous_order()

static ode.homogeneous_order(eq, *symbols)

Returns the order n if g is homogeneous and None if it is not homogeneous.

Determines if a function is homogeneous and if so of what order. A function f(x,y,...) is homogeneous of order n if f(t*x,t*y,t*...) == t**n*f(x,y,...).

If the function is of two variables, F(x, y), then f being homogeneous of any order is equivalent to being able to rewrite F(x, y) as G(x/y) or H(y/x). This fact is used to solve 1st order ordinary differential equations whose coefficients are homogeneous of the same order (see the docstrings of ode.ode_1st_homogeneous_coeff_subs_indep_div_dep() and ode.ode_1st_homogeneous_coeff_subs_indep_div_dep()

Symbols can be functions, but every argument of the function must be a symbol, and the arguments of the function that appear in the expression must match those given in the list of symbols. If a declared function appears with different arguments than given in the list of symbols, None is returned.

Examples
>>> from sympy import Function, homogeneous_order, sqrt
>>> from sympy.abc import x, y
>>> f = Function('f')
>>> homogeneous_order(f(x), f(x)) == None
True
>>> homogeneous_order(f(x,y), f(y, x), x, y) == None
True
>>> homogeneous_order(f(x), f(x), x)
1
>>> homogeneous_order(x**2*f(x)/sqrt(x**2+f(x)**2), x, f(x))
2
>>> homogeneous_order(x**2+f(x), x, f(x)) == None
True

Hint Methods

These functions are intended for internal use by dsolve() and others. Nonetheless, they contain useful information in their docstrings on the various ODE solving methods.

1st_exact

static ode.ode_1st_exact(eq, func, order, match)

Solves 1st order exact ordinary differential equations.

A 1st order differential equation is called exact if it is the total differential of a function. That is, the differential equation P(x, y)dx + Q(x, y)dy = 0 is exact if there is some function F(x, y) such that P(x, y) = dF/dx and Q(x, y) = dF/dy (d here refers to the partial derivative). It can be shown that a necessary and sufficient condition for a first order ODE to be exact is that dP/dy = dQ/dx. Then, the solution will be as given below:

>>> from sympy import Function, Eq, Integral, symbols, pprint
>>> x, y, t, x0, y0, C1= symbols('x,y,t,x0,y0,C1')
>>> P, Q, F= map(Function, ['P', 'Q', 'F'])
>>> pprint(Eq(Eq(F(x, y), Integral(P(t, y), (t, x0, x)) +
... Integral(Q(x0, t), (t, y0, y))), C1))
            x                y
            /                /
           |                |
F(x, y) =  |  P(t, y) dt +  |  Q(x0, t) dt = C1
           |                |
          /                /
          x0               y0

Where the first partials of P and Q exist and are continuous in a simply connected region.

A note: SymPy currently has no way to represent inert substitution on an expression, so the hint ‘1st_exact_Integral’ will return an integral with dy. This is supposed to represent the function that you are solving for.

Example
>>> from sympy import Function, dsolve, cos, sin
>>> from sympy.abc import x
>>> f = Function('f')
>>> dsolve(cos(f(x)) - (x*sin(f(x)) - f(x)**2)*f(x).diff(x),
... f(x), hint='1st_exact')
x*cos(f(x)) + f(x)**3/3 == C1
References

# indirect doctest

1st_homogeneous_coeff_best

static ode.ode_1st_homogeneous_coeff_best(eq, func, order, match)

Returns the best solution to an ODE from the two hints ‘1st_homogeneous_coeff_subs_dep_div_indep’ and ‘1st_homogeneous_coeff_subs_indep_div_dep’.

This is as determined by ode_sol_simplicity().

See the ode_1st_homogeneous_coeff_subs_indep_div_dep() and ode_1st_homogeneous_coeff_subs_dep_div_indep() docstrings for more information on these hints. Note that there is no ‘1st_homogeneous_coeff_best_Integral’ hint.

Example

>>> from sympy import Function, dsolve, pprint
>>> from sympy.abc import x
>>> f = Function('f')
>>> pprint(dsolve(2*x*f(x) + (x**2 + f(x)**2)*f(x).diff(x), f(x),
... hint='1st_homogeneous_coeff_best'))
      ___________
     /     2
    /   3*x
   /   ----- + 1 *f(x) = C1
3 /     2
\/     f (x)
References

# indirect doctest

1st_homogeneous_coeff_subs_dep_div_indep

static ode.ode_1st_homogeneous_coeff_subs_dep_div_indep(eq, func, order, match)

Solves a 1st order differential equation with homogeneous coefficients using the substitution u1 = <dependent variable>/<independent variable>.

This is a differential equation P(x, y) + Q(x, y)dy/dx = 0, that P and Q are homogeneous of the same order. A function F(x, y) is homogeneous of order n if F(xt, yt) = t**n*F(x, y). Equivalently, F(x, y) can be rewritten as G(y/x) or H(x/y). See also the docstring of homogeneous_order().

If the coefficients P and Q in the differential equation above are homogeneous functions of the same order, then it can be shown that the substitution y = u1*x (u1 = y/x) will turn the differential equation into an equation separable in the variables x and u. If h(u1) is the function that results from making the substitution u1 = f(x)/x on P(x, f(x)) and g(u2) is the function that results from the substitution on Q(x, f(x)) in the differential equation P(x, f(x)) + Q(x, f(x))*diff(f(x), x) = 0, then the general solution is:

>>> from sympy import Function, dsolve, pprint
>>> from sympy.abc import x
>>> f, g, h = map(Function, ['f', 'g', 'h'])
>>> genform = g(f(x)/x) + h(f(x)/x)*f(x).diff(x)
>>> pprint(genform)
 /f(x)\    /f(x)\ d
g|----| + h|----|*--(f(x))
 \ x  /    \ x  / dx
>>> pprint(dsolve(genform, f(x),
... hint='1st_homogeneous_coeff_subs_dep_div_indep_Integral'))
             f(x)
             ----
              x
               /
              |
              |       -h(u1)
log(C1*x) -   |  ---------------- d(u1) = 0
              |  u1*h(u1) + g(u1)
              |
             /

Where u1*h(u1) + g(u1) != 0 and x != 0.

See also the docstrings of ode_1st_homogeneous_coeff_best() and ode_1st_homogeneous_coeff_subs_indep_div_dep().

Example

>>> from sympy import Function, dsolve
>>> from sympy.abc import x
>>> f = Function('f')
>>> pprint(dsolve(2*x*f(x) + (x**2 + f(x)**2)*f(x).diff(x), f(x),
... hint='1st_homogeneous_coeff_subs_dep_div_indep'))
        ________________
       /           3
      /  3*f(x)   f (x)
x*   /   ------ + -----  = C1
  3 /      x         3
  \/                x
References

# indirect doctest

1st_homogeneous_coeff_subs_indep_div_dep

static ode.ode_1st_homogeneous_coeff_subs_indep_div_dep(eq, func, order, match)

Solves a 1st order differential equation with homogeneous coefficients using the substitution u2 = <independent variable>/<dependent variable>.

This is a differential equation P(x, y) + Q(x, y)dy/dx = 0, that P and Q are homogeneous of the same order. A function F(x, y) is homogeneous of order n if F(xt, yt) = t**n*F(x, y). Equivalently, F(x, y) can be rewritten as G(y/x) or H(x/y). See also the docstring of homogeneous_order().

If the coefficients P and Q in the differential equation above are homogeneous functions of the same order, then it can be shown that the substitution x = u2*y (u2 = x/y) will turn the differential equation into an equation separable in the variables y and u2. If h(u2) is the function that results from making the substitution u2 = x/f(x) on P(x, f(x)) and g(u2) is the function that results from the substitution on Q(x, f(x)) in the differential equation P(x, f(x)) + Q(x, f(x))*diff(f(x), x) = 0, then the general solution is:

>>> from sympy import Function, dsolve, pprint
>>> from sympy.abc import x
>>> f, g, h = map(Function, ['f', 'g', 'h'])
>>> genform = g(x/f(x)) + h(x/f(x))*f(x).diff(x)
>>> pprint(genform)
 / x  \    / x  \ d
g|----| + h|----|*--(f(x))
 \f(x)/    \f(x)/ dx
>>> pprint(dsolve(genform, f(x),
... hint='1st_homogeneous_coeff_subs_indep_div_dep_Integral'))
                 x
                ----
                f(x)
                  /
                 |
                 |        g(u2)
log(C1*f(x)) -   |  ----------------- d(u2) = 0
                 |  -u2*g(u2) - h(u2)
                 |
                /

Where u2*g(u2) + h(u2) != 0 and f(x) != 0.

See also the docstrings of ode_1st_homogeneous_coeff_best() and ode_1st_homogeneous_coeff_subs_dep_div_indep().

Example
>>> from sympy import Function, pprint
>>> from sympy.abc import x
>>> f = Function('f')
>>> pprint(dsolve(2*x*f(x) + (x**2 + f(x)**2)*f(x).diff(x), f(x),
... hint='1st_homogeneous_coeff_subs_indep_div_dep'))
      ___________
     /     2
    /   3*x
   /   ----- + 1 *f(x) = C1
3 /     2
\/     f (x)
References

# indirect doctest

1st_linear

static ode.ode_1st_linear(eq, func, order, match)

Solves 1st order linear differential equations.

These are differential equations of the form dy/dx _ P(x)*y = Q(x). These kinds of differential equations can be solved in a general way. The integrating factor exp(Integral(P(x), x)) will turn the equation into a separable equation. The general solution is:

>>> from sympy import Function, dsolve, Eq, pprint, diff, sin
>>> from sympy.abc import x
>>> f, P, Q = map(Function, ['f', 'P', 'Q'])
>>> genform = Eq(f(x).diff(x) + P(x)*f(x), Q(x))
>>> pprint(genform)
            d
P(x)*f(x) + --(f(x)) = Q(x)
            dx
>>> pprint(dsolve(genform, f(x), hint='1st_linear_Integral'))
       /       /                   \
       |      |                    |
       |      |         /          |     /
       |      |        |           |    |
       |      |        | P(x) dx   |  - | P(x) dx
       |      |        |           |    |
       |      |       /            |   /
f(x) = |C1 +  | Q(x)*e           dx|*e
       |      |                    |
       \     /                     /
Example
>>> f = Function('f')
>>> pprint(dsolve(Eq(x*diff(f(x), x) - f(x), x**2*sin(x)),
... f(x), '1st_linear'))
f(x) = x*(C1 - cos(x))
References

# indirect doctest

Bernoulli

static ode.ode_Bernoulli(eq, func, order, match)

Solves Bernoulli differential equations.

These are equations of the form dy/dx + P(x)*y = Q(x)*y**n, n != 1. The substitution w = 1/y**(1-n) will transform an equation of this form into one that is linear (see the docstring of ode_1st_linear()). The general solution is:

>>> from sympy import Function, dsolve, Eq, pprint
>>> from sympy.abc import x, n
>>> f, P, Q = map(Function, ['f', 'P', 'Q'])
>>> genform = Eq(f(x).diff(x) + P(x)*f(x), Q(x)*f(x)**n)
>>> pprint(genform)
            d                n
P(x)*f(x) + --(f(x)) = Q(x)*f (x)
            dx
>>> pprint(dsolve(genform, f(x), hint='Bernoulli_Integral')) 
                                                                               1
                                                                              ----
                                                                             1 - n
       //                /                            \                     \
       ||               |                             |                     |
       ||               |                  /          |             /       |
       ||               |                 |           |            |        |
       ||               |        (1 - n)* | P(x) dx   |  (-1 + n)* | P(x) dx|
       ||               |                 |           |            |        |
       ||               |                /            |           /         |
f(x) = ||C1 + (-1 + n)* | -Q(x)*e                   dx|*e                   |
       ||               |                             |                     |
       \\               /                            /                     /

Note that when n = 1, then the equation is separable (see the docstring of ode_separable()).

>>> pprint(dsolve(Eq(f(x).diff(x) + P(x)*f(x), Q(x)*f(x)), f(x),
... hint='separable_Integral'))
 f(x)
   /
  |                /
  |  1            |
  |  - dy = C1 +  | (-P(x) + Q(x)) dx
  |  y            |
  |              /
 /
Example
>>> from sympy import Function, dsolve, Eq, pprint, log
>>> from sympy.abc import x
>>> f = Function('f')
>>> pprint(dsolve(Eq(x*f(x).diff(x) + f(x), log(x)*f(x)**2),
... f(x), hint='Bernoulli'))
                1
f(x) = -------------------
         /     log(x)   1\
       x*|C1 + ------ + -|
         \       x      x/
References

# indirect doctest

Liouville

static ode.ode_Liouville(eq, func, order, match)

Solves 2nd order Liouville differential equations.

The general form of a Liouville ODE is d^2y/dx^2 + g(y)*(dy/dx)**2 + h(x)*dy/dx. The general solution is:

>>> from sympy import Function, dsolve, Eq, pprint, diff
>>> from sympy.abc import x
>>> f, g, h = map(Function, ['f', 'g', 'h'])
>>> genform = Eq(diff(f(x),x,x) + g(f(x))*diff(f(x),x)**2 +
... h(x)*diff(f(x),x), 0)
>>> pprint(genform)
                2                    2
        d                d          d
g(f(x))*--(f(x))  + h(x)*--(f(x)) + ---(f(x)) = 0
        dx               dx           2
                                    dx
>>> pprint(dsolve(genform, f(x), hint='Liouville_Integral'))
                                  f(x)
          /                     /
         |                     |
         |     /               |     /
         |    |                |    |
         |  - | h(x) dx        |    | g(y) dy
         |    |                |    |
         |   /                 |   /
C1 + C2* | e            dx +   |  e           dy = 0
         |                     |
        /                     /

Example

>>> from sympy import Function, dsolve, Eq, pprint
>>> from sympy.abc import x
>>> f = Function('f')
>>> pprint(dsolve(diff(f(x), x, x) + diff(f(x), x)**2/f(x) +
... diff(f(x), x)/x, f(x), hint='Liouville'))
           ___   ________________           ___   ________________
[f(x) = -\/ 2 *\/ C1 + C2*log(x) , f(x) = \/ 2 *\/ C1 + C2*log(x) ]
References

# indirect doctest

nth_linear_constant_coeff_homogeneous

static ode.ode_nth_linear_constant_coeff_homogeneous(eq, func, order, match, returns='sol')

Solves an nth order linear homogeneous differential equation with constant coefficients.

This is an equation of the form a_n*f(x)^(n) + a_(n-1)*f(x)^(n-1) + ... + a1*f’(x) + a0*f(x) = 0

These equations can be solved in a general manner, by taking the roots of the characteristic equation a_n*m**n + a_(n-1)*m**(n-1) + ... + a1*m + a0 = 0. The solution will then be the sum of Cn*x**i*exp(r*x) terms, for each where Cn is an arbitrary constant, r is a root of the characteristic equation and i is is one of each from 0 to the multiplicity of the root - 1 (for example, a root 3 of multiplicity 2 would create the terms C1*exp(3*x) + C2*x*exp(3*x)). The exponential is usually expanded for complex roots using Euler’s equation exp(I*x) = cos(x) + I*sin(x). Complex roots always come in conjugate pars in polynomials with real coefficients, so the two roots will be represented (after simplifying the constants) as exp(a*x)*(C1*cos(b*x) + C2*sin(b*x)).

If SymPy cannot find exact roots to the characteristic equation, a RootOf instance will be return in its stead.

>>> from sympy import Function, dsolve, Eq
>>> from sympy.abc import x
>>> f = Function('f')
>>> dsolve(f(x).diff(x, 5) + 10*f(x).diff(x) - 2*f(x), f(x),
... hint='nth_linear_constant_coeff_homogeneous')
... 
f(x) == C1*exp(x*RootOf(_x**5 + 10*_x - 2, 0)) +     C2*exp(x*RootOf(_x**5 + 10*_x - 2, 1)) +     C3*exp(x*RootOf(_x**5 + 10*_x - 2, 2)) +     C4*exp(x*RootOf(_x**5 + 10*_x - 2, 3)) +     C5*exp(x*RootOf(_x**5 + 10*_x - 2, 4))

Note that because this method does not involve integration, there is no ‘nth_linear_constant_coeff_homogeneous_Integral’ hint.

The following is for internal use:

  • returns = ‘sol’ returns the solution to the ODE.
  • returns = ‘list’ returns a list of linearly independent solutions, for use with non homogeneous solution methods like variation of parameters and undetermined coefficients. Note that, though the solutions should be linearly independent, this function does not explicitly check that. You can do “assert simplify(wronskian(sollist)) != 0” to check for linear independence. Also, “assert len(sollist) == order” will need to pass.
  • returns = ‘both’, return a dictionary {‘sol’:solution to ODE, ‘list’: list of linearly independent solutions}.
Example
>>> from sympy import Function, dsolve, pprint
>>> from sympy.abc import x
>>> f = Function('f')
>>> pprint(dsolve(f(x).diff(x, 4) + 2*f(x).diff(x, 3) -
... 2*f(x).diff(x, 2) - 6*f(x).diff(x) + 5*f(x), f(x),
... hint='nth_linear_constant_coeff_homogeneous'))
                    x                            -2*x
f(x) = (C1 + C2*x)*e  + (C3*cos(x) + C4*sin(x))*e
References

# indirect doctest

nth_linear_constant_coeff_undetermined_coefficients

static ode.ode_nth_linear_constant_coeff_undetermined_coefficients(eq, func, order, match)

Solves an nth order linear differential equation with constant coefficients using the method of undetermined coefficients.

This method works on differential equations of the form a_n*f(x)^(n) + a_(n-1)*f(x)^(n-1) + ... + a1*f’(x) + a0*f(x) = P(x), where P(x) is a function that has a finite number of linearly independent derivatives.

Functions that fit this requirement are finite sums functions of the form a*x**i*exp(b*x)*sin(c*x + d) or a*x**i*exp(b*x)*cos(c*x + d), where i is a non-negative integer and a, b, c, and d are constants. For example any polynomial in x, functions like x**2*exp(2*x), x*sin(x), and exp(x)*cos(x) can all be used. Products of sin’s and cos’s have a finite number of derivatives, because they can be expanded into sin(a*x) and cos(b*x) terms. However, SymPy currently cannot do that expansion, so you will need to manually rewrite the expression in terms of the above to use this method. So, for example, you will need to manually convert sin(x)**2 into (1 + cos(2*x))/2 to properly apply the method of undetermined coefficients on it.

This method works by creating a trial function from the expression and all of its linear independent derivatives and substituting them into the original ODE. The coefficients for each term will be a system of linear equations, which are be solved for and substituted, giving the solution. If any of the trial functions are linearly dependent on the solution to the homogeneous equation, they are multiplied by sufficient x to make them linearly independent.

Example
>>> from sympy import Function, dsolve, pprint, exp, cos
>>> from sympy.abc import x
>>> f = Function('f')
>>> pprint(dsolve(f(x).diff(x, 2) + 2*f(x).diff(x) + f(x) -
... 4*exp(-x)*x**2 + cos(2*x), f(x),
... hint='nth_linear_constant_coeff_undetermined_coefficients'))
       /             4\
       |            x |  -x   4*sin(2*x)   3*cos(2*x)
f(x) = |C1 + C2*x + --|*e   - ---------- + ----------
       \            3 /           25           25
References

# indirect doctest

nth_linear_constant_coeff_variation_of_parameters

static ode.ode_nth_linear_constant_coeff_variation_of_parameters(eq, func, order, match)

Solves an nth order linear differential equation with constant coefficients using the method of undetermined coefficients.

This method works on any differential equations of the form f(x)^(n) + a_(n-1)*f(x)^(n-1) + ... + a1*f’(x) + a0*f(x) = P(x).

This method works by assuming that the particular solution takes the form Sum(c_i(x)*y_i(x), (x, 1, n)), where y_i is the ith solution to the homogeneous equation. The solution is then solved using Wronskian’s and Cramer’s Rule. The particular solution is given by Sum(Integral(W_i(x)/W(x), x)*y_i(x), (x, 1, n)), where W(x) is the Wronskian of the fundamental system (the system of n linearly independent solutions to the homogeneous equation), and W_i(x) is the Wronskian of the fundamental system with the ith column replaced with [0, 0, ..., 0, P(x)].

This method is general enough to solve any nth order inhomogeneous linear differential equation with constant coefficients, but sometimes SymPy cannot simplify the Wronskian well enough to integrate it. If this method hangs, try using the ‘nth_linear_constant_coeff_variation_of_parameters_Integral’ hint and simplifying the integrals manually. Also, prefer using ‘nth_linear_constant_coeff_undetermined_coefficients’ when it applies, because it doesn’t use integration, making it faster and more reliable.

Warning, using simplify=False with ‘nth_linear_constant_coeff_variation_of_parameters’ in dsolve() may cause it to hang, because it will not attempt to simplify the Wronskian before integrating. It is recommended that you only use simplify=False with ‘nth_linear_constant_coeff_variation_of_parameters_Integral’ for this method, especially if the solution to the homogeneous equation has trigonometric functions in it.

Example
>>> from sympy import Function, dsolve, pprint, exp, log
>>> from sympy.abc import x
>>> f = Function('f')
>>> pprint(dsolve(f(x).diff(x, 3) - 3*f(x).diff(x, 2) +
... 3*f(x).diff(x) - f(x) - exp(x)*log(x), f(x),
... hint='nth_linear_constant_coeff_variation_of_parameters'))
       /                2    3 /log(x)   11\\  x
f(x) = |C1 + C2*x + C3*x  + x *|------ - --||*e
       \                       \  6      36//
References

# indirect doctest

separable

static ode.ode_separable(eq, func, order, match)

Solves separable 1st order differential equations.

This is any differential equation that can be written as P(y)*dy/dx = Q(x). The solution can then just be found by rearranging terms and integrating: Integral(P(y), y) = Integral(Q(x), x). This hint uses separatevars() as its back end, so if a separable equation is not caught by this solver, it is most likely the fault of that function. separatevars() is smart enough to do most expansion and factoring necessary to convert a separable equation F(x, y) into the proper form P(x)*Q(y). The general solution is:

>>> from sympy import Function, dsolve, Eq, pprint
>>> from sympy.abc import x
>>> a, b, c, d, f = map(Function, ['a', 'b', 'c', 'd', 'f'])
>>> genform = Eq(a(x)*b(f(x))*f(x).diff(x), c(x)*d(f(x)))
>>> pprint(genform)
             d
a(x)*b(f(x))*--(f(x)) = c(x)*d(f(x))
             dx
>>> pprint(dsolve(genform, f(x), hint='separable_Integral'))
     f(x)
   /                  /
  |                  |
  |  b(y)            | c(x)
  |  ---- dy = C1 +  | ---- dx
  |  d(y)            | a(x)
  |                  |
 /                  /

Example

>>> from sympy import Function, dsolve, Eq
>>> from sympy.abc import x
>>> f = Function('f')
>>> pprint(dsolve(Eq(f(x)*f(x).diff(x) + x, 3*x*f(x)**2), f(x),
... hint='separable'))
   /   2       \         2
log\3*f (x) - 1/        x
---------------- = C1 + --
       6                2
Reference
  • M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 52

# indirect doctest

Information on the ode module

This module contains dsolve() and different helper functions that it uses.

dsolve() solves ordinary differential equations. See the docstring on the various functions for their uses. Note that partial differential equations support is in pde.py. Note that ode_hint() functions have docstrings describing their various methods, but they are intended for internal use. Use dsolve(ode, func, hint=hint) to solve an ode using a specific hint. See also the docstring on dsolve().

Functions in this module

These are the user functions in this module:

  • dsolve() - Solves ODEs.
  • classify_ode() - Classifies ODEs into possible hints for dsolve().
  • checkodesol() - Checks if an equation is the solution to an ODE.
  • ode_order() - Returns the order (degree) of an ODE.
  • homogeneous_order() - Returns the homogeneous order of an expression.

These are the non-solver helper functions that are for internal use. The user should use the various options to dsolve() to obtain the functionality provided by these functions:

  • odesimp() - Does all forms of ODE simplification.
  • ode_sol_simplicity() - A key function for comparing solutions by simplicity.
  • constantsimp() - Simplifies arbitrary constants.
  • constant_renumber() - Renumber arbitrary constants
  • _handle_Integral() - Evaluate unevaluated Integrals.

See also the docstrings of these functions.

Solving methods currently implemented

The following methods are implemented for solving ordinary differential equations. See the docstrings of the various ode_hint() functions for more information on each (run help(ode)):

  • 1st order separable differential equations
  • 1st order differential equations whose coefficients or dx and dy are functions homogeneous of the same order.
  • 1st order exact differential equations.
  • 1st order linear differential equations
  • 1st order Bernoulli differential equations.
  • 2nd order Liouville differential equations.
  • nth order linear homogeneous differential equation with constant coefficients.
  • nth order linear inhomogeneous differential equation with constant coefficients using the method of undetermined coefficients.
  • nth order linear inhomogeneous differential equation with constant coefficients using the method of variation of parameters.

Philosophy behind this module

This module is designed to make it easy to add new ODE solving methods without having to mess with the solving code for other methods. The idea is that there is a classify_ode() function, which takes in an ODE and tells you what hints, if any, will solve the ODE. It does this without attempting to solve the ODE, so it is fast. Each solving method is a hint, and it has its own function, named ode_hint. That function takes in the ODE and any match expression gathered by classify_ode and returns a solved result. If this result has any integrals in it, the ode_hint function will return an unevaluated Integral class. dsolve(), which is the user wrapper function around all of this, will then call odesimp() on the result, which, among other things, will attempt to solve the equation for the dependent variable (the function we are solving for), simplify the arbitrary constants in the expression, and evaluate any integrals, if the hint allows it.

How to add new solution methods

If you have an ODE that you want dsolve() to be able to solve, try to avoid adding special case code here. Instead, try finding a general method that will solve your ODE, as well as others. This way, the ode module will become more robust, and unhindered by special case hacks. WolphramAlpha and Maple’s DETools[odeadvisor] function are two resources you can use to classify a specific ODE. It is also better for a method to work with an nth order ODE instead of only with specific orders, if possible.

To add a new method, there are a few things that you need to do. First, you need a hint name for your method. Try to name your hint so that it is unambiguous with all other methods, including ones that may not be implemented yet. If your method uses integrals, also include a “hint_Integral” hint. If there is more than one way to solve ODEs with your method, include a hint for each one, as well as a “hint_best” hint. Your ode_hint_best() function should choose the best using min with ode_sol_simplicity as the key argument. See ode_1st_homogeneous_coeff_best(), for example. The function that uses your method will be called ode_hint(), so the hint must only use characters that are allowed in a Python function name (alphanumeric characters and the underscore ‘_’ character). Include a function for every hint, except for “_Integral” hints (dsolve() takes care of those automatically). Hint names should be all lowercase, unless a word is commonly capitalized (such as Integral or Bernoulli). If you have a hint that you do not want to run with “all_Integral” that doesn’t have an “_Integral” counterpart (such as a best hint that would defeat the purpose of “all_Integral”), you will need to remove it manually in the dsolve() code. See also the classify_ode() docstring for guidelines on writing a hint name.

Determine in general how the solutions returned by your method compare with other methods that can potentially solve the same ODEs. Then, put your hints in the allhints tuple in the order that they should be called. The ordering of this tuple determines which hints are default. Note that exceptions are ok, because it is easy for the user to choose individual hints with dsolve(). In general, “_Integral” variants should go at the end of the list, and “_best” variants should go before the various hints they apply to. For example, the “undetermined_coefficients” hint comes before the “variation_of_parameters” hint because, even though variation of parameters is more general than undetermined coefficients, undetermined coefficients generally returns cleaner results for the ODEs that it can solve than variation of parameters does, and it does not require integration, so it is much faster.

Next, you need to have a match expression or a function that matches the type of the ODE, which you should put in classify_ode() (if the match function is more than just a few lines, like _undetermined_coefficients_match(), it should go outside of classify_ode()). It should match the ODE without solving for it as much as possible, so that classify_ode() remains fast and is not hindered by bugs in solving code. Be sure to consider corner cases. For example, if your solution method involves dividing by something, make sure you exclude the case where that division will be 0.

In most cases, the matching of the ODE will also give you the various parts that you need to solve it. You should put that in a dictionary (.match() will do this for you), and add that as matching_hints[‘hint’] = matchdict in the relevant part of classify_ode. classify_ode will then send this to dsolve(), which will send it to your function as the match argument. Your function should be named ode_hint(eq, func, order, match). If you need to send more information, put it in the match dictionary. For example, if you had to substitute in a dummy variable in classify_ode to match the ODE, you will need to pass it to your function using the match dict to access it. You can access the independent variable using func.args[0], and the dependent variable (the function you are trying to solve for) as func.func. If, while trying to solve the ODE, you find that you cannot, raise NotImplementedError. dsolve() will catch this error with the “all” meta-hint, rather than causing the whole routine to fail.

Add a docstring to your function that describes the method employed. Like with anything else in SymPy, you will need to add a doctest to the docstring, in addition to real tests in test_ode.py. Try to maintain consistency with the other hint functions’ docstrings. Add your method to the list at the top of this docstring. Also, add your method to ode.txt in the docs/src directory, so that the Sphinx docs will pull its docstring into the main SymPy documentation. Be sure to make the Sphinx documentation by running “make html” from within the doc directory to verify that the docstring formats correctly.

If your solution method involves integrating, use C.Integral() instead of integrate(). This allows the user to bypass hard/slow integration by using the “_Integral” variant of your hint. In most cases, calling .doit() will integrate your solution. If this is not the case, you will need to write special code in _handle_Integral(). Arbitrary constants should be symbols named C1, C2, and so on. All solution methods should return an equality instance. If you need an arbitrary number of arbitrary constants, you can use constants = numbered_symbols(prefix=’C’, cls=Symbol, start=1). If it is possible to solve for the dependent function in a general way, do so. Otherwise, do as best as you can, but do not call solve in your ode_hint() function. odesimp() will attempt to solve the solution for you, so you do not need to do that. Lastly, if your ODE has a common simplification that can be applied to your solutions, you can add a special case in odesimp() for it. For example, solutions returned from the “1st_homogeneous_coeff” hints often have many log() terms, so odesimp() calls logcombine() on them (it also helps to write the arbitrary constant as log(C1) instead of C1 in this case). Also consider common ways that you can rearrange your solution to have constantsimp() take better advantage of it. It is better to put simplification in odesimp() than in your method, because it can then be turned off with the simplify flag in dsolve(). If you have any extraneous simplification in your function, be sure to only run it using “if match.get(‘simplify’, True):”, especially if it can be slow or if it can reduce the domain of the solution.

Finally, as with every contribution to SymPy, your method will need to be tested. Add a test for each method in test_ode.py. Follow the conventions there, i.e., test the solver using dsolve(eq, f(x), hint=your_hint), and also test the solution using checkodesol (you can put these in a separate tests and skip/XFAIL if it runs too slow/doesn’t work). Be sure to call your hint specifically in dsolve, that way the test won’t be broken simply by the introduction of another matching hint. If your method works for higher order (>1) ODEs, you will need to run sol = constant_renumber(sol, ‘C’, 1, order), for each solution, where order is the order of the ODE. This is because constant_renumber renumbers the arbitrary constants by printing order, which is platform dependent. Try to test every corner case of your solver, including a range of orders if it is a nth order solver, but if your solver is slow, auch as if it involves hard integration, try to keep the test run time down.

Feel free to refactor existing hints to avoid duplicating code or creating inconsistencies. If you can show that your method exactly duplicates an existing method, including in the simplicity and speed of obtaining the solutions, then you can remove the old, less general method. The existing code is tested extensively in test_ode.py, so if anything is broken, one of those tests will surely fail.