ODE¶
Note
For a beginner-friendly guide focused on solving ODEs, refer to Solve an Ordinary Differential Equation (ODE) Algebraically.
User Functions¶
These are functions that are imported into the global namespace with from
sympy import *
. These functions (unlike Hint Functions, below) are
intended for use by ordinary users of SymPy.
- sympy.solvers.ode.dsolve(
- eq,
- func=None,
- hint='default',
- simplify=True,
- ics=None,
- xi=None,
- eta=None,
- x0=0,
- n=6,
- **kwargs,
Solves any (supported) kind of ordinary differential equation and system of ordinary differential equations.
For Single Ordinary Differential Equation
It is classified under this when number of equation in
eq
is one. Usagedsolve(eq, f(x), hint)
-> Solve ordinary differential equationeq
for functionf(x)
, using methodhint
.Details
eq
can be any supported ordinary differential equation (see theode
docstring for supported methods). This can either be anEquality
, or an expression, which is assumed to be equal to0
.f(x)
is a function of one variable whose derivatives in thatvariable make up the ordinary differential equation
eq
. In many cases it is not necessary to provide this; it will be autodetected (and an error raised if it could not be detected).hint
is the solving method that you want dsolve to use. Useclassify_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 byclassify_ode()
. See Hints below for more options that you can use for hint.simplify
enables simplification byodesimp()
. See its docstring for more information. Turn this off, for example, to disable solving of solutions forfunc
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.xi
andeta
are the infinitesimal functions of an ordinarydifferential equation. They are the infinitesimals of the Lie group of point transformations for which the differential equation is invariant. The user can specify values for the infinitesimals. If nothing is specified,
xi
andeta
are calculated usinginfinitesimals()
with the help of various heuristics.ics
is the set of initial/boundary conditions for the differential equation.It should be given in the form of
{f(x0): x1, f(x).diff(x).subs(x, x2): x3}
and so on. For power series solutions, if no initial conditions are specifiedf(0)
is assumed to beC0
and the power series solution is calculated about 0.x0
is the point about which the power series solution of a differentialequation is to be evaluated.
n
gives the exponent of the dependent variable up to which the power seriessolution of a differential equation is to be evaluated.
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 todsolve()
.all
:To make
dsolve()
apply all relevant classification hints, usedsolve(ODE, func, hint="all")
. This will return a dictionary ofhint:solution
terms. If a hint causes dsolve to raise theNotImplementedError
, 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 alsoode_order()
indeutils.py
.best
: The simplest hint; what would be returned bybest
below.best_hint
: The hint that would produce the solution given bybest
. If more than one hint produces the best solution, the first one in the tuple returned byclassify_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 byclassify_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 ifall
causesdsolve()
to hang because of a difficult or impossible integral. This meta-hint will also be much faster thanall
, becauseintegrate()
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 theode
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 usedsolve()
.dsolve()
always returns anEquality
class (except for the case when the hint isall
orall_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
.
For System Of Ordinary Differential Equations
- Usage
dsolve(eq, func)
-> Solve a system of ordinary differential equationseq
forfunc
being list of functions including \(x(t)\), \(y(t)\), \(z(t)\) where number of functions in the list depends upon the number of equations provided ineq
.
Details
eq
can be any supported system of ordinary differential equations This can either be anEquality
, or an expression, which is assumed to be equal to0
.func
holdsx(t)
andy(t)
being functions of one variable which together with some of their derivatives make up the system of ordinary differential equationeq
. It is not necessary to provide this; it will be autodetected (and an error raised if it could not be detected).Hints
The hints are formed by parameters returned by classify_sysode, combining them give hints name used later for forming method name.
Examples
>>> from sympy import Function, dsolve, Eq, Derivative, sin, cos, symbols >>> from sympy.abc import x >>> f = Function('f') >>> dsolve(Derivative(f(x), x, x) + 9*f(x), f(x)) Eq(f(x), C1*sin(3*x) + C2*cos(3*x))
>>> eq = sin(x)*cos(f(x)) + cos(x)*sin(f(x))*f(x).diff(x) >>> dsolve(eq, hint='1st_exact') [Eq(f(x), -acos(C1/cos(x)) + 2*pi), Eq(f(x), acos(C1/cos(x)))] >>> dsolve(eq, hint='almost_linear') [Eq(f(x), -acos(C1/cos(x)) + 2*pi), Eq(f(x), acos(C1/cos(x)))] >>> t = symbols('t') >>> x, y = symbols('x, y', cls=Function) >>> eq = (Eq(Derivative(x(t),t), 12*t*x(t) + 8*y(t)), Eq(Derivative(y(t),t), 21*x(t) + 7*t*y(t))) >>> dsolve(eq) [Eq(x(t), C1*x0(t) + C2*x0(t)*Integral(8*exp(Integral(7*t, t))*exp(Integral(12*t, t))/x0(t)**2, t)), Eq(y(t), C1*y0(t) + C2*(y0(t)*Integral(8*exp(Integral(7*t, t))*exp(Integral(12*t, t))/x0(t)**2, t) + exp(Integral(7*t, t))*exp(Integral(12*t, t))/x0(t)))] >>> eq = (Eq(Derivative(x(t),t),x(t)*y(t)*sin(t)), Eq(Derivative(y(t),t),y(t)**2*sin(t))) >>> dsolve(eq) {Eq(x(t), -exp(C1)/(C2*exp(C1) - cos(t))), Eq(y(t), -1/(C1 - cos(t)))}
- sympy.solvers.ode.systems.dsolve_system(
- eqs,
- funcs=None,
- t=None,
- ics=None,
- doit=False,
- simplify=True,
Solves any(supported) system of Ordinary Differential Equations
- Parameters:
eqs : List
system of ODEs to be solved
funcs : List or None
List of dependent variables that make up the system of ODEs
t : Symbol or None
Independent variable in the system of ODEs
ics : Dict or None
Set of initial boundary/conditions for the system of ODEs
doit : Boolean
Evaluate the solutions if True. Default value is True. Can be set to false if the integral evaluation takes too much time and/or is not required.
simplify: Boolean
Simplify the solutions for the systems. Default value is True. Can be set to false if simplification takes too much time and/or is not required.
- Returns:
List of List of Equations
- Raises:
NotImplementedError
When the system of ODEs is not solvable by this function.
ValueError
When the parameters passed are not in the required form.
Explanation
This function takes a system of ODEs as an input, determines if the it is solvable by this function, and returns the solution if found any.
This function can handle: 1. Linear, First Order, Constant coefficient homogeneous system of ODEs 2. Linear, First Order, Constant coefficient non-homogeneous system of ODEs 3. Linear, First Order, non-constant coefficient homogeneous system of ODEs 4. Linear, First Order, non-constant coefficient non-homogeneous system of ODEs 5. Any implicit system which can be divided into system of ODEs which is of the above 4 forms 6. Any higher order linear system of ODEs that can be reduced to one of the 5 forms of systems described above.
The types of systems described above are not limited by the number of equations, i.e. this function can solve the above types irrespective of the number of equations in the system passed. But, the bigger the system, the more time it will take to solve the system.
This function returns a list of solutions. Each solution is a list of equations where LHS is the dependent variable and RHS is an expression in terms of the independent variable.
Among the non constant coefficient types, not all the systems are solvable by this function. Only those which have either a coefficient matrix with a commutative antiderivative or those systems which may be divided further so that the divided systems may have coefficient matrix with commutative antiderivative.
Examples
>>> from sympy import symbols, Eq, Function >>> from sympy.solvers.ode.systems import dsolve_system >>> f, g = symbols("f g", cls=Function) >>> x = symbols("x")
>>> eqs = [Eq(f(x).diff(x), g(x)), Eq(g(x).diff(x), f(x))] >>> dsolve_system(eqs) [[Eq(f(x), -C1*exp(-x) + C2*exp(x)), Eq(g(x), C1*exp(-x) + C2*exp(x))]]
You can also pass the initial conditions for the system of ODEs:
>>> dsolve_system(eqs, ics={f(0): 1, g(0): 0}) [[Eq(f(x), exp(x)/2 + exp(-x)/2), Eq(g(x), exp(x)/2 - exp(-x)/2)]]
Optionally, you can pass the dependent variables and the independent variable for which the system is to be solved:
>>> funcs = [f(x), g(x)] >>> dsolve_system(eqs, funcs=funcs, t=x) [[Eq(f(x), -C1*exp(-x) + C2*exp(x)), Eq(g(x), C1*exp(-x) + C2*exp(x))]]
Lets look at an implicit system of ODEs:
>>> eqs = [Eq(f(x).diff(x)**2, g(x)**2), Eq(g(x).diff(x), g(x))] >>> dsolve_system(eqs) [[Eq(f(x), C1 - C2*exp(x)), Eq(g(x), C2*exp(x))], [Eq(f(x), C1 + C2*exp(x)), Eq(g(x), C2*exp(x))]]
- sympy.solvers.ode.classify_ode(
- eq,
- func=None,
- dict=False,
- ics=None,
- *,
- prep=True,
- xi=None,
- eta=None,
- n=None,
- **kwargs,
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 makedsolve()
use a different classification, usedsolve(ODE, func, hint=<classification>)
. See also thedsolve()
docstring for different meta-hints you can use.If
dict
is true,classify_ode()
will return a dictionary ofhint:match
expression terms. This is intended for internal use bydsolve()
. 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 executing
help(ode.ode_hintname)
, wherehintname
is the name of the hint without_Integral
.See
allhints
or theode
docstring for a list of all supported hints that can be returned fromclassify_ode()
.Notes
These are remarks on hint names.
_Integral
If a classification has
_Integral
at the end, it will return the expression with an unevaluatedIntegral
class in it. Note that a hint may do this anyway ifintegrate()
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
becauseintegrate()
is an expensive routine. Ifdsolve()
hangs, it is probably becauseintegrate()
is hanging on a tough or impossible integral. Try using an_Integral
hint orall_Integral
to get it return something.Note that some hints do not have
_Integral
counterparts. This is becauseintegrate()
is not used in solving the ODE for those method. For example, \(n\)th order linear homogeneous ODEs with constant coefficients do not require integration to solve, so there is nonth_linear_homogeneous_constant_coeff_Integrate
hint. You can easily evaluate any unevaluatedIntegral
s in an expression by doingexpr.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 hasnth
in it, such as thenth_linear
hints, this means that the method used to applies to ODEs of any order.indep
anddep
Some hints contain the words
indep
ordep
. These reference the independent variable and the dependent function, respectively. For example, if an ODE is in terms of \(f(x)\), thenindep
will refer to \(x\) anddep
will refer to \(f\).subs
If a hints has the word
subs
in it, it means that the ODE is solved by substituting the expression given after the wordsubs
for a single dummy variable. This is usually in terms ofindep
anddep
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 asindep_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 tocoefficients
, as inundetermined_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 normalbest
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)) ('nth_algebraic', 'separable', '1st_exact', '1st_linear', 'Bernoulli', '1st_homogeneous_coeff_best', '1st_homogeneous_coeff_subs_indep_div_dep', '1st_homogeneous_coeff_subs_dep_div_indep', '1st_power_series', 'lie_group', 'nth_linear_constant_coeff_homogeneous', 'nth_linear_euler_eq_homogeneous', 'nth_algebraic_Integral', 'separable_Integral', '1st_exact_Integral', '1st_linear_Integral', 'Bernoulli_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) ('factorable', 'nth_linear_constant_coeff_undetermined_coefficients', 'nth_linear_constant_coeff_variation_of_parameters', 'nth_linear_constant_coeff_variation_of_parameters_Integral')
- sympy.solvers.ode.checkodesol(
- ode,
- sol,
- func=None,
- order='auto',
- solve_for_func=True,
Substitutes
sol
intoode
and checks that the result is0
.This works when
func
is one function, like \(f(x)\) or a list of functions like \([f(x), g(x)]\) when \(ode\) is a system of ODEs.sol
can be a single solution or a list of solutions. Each solution may be anEquality
that the solution satisfies, e.g.Eq(f(x), C1), Eq(f(x) + C1, 0)
; or simply anExpr
, e.g.f(x) - C1
. In most cases it will not be necessary to explicitly identify the function, but if the function cannot be inferred from the original equation it can be supplied through thefunc
argument.If a sequence of solutions is passed, the same sort of container will be used to return the result for each solution.
It tries the following methods, in order, until it finds zero equivalence:
Substitute the solution for \(f\) in the original equation. This only works if
ode
is solved for \(f\). It will attempt to solve it first unlesssolve_for_func == False
.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.Take the 1st, 2nd, …, \(n\)th 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 in0
, andFalse
otherwise. The second item in the tuple is what the substitution results in. It should always be0
if the first item isTrue
. Sometimes this function will returnFalse
even when an expression is identically equal to0
. This happens whensimplify()
does not reduce the expression to0
. If an expression returned by this function vanishes identically, thensol
really is a solution to theode
.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, ... Derivative, exp) >>> x, C1, C2 = symbols('x,C1,C2') >>> f, g = symbols('f g', cls=Function) >>> checkodesol(f(x).diff(x), Eq(f(x), C1)) (True, 0) >>> assert checkodesol(f(x).diff(x), C1)[0] >>> assert not checkodesol(f(x).diff(x), x)[0] >>> checkodesol(f(x).diff(x, 2), x**2) (False, 2)
>>> eqs = [Eq(Derivative(f(x), x), f(x)), Eq(Derivative(g(x), x), g(x))] >>> sol = [Eq(f(x), C1*exp(x)), Eq(g(x), C2*exp(x))] >>> checkodesol(eqs, sol) (True, [0, 0])
- sympy.solvers.ode.homogeneous_order(eq, *symbols)[source]¶
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, \cdots)\) is homogeneous of order \(n\) if \(f(t x, t y, \cdots) = t^n f(x, y, \cdots)\).
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
HomogeneousCoeffSubsDepDivIndep
andHomogeneousCoeffSubsIndepDivDep
).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)) is None True >>> homogeneous_order(f(x,y), f(y, x), x, y) is 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)) is None True
- sympy.solvers.ode.infinitesimals(
- eq,
- func=None,
- order=None,
- hint='default',
- match=None,
The infinitesimal functions of an ordinary differential equation, \(\xi(x,y)\) and \(\eta(x,y)\), are the infinitesimals of the Lie group of point transformations for which the differential equation is invariant. So, the ODE \(y'=f(x,y)\) would admit a Lie group \(x^*=X(x,y;\varepsilon)=x+\varepsilon\xi(x,y)\), \(y^*=Y(x,y;\varepsilon)=y+\varepsilon\eta(x,y)\) such that \((y^*)'=f(x^*, y^*)\). A change of coordinates, to \(r(x,y)\) and \(s(x,y)\), can be performed so this Lie group becomes the translation group, \(r^*=r\) and \(s^*=s+\varepsilon\). They are tangents to the coordinate curves of the new system.
Consider the transformation \((x, y) \to (X, Y)\) such that the differential equation remains invariant. \(\xi\) and \(\eta\) are the tangents to the transformed coordinates \(X\) and \(Y\), at \(\varepsilon=0\).
\[\left(\frac{\partial X(x,y;\varepsilon)}{\partial\varepsilon }\right)|_{\varepsilon=0} = \xi, \left(\frac{\partial Y(x,y;\varepsilon)}{\partial\varepsilon }\right)|_{\varepsilon=0} = \eta,\]The infinitesimals can be found by solving the following PDE:
>>> from sympy import Function, Eq, pprint >>> from sympy.abc import x, y >>> xi, eta, h = map(Function, ['xi', 'eta', 'h']) >>> h = h(x, y) # dy/dx = h >>> eta = eta(x, y) >>> xi = xi(x, y) >>> genform = Eq(eta.diff(x) + (eta.diff(y) - xi.diff(x))*h ... - (xi.diff(y))*h**2 - xi*(h.diff(x)) - eta*(h.diff(y)), 0) >>> pprint(genform) /d d \ d 2 d d d |--(eta(x, y)) - --(xi(x, y))|*h(x, y) - eta(x, y)*--(h(x, y)) - h (x, y)*--(xi(x, y)) - xi(x, y)*--(h(x, y)) + --(eta(x, y)) = 0 \dy dx / dy dy dx dx
Solving the above mentioned PDE is not trivial, and can be solved only by making intelligent assumptions for \(\xi\) and \(\eta\) (heuristics). Once an infinitesimal is found, the attempt to find more heuristics stops. This is done to optimise the speed of solving the differential equation. If a list of all the infinitesimals is needed,
hint
should be flagged asall
, which gives the complete list of infinitesimals. If the infinitesimals for a particular heuristic needs to be found, it can be passed as a flag tohint
.Examples
>>> from sympy import Function >>> from sympy.solvers.ode.lie_group import infinitesimals >>> from sympy.abc import x >>> f = Function('f') >>> eq = f(x).diff(x) - x**2*f(x) >>> infinitesimals(eq) [{eta(x, f(x)): exp(x**3/3), xi(x, f(x)): 0}]
References
Solving differential equations by Symmetry Groups, John Starrett, pp. 1 - pp. 14
- sympy.solvers.ode.checkinfsol(
- eq,
- infinitesimals,
- func=None,
- order=None,
This function is used to check if the given infinitesimals are the actual infinitesimals of the given first order differential equation. This method is specific to the Lie Group Solver of ODEs.
As of now, it simply checks, by substituting the infinitesimals in the partial differential equation.
\[\frac{\partial \eta}{\partial x} + \left(\frac{\partial \eta}{\partial y} - \frac{\partial \xi}{\partial x}\right)*h - \frac{\partial \xi}{\partial y}*h^{2} - \xi\frac{\partial h}{\partial x} - \eta\frac{\partial h}{\partial y} = 0\]where \(\eta\), and \(\xi\) are the infinitesimals and \(h(x,y) = \frac{dy}{dx}\)
The infinitesimals should be given in the form of a list of dicts
[{xi(x, y): inf, eta(x, y): inf}]
, corresponding to the output of the function infinitesimals. It returns a list of values of the form[(True/False, sol)]
wheresol
is the value obtained after substituting the infinitesimals in the PDE. If it isTrue
, thensol
would be 0.
- sympy.solvers.ode.constantsimp(expr, constants)[source]¶
Simplifies an expression with arbitrary constants in it.
This function is written specifically to work with
dsolve()
, and is not intended for general use.Simplification is done by “absorbing” the arbitrary constants into other arbitrary constants, numbers, and symbols that they are not independent of.
The symbols must all have the same name with numbers after it, for example,
C1
,C2
,C3
. Thesymbolname
here would be ‘C
’, thestartnumber
would be 1, and theendnumber
would be 3. If the arbitrary constants are independent of the variablex
, then the independent symbol would bex
. There is no need to specify the dependent function, such asf(x)
, because it already has the independent symbol,x
, in it.Because terms are “absorbed” into arbitrary constants and because constants are renumbered after simplifying, the arbitrary constants in expr are not necessarily equal to the ones of the same name in the returned result.
If two or more arbitrary constants are added, multiplied, or raised to the power of each other, they are first absorbed together into a single arbitrary constant. Then the new constant is combined into other terms if necessary.
Absorption of constants is done with limited assistance:
terms of
Add
s are collected to try join constants so \(e^x (C_1 \cos(x) + C_2 \cos(x))\) will simplify to \(e^x C_1 \cos(x)\);powers with exponents that are
Add
s are expanded so \(e^{C_1 + x}\) will be simplified to \(C_1 e^x\).
Use
constant_renumber()
to renumber constants after simplification or else arbitrary numbers on constants may appear, e.g. \(C_1 + C_3 x\).In rare cases, a single constant can be “simplified” into two constants. Every differential equation solution should have as many arbitrary constants as the order of the differential equation. The result here will be technically correct, but it may, for example, have \(C_1\) and \(C_2\) in an expression, when \(C_1\) is actually equal to \(C_2\). Use your discretion in such situations, and also take advantage of the ability to use hints in
dsolve()
.Examples
>>> from sympy import symbols >>> from sympy.solvers.ode.ode import constantsimp >>> C1, C2, C3, x, y = symbols('C1, C2, C3, x, y') >>> constantsimp(2*C1*x, {C1, C2, C3}) C1*x >>> constantsimp(C1 + 2 + x, {C1, C2, C3}) C1 + x >>> constantsimp(C1*C2 + 2 + C2 + C3*x, {C1, C2, C3}) C1 + C3*x
Hint Functions¶
These functions are intended for internal use by
dsolve()
and others. Unlike User Functions,
above, these are not intended for every-day use by ordinary SymPy users.
Instead, functions such as dsolve()
should be used.
Nonetheless, these functions contain useful information in their docstrings on
the various ODE solving methods. For this reason, they are documented here.
- sympy.solvers.ode.allhints = ('factorable', 'nth_algebraic', 'separable', '1st_exact', '1st_linear', 'Bernoulli', '1st_rational_riccati', 'Riccati_special_minus2', '1st_homogeneous_coeff_best', '1st_homogeneous_coeff_subs_indep_div_dep', '1st_homogeneous_coeff_subs_dep_div_indep', 'almost_linear', 'linear_coefficients', 'separable_reduced', '1st_power_series', 'lie_group', 'nth_linear_constant_coeff_homogeneous', 'nth_linear_euler_eq_homogeneous', 'nth_linear_constant_coeff_undetermined_coefficients', 'nth_linear_euler_eq_nonhomogeneous_undetermined_coefficients', 'nth_linear_constant_coeff_variation_of_parameters', 'nth_linear_euler_eq_nonhomogeneous_variation_of_parameters', 'Liouville', '2nd_linear_airy', '2nd_linear_bessel', '2nd_hypergeometric', '2nd_hypergeometric_Integral', 'nth_order_reducible', '2nd_power_series_ordinary', '2nd_power_series_regular', 'nth_algebraic_Integral', 'separable_Integral', '1st_exact_Integral', '1st_linear_Integral', 'Bernoulli_Integral', '1st_homogeneous_coeff_subs_indep_div_dep_Integral', '1st_homogeneous_coeff_subs_dep_div_indep_Integral', 'almost_linear_Integral', 'linear_coefficients_Integral', 'separable_reduced_Integral', 'nth_linear_constant_coeff_variation_of_parameters_Integral', 'nth_linear_euler_eq_nonhomogeneous_variation_of_parameters_Integral', 'Liouville_Integral', '2nd_nonlinear_autonomous_conserved', '2nd_nonlinear_autonomous_conserved_Integral')¶
Built-in immutable sequence.
If no argument is given, the constructor returns an empty tuple. If iterable is specified the tuple is initialized from iterable’s items.
If the argument is a tuple, the return value is the same object.
- sympy.solvers.ode.ode.odesimp(ode, eq, func, hint)[source]¶
Simplifies solutions of ODEs, including trying to solve for
func
and runningconstantsimp()
.It may use knowledge of the type of solution that the hint returns to apply additional simplifications.
It also attempts to integrate any
Integral
s in the expression, if the hint is not an_Integral
hint.This function should have no effect on expressions returned by
dsolve()
, asdsolve()
already callsodesimp()
, but the individual hint functions do not callodesimp()
(because thedsolve()
wrapper does). Therefore, this function is designed for mainly internal use.Examples
>>> from sympy import sin, symbols, dsolve, pprint, Function >>> from sympy.solvers.ode.ode import odesimp >>> x, u2, C1= symbols('x,u2,C1') >>> f = Function('f')
>>> eq = dsolve(x*f(x).diff(x) - f(x) - x*sin(f(x)/x), f(x), ... hint='1st_homogeneous_coeff_subs_indep_div_dep_Integral', ... simplify=False) >>> pprint(eq, wrap_line=False) x ---- f(x) / | | / 1 \ | -|u1 + -------| | | /1 \| | | sin|--|| | \ \u1// log(f(x)) = log(C1) + | ---------------- d(u1) | 2 | u1 | /
>>> pprint(odesimp(eq, f(x), 1, {C1}, ... hint='1st_homogeneous_coeff_subs_indep_div_dep' ... )) x --------- = C1 /f(x)\ tan|----| \2*x /
- sympy.solvers.ode.ode.constant_renumber(
- expr,
- variables=None,
- newconstants=None,
Renumber arbitrary constants in
expr
to use the symbol names as given innewconstants
. In the process, this reorders expression terms in a standard way.If
newconstants
is not provided then the new constant names will beC1
,C2
etc. Otherwisenewconstants
should be an iterable giving the new symbols to use for the constants in order.The
variables
argument is a list of non-constant symbols. All other free symbols found inexpr
are assumed to be constants and will be renumbered. Ifvariables
is not given then any numbered symbol beginning withC
(e.g.C1
) is assumed to be a constant.Symbols are renumbered based on
.sort_key()
, so they should be numbered roughly in the order that they appear in the final, printed expression. Note that this ordering is based in part on hashes, so it can produce different results on different machines.The structure of this function is very similar to that of
constantsimp()
.Examples
>>> from sympy import symbols >>> from sympy.solvers.ode.ode import constant_renumber >>> x, C1, C2, C3 = symbols('x,C1:4') >>> expr = C3 + C2*x + C1*x**2 >>> expr C1*x**2 + C2*x + C3 >>> constant_renumber(expr) C1 + C2*x + C3*x**2
The
variables
argument specifies which are constants so that the other symbols will not be renumbered:>>> constant_renumber(expr, [C1, x]) C1*x**2 + C2 + C3*x
The
newconstants
argument is used to specify what symbols to use when replacing the constants:>>> constant_renumber(expr, [x], newconstants=symbols('E1:4')) E1 + E2*x + E3*x**2
- sympy.solvers.ode.ode.ode_sol_simplicity(sol, func, trysolving=True)[source]¶
Returns an extended integer representing how simple a solution to an ODE is.
The following things are considered, in order from most simple to least:
sol
is solved forfunc
.sol
is not solved forfunc
, but can be if passed to solve (e.g., a solution returned bydsolve(ode, func, simplify=False
).If
sol
is not solved forfunc
, then base the result on the length ofsol
, as computed bylen(str(sol))
.If
sol
has any unevaluatedIntegral
s, this will automatically be considered less simple than any of the above.
This function returns an integer such that if solution A is simpler than solution B by above metric, then
ode_sol_simplicity(sola, func) < ode_sol_simplicity(solb, func)
.Currently, the following are the numbers returned, but if the heuristic is ever improved, this may change. Only the ordering is guaranteed.
Simplicity
Return
sol
solved forfunc
-2
sol
not solved forfunc
but can be-1
sol
is not solved nor solvable forfunc
len(str(sol))
sol
contains anIntegral
oo
oo
here means the SymPy infinity, which should compare greater than any integer.If you already know
solve()
cannot solvesol
, you can usetrysolving=False
to skip that step, which is the only potentially slow step. For example,dsolve()
with thesimplify=False
flag should do this.If
sol
is a list of solutions, if the worst solution in the list returnsoo
it returns that, otherwise it returnslen(str(sol))
, that is, the length of the string representation of the whole list.Examples
This function is designed to be passed to
min
as the key argument, such asmin(listofsolutions, key=lambda i: ode_sol_simplicity(i, f(x)))
.>>> from sympy import symbols, Function, Eq, tan, Integral >>> from sympy.solvers.ode.ode import ode_sol_simplicity >>> x, C1, C2 = symbols('x, C1, C2') >>> f = Function('f')
>>> ode_sol_simplicity(Eq(f(x), C1*x**2), f(x)) -2 >>> ode_sol_simplicity(Eq(x**2 + f(x), C1), f(x)) -1 >>> ode_sol_simplicity(Eq(f(x), C1*Integral(2*x, x)), f(x)) oo >>> eq1 = Eq(f(x)/tan(f(x)/(2*x)), C1) >>> eq2 = Eq(f(x)/tan(f(x)/(2*x) + f(x)), C2) >>> [ode_sol_simplicity(eq, f(x)) for eq in [eq1, eq2]] [28, 35] >>> min([eq1, eq2], key=lambda i: ode_sol_simplicity(i, f(x))) Eq(f(x)/tan(f(x)/(2*x)), C1)
- class sympy.solvers.ode.single.Factorable(ode_problem)[source]¶
Solves equations having a solvable factor.
This function is used to solve the equation having factors. Factors may be of type algebraic or ode. It will try to solve each factor independently. Factors will be solved by calling dsolve. We will return the list of solutions.
Examples
>>> from sympy import Function, dsolve, pprint >>> from sympy.abc import x >>> f = Function('f') >>> eq = (f(x)**2-4)*(f(x).diff(x)+f(x)) >>> pprint(dsolve(eq, f(x))) -x [f(x) = 2, f(x) = -2, f(x) = C1*e ]
- class sympy.solvers.ode.single.FirstExact(ode_problem)[source]¶
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) \,\partial{}x + Q(x, y) \,\partial{}y = 0\]is exact if there is some function \(F(x, y)\) such that \(P(x, y) = \partial{}F/\partial{}x\) and \(Q(x, y) = \partial{}F/\partial{}y\). It can be shown that a necessary and sufficient condition for a first order ODE to be exact is that \(\partial{}P/\partial{}y = \partial{}Q/\partial{}x\). 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.Examples
>>> 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') Eq(x*cos(f(x)) + f(x)**3/3, C1)
References
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 73
# indirect doctest
- class sympy.solvers.ode.single.HomogeneousCoeffBest(ode_problem)[source]¶
Returns the best solution to an ODE from the two hints
1st_homogeneous_coeff_subs_dep_div_indep
and1st_homogeneous_coeff_subs_indep_div_dep
.This is as determined by
ode_sol_simplicity()
.See the
HomogeneousCoeffSubsIndepDivDep
andHomogeneousCoeffSubsDepDivIndep
docstrings for more information on these hints. Note that there is noode_1st_homogeneous_coeff_best_Integral
hint.Examples
>>> 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', simplify=False)) / 2 \ |3*x | log|----- + 1| | 2 | \f (x) / log(f(x)) = log(C1) - -------------- 3
References
https://en.wikipedia.org/wiki/Homogeneous_differential_equation
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 59
# indirect doctest
- class sympy.solvers.ode.single.HomogeneousCoeffSubsDepDivIndep(ode_problem)[source]¶
Solves a 1st order differential equation with homogeneous coefficients using the substitution \(u_1 = \frac{\text{<dependent variable>}}{\text{<independent variable>}}\).
This is a differential equation
\[P(x, y) + Q(x, y) dy/dx = 0\]such that \(P\) and \(Q\) are homogeneous and of the same order. A function \(F(x, y)\) is homogeneous of order \(n\) if \(F(x t, y t) = 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 = u_1 x\) (i.e. \(u_1 = y/x\)) will turn the differential equation into an equation separable in the variables \(x\) and \(u\). If \(h(u_1)\) is the function that results from making the substitution \(u_1 = f(x)/x\) on \(P(x, f(x))\) and \(g(u_2)\) 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)) f'(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(x) = C1 + | ---------------- d(u1) | u1*h(u1) + g(u1) | /
Where \(u_1 h(u_1) + g(u_1) \ne 0\) and \(x \ne 0\).
See also the docstrings of
HomogeneousCoeffBest
andHomogeneousCoeffSubsIndepDivDep
.Examples
>>> 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', simplify=False)) / 3 \ |3*f(x) f (x)| log|------ + -----| | x 3 | \ x / log(x) = log(C1) - ------------------- 3
References
https://en.wikipedia.org/wiki/Homogeneous_differential_equation
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 59
# indirect doctest
- class sympy.solvers.ode.single.HomogeneousCoeffSubsIndepDivDep(ode_problem)[source]¶
Solves a 1st order differential equation with homogeneous coefficients using the substitution \(u_2 = \frac{\text{<independent variable>}}{\text{<dependent variable>}}\).
This is a differential equation
\[P(x, y) + Q(x, y) dy/dx = 0\]such that \(P\) and \(Q\) are homogeneous and of the same order. A function \(F(x, y)\) is homogeneous of order \(n\) if \(F(x t, y t) = 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 = u_2 y\) (i.e. \(u_2 = x/y\)) will turn the differential equation into an equation separable in the variables \(y\) and \(u_2\). If \(h(u_2)\) is the function that results from making the substitution \(u_2 = x/f(x)\) on \(P(x, f(x))\) and \(g(u_2)\) 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)) f'(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(u1) | ---------------- d(u1) | u1*g(u1) + h(u1) | / f(x) = C1*e
Where \(u_1 g(u_1) + h(u_1) \ne 0\) and \(f(x) \ne 0\).
See also the docstrings of
HomogeneousCoeffBest
andHomogeneousCoeffSubsDepDivIndep
.Examples
>>> from sympy import Function, pprint, 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_indep_div_dep', ... simplify=False)) / 2 \ |3*x | log|----- + 1| | 2 | \f (x) / log(f(x)) = log(C1) - -------------- 3
References
https://en.wikipedia.org/wiki/Homogeneous_differential_equation
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 59
# indirect doctest
- class sympy.solvers.ode.single.FirstLinear(ode_problem)[source]¶
Solves 1st order linear differential equations.
These are differential equations of the form
\[dy/dx + P(x) y = Q(x)\text{.}\]These kinds of differential equations can be solved in a general way. The integrating factor \(e^{\int P(x) \,dx}\) 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 | | | \ / /
Examples
>>> 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
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 92
# indirect doctest
- class sympy.solvers.ode.single.RationalRiccati(ode_problem)[source]¶
Gives general solutions to the first order Riccati differential equations that have atleast one rational particular solution.
\[y' = b_0(x) + b_1(x) y + b_2(x) y^2\]where \(b_0\), \(b_1\) and \(b_2\) are rational functions of \(x\) with \(b_2 \ne 0\) (\(b_2 = 0\) would make it a Bernoulli equation).
Examples
>>> from sympy import Symbol, Function, dsolve, checkodesol >>> f = Function('f') >>> x = Symbol('x')
>>> eq = -x**4*f(x)**2 + x**3*f(x).diff(x) + x**2*f(x) + 20 >>> sol = dsolve(eq, hint="1st_rational_riccati") >>> sol Eq(f(x), (4*C1 - 5*x**9 - 4)/(x**2*(C1 + x**9 - 1))) >>> checkodesol(eq, sol) (True, 0)
References
Riccati ODE: https://en.wikipedia.org/wiki/Riccati_equation
N. Thieu Vo - Rational and Algebraic Solutions of First-Order Algebraic ODEs: Algorithm 11, pp. 78 - https://www3.risc.jku.at/publications/download/risc_5387/PhDThesisThieu.pdf
- class sympy.solvers.ode.single.SecondLinearAiry(ode_problem)[source]¶
Gives solution of the Airy differential equation
\[\frac{d^2y}{dx^2} + (a + b x) y(x) = 0\]in terms of Airy special functions airyai and airybi.
Examples
>>> from sympy import dsolve, Function >>> from sympy.abc import x >>> f = Function("f") >>> eq = f(x).diff(x, 2) - x*f(x) >>> dsolve(eq) Eq(f(x), C1*airyai(x) + C2*airybi(x))
- class sympy.solvers.ode.single.SecondLinearBessel(ode_problem)[source]¶
Gives solution of the Bessel differential equation
\[x^2 \frac{d^2y}{dx^2} + x \frac{dy}{dx} y(x) + (x^2-n^2) y(x)\]if \(n\) is integer then the solution is of the form
Eq(f(x), C0 besselj(n,x) + C1 bessely(n,x))
as both the solutions are linearly independent else if \(n\) is a fraction then the solution is of the formEq(f(x), C0 besselj(n,x) + C1 besselj(-n,x))
which can also transform intoEq(f(x), C0 besselj(n,x) + C1 bessely(n,x))
.Examples
>>> from sympy.abc import x >>> from sympy import Symbol >>> v = Symbol('v', positive=True) >>> from sympy import dsolve, Function >>> f = Function('f') >>> y = f(x) >>> genform = x**2*y.diff(x, 2) + x*y.diff(x) + (x**2 - v**2)*y >>> dsolve(genform) Eq(f(x), C1*besselj(v, x) + C2*bessely(v, x))
References
- class sympy.solvers.ode.single.Bernoulli(ode_problem)[source]¶
Solves Bernoulli differential equations.
These are equations of the form
\[dy/dx + P(x) y = Q(x) y^n\text{, }n \ne 1`\text{.}\]The substitution \(w = 1/y^{1-n}\) will transform an equation of this form into one that is linear (see the docstring of
FirstLinear
). 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'), num_columns=110) -1 ----- n - 1 // / / \ \ || | | | | || | / | / | / | || | | | | | | | || | -(n - 1)* | P(x) dx | -(n - 1)* | P(x) dx | (n - 1)* | P(x) dx| || | | | | | | | || | / | / | / | f(x) = ||C1 - n* | Q(x)*e dx + | Q(x)*e dx|*e | || | | | | \\ / / / /
Note that the equation is separable when \(n = 1\) (see the docstring of
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 | | / /
Examples
>>> 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) = ----------------- C1*x + log(x) + 1
References
https://en.wikipedia.org/wiki/Bernoulli_differential_equation
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 95
# indirect doctest
- class sympy.solvers.ode.single.Liouville(ode_problem)[source]¶
Solves 2nd order Liouville differential equations.
The general form of a Liouville ODE is
\[\frac{d^2 y}{dx^2} + g(y) \left(\! \frac{dy}{dx}\!\right)^2 + h(x) \frac{dy}{dx}\text{.}\]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 | | / /
Examples
>>> 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) = -\/ C1 + C2*log(x) , f(x) = \/ C1 + C2*log(x) ]
References
Goldstein and Braun, “Advanced Methods for the Solution of Differential Equations”, pp. 98
https://www.maplesoft.com/support/help/Maple/view.aspx?path=odeadvisor/Liouville
# indirect doctest
- class sympy.solvers.ode.single.RiccatiSpecial(ode_problem)[source]¶
The general Riccati equation has the form
\[dy/dx = f(x) y^2 + g(x) y + h(x)\text{.}\]While it does not have a general solution [1], the “special” form, \(dy/dx = a y^2 - b x^c\), does have solutions in many cases [2]. This routine returns a solution for \(a(dy/dx) = b y^2 + c y/x + d/x^2\) that is obtained by using a suitable change of variables to reduce it to the special form and is valid when neither \(a\) nor \(b\) are zero and either \(c\) or \(d\) is zero.
>>> from sympy.abc import x, a, b, c, d >>> from sympy import dsolve, checkodesol, pprint, Function >>> f = Function('f') >>> y = f(x) >>> genform = a*y.diff(x) - (b*y**2 + c*y/x + d/x**2) >>> sol = dsolve(genform, y, hint="Riccati_special_minus2") >>> pprint(sol, wrap_line=False) / / __________________ \\ | __________________ | / 2 || | / 2 | \/ 4*b*d - (a + c) *log(x)|| -|a + c - \/ 4*b*d - (a + c) *tan|C1 + ----------------------------|| \ \ 2*a // f(x) = ------------------------------------------------------------------------ 2*b*x
>>> checkodesol(genform, sol, order=1)[0] True
References
- class sympy.solvers.ode.single.NthLinearConstantCoeffHomogeneous(ode_problem)[source]¶
Solves an \(n\)th order linear homogeneous differential equation with constant coefficients.
This is an equation of the form
\[a_n f^{(n)}(x) + a_{n-1} f^{(n-1)}(x) + \cdots + a_1 f'(x) + a_0 f(x) = 0\text{.}\]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} + \cdots + a_1 m + a_0 = 0\). The solution will then be the sum of \(C_n x^i e^{r x}\) terms, for each where \(C_n\) is an arbitrary constant, \(r\) is a root of the characteristic equation and \(i\) 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 \(C_1 e^{3 x} + C_2 x e^{3 x}\)). The exponential is usually expanded for complex roots using Euler’s equation \(e^{I x} = \cos(x) + I \sin(x)\). Complex roots always come in conjugate pairs in polynomials with real coefficients, so the two roots will be represented (after simplifying the constants) as \(e^{a x} \left(C_1 \cos(b x) + C_2 \sin(b x)\right)\).
If SymPy cannot find exact roots to the characteristic equation, a
ComplexRootOf
instance will be return instead.>>> from sympy import Function, dsolve >>> 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') ... Eq(f(x), C5*exp(x*CRootOf(_x**5 + 10*_x - 2, 0)) + (C1*sin(x*im(CRootOf(_x**5 + 10*_x - 2, 1))) + C2*cos(x*im(CRootOf(_x**5 + 10*_x - 2, 1))))*exp(x*re(CRootOf(_x**5 + 10*_x - 2, 1))) + (C3*sin(x*im(CRootOf(_x**5 + 10*_x - 2, 3))) + C4*cos(x*im(CRootOf(_x**5 + 10*_x - 2, 3))))*exp(x*re(CRootOf(_x**5 + 10*_x - 2, 3))))
Note that because this method does not involve integration, there is no
nth_linear_constant_coeff_homogeneous_Integral
hint.Examples
>>> 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*sin(x) + C4*cos(x))*e
References
https://en.wikipedia.org/wiki/Linear_differential_equation section: Nonhomogeneous_equation_with_constant_coefficients
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 211
# indirect doctest
- class sympy.solvers.ode.single.NthLinearConstantCoeffUndeterminedCoefficients(
- ode_problem,
Solves an \(n\)th 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^{(n)}(x) + a_{n-1} f^{(n-1)}(x) + \cdots + a_1 f'(x) + a_0 f(x) = P(x)\text{,}\]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 e^{b x} \sin(c x + d)\) or \(a x^i e^{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 e^{2 x}\), \(x \sin(x)\), and \(e^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^2(x)\) 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.
Examples
>>> 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')) / / 3\\ | | x || -x 4*sin(2*x) 3*cos(2*x) f(x) = |C1 + x*|C2 + --||*e - ---------- + ---------- \ \ 3 // 25 25
References
https://en.wikipedia.org/wiki/Method_of_undetermined_coefficients
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 221
# indirect doctest
- class sympy.solvers.ode.single.NthLinearConstantCoeffVariationOfParameters(
- ode_problem,
Solves an \(n\)th order linear differential equation with constant coefficients using the method of variation of parameters.
This method works on any differential equations of the form
\[f^{(n)}(x) + a_{n-1} f^{(n-1)}(x) + \cdots + a_1 f'(x) + a_0 f(x) = P(x)\text{.}\]This method works by assuming that the particular solution takes the form
\[\sum_{x=1}^{n} c_i(x) y_i(x)\text{,}\]where \(y_i\) is the \(i\)th solution to the homogeneous equation. The solution is then solved using Wronskian’s and Cramer’s Rule. The particular solution is given by
\[\sum_{x=1}^n \left( \int \frac{W_i(x)}{W(x)} \,dx \right) y_i(x) \text{,}\]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 \(i\)th column replaced with \([0, 0, \cdots, 0, P(x)]\).
This method is general enough to solve any \(n\)th 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 usingnth_linear_constant_coeff_undetermined_coefficients
when it applies, because it does not 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.Examples
>>> 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')) / / / x*log(x) 11*x\\\ x f(x) = |C1 + x*|C2 + x*|C3 + -------- - ----|||*e \ \ \ 6 36 ///
References
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 233
# indirect doctest
- class sympy.solvers.ode.single.NthLinearEulerEqHomogeneous(ode_problem)[source]¶
Solves an \(n\)th order linear homogeneous variable-coefficient Cauchy-Euler equidimensional ordinary differential equation.
This is an equation with form \(0 = a_0 f(x) + a_1 x f'(x) + a_2 x^2 f''(x) \cdots\).
These equations can be solved in a general manner, by substituting solutions of the form \(f(x) = x^r\), and deriving a characteristic equation for \(r\). When there are repeated roots, we include extra terms of the form \(C_{r k} \ln^k(x) x^r\), where \(C_{r k}\) is an arbitrary integration constant, \(r\) is a root of the characteristic equation, and \(k\) ranges over the multiplicity of \(r\). In the cases where the roots are complex, solutions of the form \(C_1 x^a \sin(b \log(x)) + C_2 x^a \cos(b \log(x))\) are returned, based on expansions with Euler’s formula. The general solution is the sum of the terms found. If SymPy cannot find exact roots to the characteristic equation, a
ComplexRootOf
instance will be returned instead.>>> from sympy import Function, dsolve >>> from sympy.abc import x >>> f = Function('f') >>> dsolve(4*x**2*f(x).diff(x, 2) + f(x), f(x), ... hint='nth_linear_euler_eq_homogeneous') ... Eq(f(x), sqrt(x)*(C1 + C2*log(x)))
Note that because this method does not involve integration, there is no
nth_linear_euler_eq_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, corresponding to the fundamental solution set, 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 doassert 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>}
.
Examples
>>> from sympy import Function, dsolve, pprint >>> from sympy.abc import x >>> f = Function('f') >>> eq = f(x).diff(x, 2)*x**2 - 4*f(x).diff(x)*x + 6*f(x) >>> pprint(dsolve(eq, f(x), ... hint='nth_linear_euler_eq_homogeneous')) 2 f(x) = x *(C1 + C2*x)
References
C. Bender & S. Orszag, “Advanced Mathematical Methods for Scientists and Engineers”, Springer 1999, pp. 12
# indirect doctest
- class sympy.solvers.ode.single.NthLinearEulerEqNonhomogeneousVariationOfParameters(
- ode_problem,
Solves an \(n\)th order linear non homogeneous Cauchy-Euler equidimensional ordinary differential equation using variation of parameters.
This is an equation with form \(g(x) = a_0 f(x) + a_1 x f'(x) + a_2 x^2 f''(x) \cdots\).
This method works by assuming that the particular solution takes the form
\[\sum_{x=1}^{n} c_i(x) y_i(x) {a_n} {x^n} \text{, }\]where \(y_i\) is the \(i\)th solution to the homogeneous equation. The solution is then solved using Wronskian’s and Cramer’s Rule. The particular solution is given by multiplying eq given below with \(a_n x^{n}\)
\[\sum_{x=1}^n \left( \int \frac{W_i(x)}{W(x)} \, dx \right) y_i(x) \text{, }\]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 \(i\)th column replaced with \([0, 0, \cdots, 0, \frac{x^{- n}}{a_n} g{\left(x \right)}]\).
This method is general enough to solve any \(n\)th order inhomogeneous linear differential equation, 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 usingnth_linear_constant_coeff_undetermined_coefficients
when it applies, because it does not 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.Examples
>>> from sympy import Function, dsolve, Derivative >>> from sympy.abc import x >>> f = Function('f') >>> eq = x**2*Derivative(f(x), x, x) - 2*x*Derivative(f(x), x) + 2*f(x) - x**4 >>> dsolve(eq, f(x), ... hint='nth_linear_euler_eq_nonhomogeneous_variation_of_parameters').expand() Eq(f(x), C1*x + C2*x**2 + x**4/6)
- class sympy.solvers.ode.single.NthLinearEulerEqNonhomogeneousUndeterminedCoefficients(
- ode_problem,
Solves an \(n\)th order linear non homogeneous Cauchy-Euler equidimensional ordinary differential equation using undetermined coefficients.
This is an equation with form \(g(x) = a_0 f(x) + a_1 x f'(x) + a_2 x^2 f''(x) \cdots\).
These equations can be solved in a general manner, by substituting solutions of the form \(x = exp(t)\), and deriving a characteristic equation of form \(g(exp(t)) = b_0 f(t) + b_1 f'(t) + b_2 f''(t) \cdots\) which can be then solved by nth_linear_constant_coeff_undetermined_coefficients if g(exp(t)) has finite number of linearly independent derivatives.
Functions that fit this requirement are finite sums functions of the form \(a x^i e^{b x} \sin(c x + d)\) or \(a x^i e^{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 e^{2 x}\), \(x \sin(x)\), and \(e^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^2(x)\) into \((1 + \cos(2 x))/2\) to properly apply the method of undetermined coefficients on it.
After replacement of x by exp(t), 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.
Examples
>>> from sympy import dsolve, Function, Derivative, log >>> from sympy.abc import x >>> f = Function('f') >>> eq = x**2*Derivative(f(x), x, x) - 2*x*Derivative(f(x), x) + 2*f(x) - log(x) >>> dsolve(eq, f(x), ... hint='nth_linear_euler_eq_nonhomogeneous_undetermined_coefficients').expand() Eq(f(x), C1*x + C2*x**2 + log(x)/2 + 3/4)
- class sympy.solvers.ode.single.NthAlgebraic(ode_problem)[source]¶
Solves an \(n\)th order ordinary differential equation using algebra and integrals.
There is no general form for the kind of equation that this can solve. The the equation is solved algebraically treating differentiation as an invertible algebraic function.
Examples
>>> from sympy import Function, dsolve, Eq >>> from sympy.abc import x >>> f = Function('f') >>> eq = Eq(f(x) * (f(x).diff(x)**2 - 1), 0) >>> dsolve(eq, f(x), hint='nth_algebraic') [Eq(f(x), 0), Eq(f(x), C1 - x), Eq(f(x), C1 + x)]
Note that this solver can return algebraic solutions that do not have any integration constants (f(x) = 0 in the above example).
- class sympy.solvers.ode.single.NthOrderReducible(ode_problem)[source]¶
Solves ODEs that only involve derivatives of the dependent variable using a substitution of the form \(f^n(x) = g(x)\).
For example any second order ODE of the form \(f''(x) = h(f'(x), x)\) can be transformed into a pair of 1st order ODEs \(g'(x) = h(g(x), x)\) and \(f'(x) = g(x)\). Usually the 1st order ODE for \(g\) is easier to solve. If that gives an explicit solution for \(g\) then \(f\) is found simply by integration.
Examples
>>> from sympy import Function, dsolve, Eq >>> from sympy.abc import x >>> f = Function('f') >>> eq = Eq(x*f(x).diff(x)**2 + f(x).diff(x, 2), 0) >>> dsolve(eq, f(x), hint='nth_order_reducible') ... Eq(f(x), C1 - sqrt(-1/C2)*log(-C2*sqrt(-1/C2) + x) + sqrt(-1/C2)*log(C2*sqrt(-1/C2) + x))
- class sympy.solvers.ode.single.Separable(ode_problem)[source]¶
Solves separable 1st order differential equations.
This is any differential equation that can be written as \(P(y) \tfrac{dy}{dx} = Q(x)\). The solution can then just be found by rearranging terms and integrating: \(\int P(y) \,dy = \int Q(x) \,dx\). This hint uses
sympy.simplify.simplify.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)\cdot{}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) | | / /
Examples
>>> 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', simplify=False)) / 2 \ 2 log\3*f (x) - 1/ x ---------------- = C1 + -- 6 2
References
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 52
# indirect doctest
- class sympy.solvers.ode.single.AlmostLinear(ode_problem)[source]¶
Solves an almost-linear differential equation.
The general form of an almost linear differential equation is
\[a(x) g'(f(x)) f'(x) + b(x) g(f(x)) + c(x)\]Here \(f(x)\) is the function to be solved for (the dependent variable). The substitution \(g(f(x)) = u(x)\) leads to a linear differential equation for \(u(x)\) of the form \(a(x) u' + b(x) u + c(x) = 0\). This can be solved for \(u(x)\) by the \(first_linear\) hint and then \(f(x)\) is found by solving \(g(f(x)) = u(x)\).
Examples
>>> from sympy import dsolve, Function, pprint, sin, cos >>> from sympy.abc import x >>> f = Function('f') >>> d = f(x).diff(x) >>> eq = x*d + x*f(x) + 1 >>> dsolve(eq, f(x), hint='almost_linear') Eq(f(x), (C1 - Ei(x))*exp(-x)) >>> pprint(dsolve(eq, f(x), hint='almost_linear')) -x f(x) = (C1 - Ei(x))*e >>> example = cos(f(x))*f(x).diff(x) + sin(f(x)) + 1 >>> pprint(example) d sin(f(x)) + cos(f(x))*--(f(x)) + 1 dx >>> pprint(dsolve(example, f(x), hint='almost_linear')) / -x \ / -x \ [f(x) = pi - asin\C1*e - 1/, f(x) = asin\C1*e - 1/]
See also
References
Joel Moses, “Symbolic Integration - The Stormy Decade”, Communications of the ACM, Volume 14, Number 8, August 1971, pp. 558
- class sympy.solvers.ode.single.LinearCoefficients(ode_problem)[source]¶
Solves a differential equation with linear coefficients.
The general form of a differential equation with linear coefficients is
\[y' + F\left(\!\frac{a_1 x + b_1 y + c_1}{a_2 x + b_2 y + c_2}\!\right) = 0\text{,}\]where \(a_1\), \(b_1\), \(c_1\), \(a_2\), \(b_2\), \(c_2\) are constants and \(a_1 b_2 - a_2 b_1 \ne 0\).
This can be solved by substituting:
\[ \begin{align}\begin{aligned}x = x' + \frac{b_2 c_1 - b_1 c_2}{a_2 b_1 - a_1 b_2}\\y = y' + \frac{a_1 c_2 - a_2 c_1}{a_2 b_1 - a_1 b_2}\text{.}\end{aligned}\end{align} \]This substitution reduces the equation to a homogeneous differential equation.
Examples
>>> from sympy import dsolve, Function, pprint >>> from sympy.abc import x >>> f = Function('f') >>> df = f(x).diff(x) >>> eq = (x + f(x) + 1)*df + (f(x) - 6*x + 1) >>> dsolve(eq, hint='linear_coefficients') [Eq(f(x), -x - sqrt(C1 + 7*x**2) - 1), Eq(f(x), -x + sqrt(C1 + 7*x**2) - 1)] >>> pprint(dsolve(eq, hint='linear_coefficients')) ___________ ___________ / 2 / 2 [f(x) = -x - \/ C1 + 7*x - 1, f(x) = -x + \/ C1 + 7*x - 1]
See also
sympy.solvers.ode.single.HomogeneousCoeffBest
,sympy.solvers.ode.single.HomogeneousCoeffSubsIndepDivDep
,sympy.solvers.ode.single.HomogeneousCoeffSubsDepDivIndep
References
Joel Moses, “Symbolic Integration - The Stormy Decade”, Communications of the ACM, Volume 14, Number 8, August 1971, pp. 558
- class sympy.solvers.ode.single.SeparableReduced(ode_problem)[source]¶
Solves a differential equation that can be reduced to the separable form.
The general form of this equation is
\[y' + (y/x) H(x^n y) = 0\text{}.\]This can be solved by substituting \(u(y) = x^n y\). The equation then reduces to the separable form \(\frac{u'}{u (\mathrm{power} - H(u))} - \frac{1}{x} = 0\).
The general solution is:
>>> from sympy import Function, dsolve, pprint >>> from sympy.abc import x, n >>> f, g = map(Function, ['f', 'g']) >>> genform = f(x).diff(x) + (f(x)/x)*g(x**n*f(x)) >>> pprint(genform) / n \ d f(x)*g\x *f(x)/ --(f(x)) + --------------- dx x >>> pprint(dsolve(genform, hint='separable_reduced')) n x *f(x) / | | 1 | ------------ dy = C1 + log(x) | y*(n - g(y)) | /
Examples
>>> from sympy import dsolve, Function, pprint >>> from sympy.abc import x >>> f = Function('f') >>> d = f(x).diff(x) >>> eq = (x - x**2*f(x))*d - f(x) >>> dsolve(eq, hint='separable_reduced') [Eq(f(x), (1 - sqrt(C1*x**2 + 1))/x), Eq(f(x), (sqrt(C1*x**2 + 1) + 1)/x)] >>> pprint(dsolve(eq, hint='separable_reduced')) ___________ ___________ / 2 / 2 1 - \/ C1*x + 1 \/ C1*x + 1 + 1 [f(x) = ------------------, f(x) = ------------------] x x
See also
References
Joel Moses, “Symbolic Integration - The Stormy Decade”, Communications of the ACM, Volume 14, Number 8, August 1971, pp. 558
- class sympy.solvers.ode.single.LieGroup(ode_problem)[source]¶
This hint implements the Lie group method of solving first order differential equations. The aim is to convert the given differential equation from the given coordinate system into another coordinate system where it becomes invariant under the one-parameter Lie group of translations. The converted ODE can be easily solved by quadrature. It makes use of the
sympy.solvers.ode.infinitesimals()
function which returns the infinitesimals of the transformation.The coordinates \(r\) and \(s\) can be found by solving the following Partial Differential Equations.
\[\xi\frac{\partial r}{\partial x} + \eta\frac{\partial r}{\partial y} = 0\]\[\xi\frac{\partial s}{\partial x} + \eta\frac{\partial s}{\partial y} = 1\]The differential equation becomes separable in the new coordinate system
\[\frac{ds}{dr} = \frac{\frac{\partial s}{\partial x} + h(x, y)\frac{\partial s}{\partial y}}{ \frac{\partial r}{\partial x} + h(x, y)\frac{\partial r}{\partial y}}\]After finding the solution by integration, it is then converted back to the original coordinate system by substituting \(r\) and \(s\) in terms of \(x\) and \(y\) again.
Examples
>>> from sympy import Function, dsolve, exp, pprint >>> from sympy.abc import x >>> f = Function('f') >>> pprint(dsolve(f(x).diff(x) + 2*x*f(x) - x*exp(-x**2), f(x), ... hint='lie_group')) / 2\ 2 | x | -x f(x) = |C1 + --|*e \ 2 /
References
Solving differential equations by Symmetry Groups, John Starrett, pp. 1 - pp. 14
- class sympy.solvers.ode.single.SecondHypergeometric(ode_problem)[source]¶
Solves 2nd order linear differential equations.
It computes special function solutions which can be expressed using the 2F1, 1F1 or 0F1 hypergeometric functions.
\[y'' + A(x) y' + B(x) y = 0\text{,}\]where \(A\) and \(B\) are rational functions.
These kinds of differential equations have solution of non-Liouvillian form.
Given linear ODE can be obtained from 2F1 given by
\[(x^2 - x) y'' + ((a + b + 1) x - c) y' + b a y = 0\text{,}\]where {a, b, c} are arbitrary constants.
Notes
The algorithm should find any solution of the form
\[y = P(x) _pF_q(..; ..;\frac{\alpha x^k + \beta}{\gamma x^k + \delta})\text{,}\]where pFq is any of 2F1, 1F1 or 0F1 and \(P\) is an “arbitrary function”. Currently only the 2F1 case is implemented in SymPy but the other cases are described in the paper and could be implemented in future (contributions welcome!).
Examples
>>> from sympy import Function, dsolve, pprint >>> from sympy.abc import x >>> f = Function('f') >>> eq = (x*x - x)*f(x).diff(x,2) + (5*x - 1)*f(x).diff(x) + 4*f(x) >>> pprint(dsolve(eq, f(x), '2nd_hypergeometric')) _ / / 4 \\ |_ /-1, -1 | \ |C1 + C2*|log(x) + -----||* | | | x| \ \ x + 1// 2 1 \ 1 | / f(x) = -------------------------------------------- 3 (x - 1)
References
“Non-Liouvillian solutions for second order linear ODEs” by L. Chan, E.S. Cheb-Terrab
- sympy.solvers.ode.ode.ode_1st_power_series(eq, func, order, match)[source]¶
The power series solution is a method which gives the Taylor series expansion to the solution of a differential equation.
For a first order differential equation \(\frac{dy}{dx} = h(x, y)\), a power series solution exists at a point \(x = x_{0}\) if \(h(x, y)\) is analytic at \(x_{0}\). The solution is given by
\[y(x) = y(x_{0}) + \sum_{n = 1}^{\infty} \frac{F_{n}(x_{0},b)(x - x_{0})^n}{n!},\]where \(y(x_{0}) = b\) is the value of y at the initial value of \(x_{0}\). To compute the values of the \(F_{n}(x_{0},b)\) the following algorithm is followed, until the required number of terms are generated.
\(F_1 = h(x_{0}, b)\)
\(F_{n+1} = \frac{\partial F_{n}}{\partial x} + \frac{\partial F_{n}}{\partial y}F_{1}\)
Examples
>>> from sympy import Function, pprint, exp, dsolve >>> from sympy.abc import x >>> f = Function('f') >>> eq = exp(x)*(f(x).diff(x)) - f(x) >>> pprint(dsolve(eq, hint='1st_power_series')) 3 4 5 C1*x C1*x C1*x / 6\ f(x) = C1 + C1*x - ----- + ----- + ----- + O\x / 6 24 60
References
Travis W. Walker, Analytic power series technique for solving first-order differential equations, p.p 17, 18
- sympy.solvers.ode.ode.ode_2nd_power_series_ordinary(
- eq,
- func,
- order,
- match,
Gives a power series solution to a second order homogeneous differential equation with polynomial coefficients at an ordinary point. A homogeneous differential equation is of the form
\[P(x)\frac{d^2y}{dx^2} + Q(x)\frac{dy}{dx} + R(x) y(x) = 0\]For simplicity it is assumed that \(P(x)\), \(Q(x)\) and \(R(x)\) are polynomials, it is sufficient that \(\frac{Q(x)}{P(x)}\) and \(\frac{R(x)}{P(x)}\) exists at \(x_{0}\). A recurrence relation is obtained by substituting \(y\) as \(\sum_{n=0}^\infty a_{n}x^{n}\), in the differential equation, and equating the nth term. Using this relation various terms can be generated.
Examples
>>> from sympy import dsolve, Function, pprint >>> from sympy.abc import x >>> f = Function("f") >>> eq = f(x).diff(x, 2) + f(x) >>> pprint(dsolve(eq, hint='2nd_power_series_ordinary')) / 4 2 \ / 2\ |x x | | x | / 6\ f(x) = C2*|-- - -- + 1| + C1*x*|1 - --| + O\x / \24 2 / \ 6 /
References
https://tutorial.math.lamar.edu/Classes/DE/SeriesSolutions.aspx
George E. Simmons, “Differential Equations with Applications and Historical Notes”, p.p 176 - 184
- sympy.solvers.ode.ode.ode_2nd_power_series_regular(
- eq,
- func,
- order,
- match,
Gives a power series solution to a second order homogeneous differential equation with polynomial coefficients at a regular point. A second order homogeneous differential equation is of the form
\[P(x)\frac{d^2y}{dx^2} + Q(x)\frac{dy}{dx} + R(x) y(x) = 0\]A point is said to regular singular at \(x0\) if \(x - x0\frac{Q(x)}{P(x)}\) and \((x - x0)^{2}\frac{R(x)}{P(x)}\) are analytic at \(x0\). For simplicity \(P(x)\), \(Q(x)\) and \(R(x)\) are assumed to be polynomials. The algorithm for finding the power series solutions is:
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.
Solve the indicial equation \(f(m) = m(m - 1) + m*p0 + q0\), to obtain the roots \(m1\) and \(m2\) of the indicial equation.
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.
Examples
>>> from sympy import dsolve, Function, pprint >>> from sympy.abc import x >>> f = Function("f") >>> eq = x*(f(x).diff(x, 2)) + 2*(f(x).diff(x)) + x*f(x) >>> pprint(dsolve(eq, hint='2nd_power_series_regular')) / 6 4 2 \ | x x x | / 4 2 \ C1*|- --- + -- - -- + 1| |x x | \ 720 24 2 / / 6\ f(x) = C2*|--- - -- + 1| + ------------------------ + O\x / \120 6 / x
References
George E. Simmons, “Differential Equations with Applications and Historical Notes”, p.p 176 - 184
Lie heuristics¶
These functions are intended for internal use of the Lie Group Solver. Nonetheless, they contain useful information in their docstrings on the algorithms implemented for the various heuristics.
- sympy.solvers.ode.lie_group.lie_heuristic_abaco1_simple(match, comp=False)[source]¶
The first heuristic uses the following four sets of assumptions on \(\xi\) and \(\eta\)
\[\xi = 0, \eta = f(x)\]\[\xi = 0, \eta = f(y)\]\[\xi = f(x), \eta = 0\]\[\xi = f(y), \eta = 0\]The success of this heuristic is determined by algebraic factorisation. For the first assumption \(\xi = 0\) and \(\eta\) to be a function of \(x\), the PDE
\[\frac{\partial \eta}{\partial x} + (\frac{\partial \eta}{\partial y} - \frac{\partial \xi}{\partial x})*h - \frac{\partial \xi}{\partial y}*h^{2} - \xi*\frac{\partial h}{\partial x} - \eta*\frac{\partial h}{\partial y} = 0\]reduces to \(f'(x) - f\frac{\partial h}{\partial y} = 0\) If \(\frac{\partial h}{\partial y}\) is a function of \(x\), then this can usually be integrated easily. A similar idea is applied to the other 3 assumptions as well.
References
E.S Cheb-Terrab, L.G.S Duarte and L.A,C.P da Mota, Computer Algebra Solving of First Order ODEs Using Symmetry Methods, pp. 8
- sympy.solvers.ode.lie_group.lie_heuristic_abaco1_product(match, comp=False)[source]¶
The second heuristic uses the following two assumptions on \(\xi\) and \(\eta\)
\[\eta = 0, \xi = f(x)*g(y)\]\[\eta = f(x)*g(y), \xi = 0\]The first assumption of this heuristic holds good if \(\frac{1}{h^{2}}\frac{\partial^2}{\partial x \partial y}\log(h)\) is separable in \(x\) and \(y\), then the separated factors containing \(x\) is \(f(x)\), and \(g(y)\) is obtained by
\[e^{\int f\frac{\partial}{\partial x}\left(\frac{1}{f*h}\right)\,dy}\]provided \(f\frac{\partial}{\partial x}\left(\frac{1}{f*h}\right)\) is a function of \(y\) only.
The second assumption holds good if \(\frac{dy}{dx} = h(x, y)\) is rewritten as \(\frac{dy}{dx} = \frac{1}{h(y, x)}\) and the same properties of the first assumption satisfies. After obtaining \(f(x)\) and \(g(y)\), the coordinates are again interchanged, to get \(\eta\) as \(f(x)*g(y)\)
References
E.S. Cheb-Terrab, A.D. Roche, Symmetries and First Order ODE Patterns, pp. 7 - pp. 8
- sympy.solvers.ode.lie_group.lie_heuristic_bivariate(match, comp=False)[source]¶
The third heuristic assumes the infinitesimals \(\xi\) and \(\eta\) to be bi-variate polynomials in \(x\) and \(y\). The assumption made here for the logic below is that \(h\) is a rational function in \(x\) and \(y\) though that may not be necessary for the infinitesimals to be bivariate polynomials. The coefficients of the infinitesimals are found out by substituting them in the PDE and grouping similar terms that are polynomials and since they form a linear system, solve and check for non trivial solutions. The degree of the assumed bivariates are increased till a certain maximum value.
References
Lie Groups and Differential Equations pp. 327 - pp. 329
- sympy.solvers.ode.lie_group.lie_heuristic_chi(match, comp=False)[source]¶
The aim of the fourth heuristic is to find the function \(\chi(x, y)\) that satisfies the PDE \(\frac{d\chi}{dx} + h\frac{d\chi}{dx} - \frac{\partial h}{\partial y}\chi = 0\).
This assumes \(\chi\) to be a bivariate polynomial in \(x\) and \(y\). By intuition, \(h\) should be a rational function in \(x\) and \(y\). The method used here is to substitute a general binomial for \(\chi\) up to a certain maximum degree is reached. The coefficients of the polynomials, are calculated by by collecting terms of the same order in \(x\) and \(y\).
After finding \(\chi\), the next step is to use \(\eta = \xi*h + \chi\), to determine \(\xi\) and \(\eta\). This can be done by dividing \(\chi\) by \(h\) which would give \(-\xi\) as the quotient and \(\eta\) as the remainder.
References
E.S Cheb-Terrab, L.G.S Duarte and L.A,C.P da Mota, Computer Algebra Solving of First Order ODEs Using Symmetry Methods, pp. 8
- sympy.solvers.ode.lie_group.lie_heuristic_abaco2_similar(match, comp=False)[source]¶
This heuristic uses the following two assumptions on \(\xi\) and \(\eta\)
\[\eta = g(x), \xi = f(x)\]\[\eta = f(y), \xi = g(y)\]For the first assumption,
First \(\frac{\frac{\partial h}{\partial y}}{\frac{\partial^{2} h}{ \partial yy}}\) is calculated. Let us say this value is A
If this is constant, then \(h\) is matched to the form \(A(x) + B(x)e^{ \frac{y}{C}}\) then, \(\frac{e^{\int \frac{A(x)}{C} \,dx}}{B(x)}\) gives \(f(x)\) and \(A(x)*f(x)\) gives \(g(x)\)
Otherwise \(\frac{\frac{\partial A}{\partial X}}{\frac{\partial A}{ \partial Y}} = \gamma\) is calculated. If
a] \(\gamma\) is a function of \(x\) alone
b] \(\frac{\gamma\frac{\partial h}{\partial y} - \gamma'(x) - \frac{ \partial h}{\partial x}}{h + \gamma} = G\) is a function of \(x\) alone. then, \(e^{\int G \,dx}\) gives \(f(x)\) and \(-\gamma*f(x)\) gives \(g(x)\)
The second assumption holds good if \(\frac{dy}{dx} = h(x, y)\) is rewritten as \(\frac{dy}{dx} = \frac{1}{h(y, x)}\) and the same properties of the first assumption satisfies. After obtaining \(f(x)\) and \(g(x)\), the coordinates are again interchanged, to get \(\xi\) as \(f(x^*)\) and \(\eta\) as \(g(y^*)\)
References
E.S. Cheb-Terrab, A.D. Roche, Symmetries and First Order ODE Patterns, pp. 10 - pp. 12
- sympy.solvers.ode.lie_group.lie_heuristic_function_sum(match, comp=False)[source]¶
This heuristic uses the following two assumptions on \(\xi\) and \(\eta\)
\[\eta = 0, \xi = f(x) + g(y)\]\[\eta = f(x) + g(y), \xi = 0\]The first assumption of this heuristic holds good if
\[\frac{\partial}{\partial y}[(h\frac{\partial^{2}}{ \partial x^{2}}(h^{-1}))^{-1}]\]is separable in \(x\) and \(y\),
The separated factors containing \(y\) is \(\frac{\partial g}{\partial y}\). From this \(g(y)\) can be determined.
The separated factors containing \(x\) is \(f''(x)\).
\(h\frac{\partial^{2}}{\partial x^{2}}(h^{-1})\) equals \(\frac{f''(x)}{f(x) + g(y)}\). From this \(f(x)\) can be determined.
The second assumption holds good if \(\frac{dy}{dx} = h(x, y)\) is rewritten as \(\frac{dy}{dx} = \frac{1}{h(y, x)}\) and the same properties of the first assumption satisfies. After obtaining \(f(x)\) and \(g(y)\), the coordinates are again interchanged, to get \(\eta\) as \(f(x) + g(y)\).
For both assumptions, the constant factors are separated among \(g(y)\) and \(f''(x)\), such that \(f''(x)\) obtained from 3] is the same as that obtained from 2]. If not possible, then this heuristic fails.
References
E.S. Cheb-Terrab, A.D. Roche, Symmetries and First Order ODE Patterns, pp. 7 - pp. 8
- sympy.solvers.ode.lie_group.lie_heuristic_abaco2_unique_unknown(
- match,
- comp=False,
This heuristic assumes the presence of unknown functions or known functions with non-integer powers.
A list of all functions and non-integer powers containing x and y
Loop over each element \(f\) in the list, find \(\frac{\frac{\partial f}{\partial x}}{ \frac{\partial f}{\partial x}} = R\)
If it is separable in \(x\) and \(y\), let \(X\) be the factors containing \(x\). Then
- a] Check if \(\xi = X\) and \(\eta = -\frac{X}{R}\) satisfy the PDE. If yes, then return
\(\xi\) and \(\eta\)
- b] Check if \(\xi = \frac{-R}{X}\) and \(\eta = -\frac{1}{X}\) satisfy the PDE.
If yes, then return \(\xi\) and \(\eta\)
If not, then check if
a] \(\xi = -R,\eta = 1\)
b] \(\xi = 1, \eta = -\frac{1}{R}\)
are solutions.
References
E.S. Cheb-Terrab, A.D. Roche, Symmetries and First Order ODE Patterns, pp. 10 - pp. 12
- sympy.solvers.ode.lie_group.lie_heuristic_abaco2_unique_general(
- match,
- comp=False,
This heuristic finds if infinitesimals of the form \(\eta = f(x)\), \(\xi = g(y)\) without making any assumptions on \(h\).
The complete sequence of steps is given in the paper mentioned below.
References
E.S. Cheb-Terrab, A.D. Roche, Symmetries and First Order ODE Patterns, pp. 10 - pp. 12
- sympy.solvers.ode.lie_group.lie_heuristic_linear(match, comp=False)[source]¶
This heuristic assumes
\(\xi = ax + by + c\) and
\(\eta = fx + gy + h\)
After substituting the following assumptions in the determining PDE, it reduces to
\[f + (g - a)h - bh^{2} - (ax + by + c)\frac{\partial h}{\partial x} - (fx + gy + c)\frac{\partial h}{\partial y}\]Solving the reduced PDE obtained, using the method of characteristics, becomes impractical. The method followed is grouping similar terms and solving the system of linear equations obtained. The difference between the bivariate heuristic is that \(h\) need not be a rational function in this case.
References
E.S. Cheb-Terrab, A.D. Roche, Symmetries and First Order ODE Patterns, pp. 10 - pp. 12
Rational Riccati Solver¶
These functions are intended for internal use to solve a first order Riccati differential equation with atleast one rational particular solution.
- sympy.solvers.ode.riccati.riccati_normal(w, x, b1, b2)[source]¶
Given a solution \(w(x)\) to the equation
\[w'(x) = b_0(x) + b_1(x)*w(x) + b_2(x)*w(x)^2\]and rational function coefficients \(b_1(x)\) and \(b_2(x)\), this function transforms the solution to give a solution \(y(x)\) for its corresponding normal Riccati ODE
\[y'(x) + y(x)^2 = a(x)\]using the transformation
\[y(x) = -b_2(x)*w(x) - b'_2(x)/(2*b_2(x)) - b_1(x)/2\]
- sympy.solvers.ode.riccati.riccati_inverse_normal(y, x, b1, b2, bp=None)[source]¶
Inverse transforming the solution to the normal Riccati ODE to get the solution to the Riccati ODE.
- sympy.solvers.ode.riccati.riccati_reduced(eq, f, x)[source]¶
Convert a Riccati ODE into its corresponding normal Riccati ODE.
- sympy.solvers.ode.riccati.construct_c(num, den, x, poles, muls)[source]¶
Helper function to calculate the coefficients in the c-vector for each pole.
- sympy.solvers.ode.riccati.construct_d(num, den, x, val_inf)[source]¶
Helper function to calculate the coefficients in the d-vector based on the valuation of the function at oo.
- sympy.solvers.ode.riccati.rational_laurent_series(num, den, x, r, m, n)[source]¶
The function computes the Laurent series coefficients of a rational function.
- Parameters:
num: A Poly object that is the numerator of `f(x)`.
den: A Poly object that is the denominator of `f(x)`.
x: The variable of expansion of the series.
r: The point of expansion of the series.
m: Multiplicity of r if r is a pole of `f(x)`. Should
be zero otherwise.
n: Order of the term upto which the series is expanded.
- Returns:
series: A dictionary that has power of the term as key
and coefficient of that term as value.
Below is a basic outline of how the Laurent series of a
rational function \(f(x)\) about \(x_0\) is being calculated -
Substitute \(x + x_0\) in place of \(x\). If \(x_0\)
is a pole of \(f(x)\), multiply the expression by \(x^m\)
where \(m\) is the multiplicity of \(x_0\). Denote the
the resulting expression as g(x). We do this substitution
so that we can now find the Laurent series of g(x) about
\(x = 0\).
We can then assume that the Laurent series of \(g(x)\)
takes the following form -
\[g(x) = \frac{num(x)}{den(x)} = \sum_{m = 0}^{\infty} a_m x^m\]where \(a_m\) denotes the Laurent series coefficients.
Multiply the denominator to the RHS of the equation
and form a recurrence relation for the coefficients \(a_m\).
- sympy.solvers.ode.riccati.compute_m_ybar(x, poles, choice, N)[source]¶
Helper function to calculate -
1. m - The degree bound for the polynomial solution that must be found for the auxiliary differential equation.
2. ybar - Part of the solution which can be computed using the poles, c and d vectors.
- sympy.solvers.ode.riccati.solve_aux_eq(numa, dena, numy, deny, x, m)[source]¶
Helper function to find a polynomial solution of degree m for the auxiliary differential equation.
- sympy.solvers.ode.riccati.remove_redundant_sols(sol1, sol2, x)[source]¶
Helper function to remove redundant solutions to the differential equation.
- sympy.solvers.ode.riccati.get_gen_sol_from_part_sol(part_sols, a, x)[source]¶
” Helper function which computes the general solution for a Riccati ODE from its particular solutions.
There are 3 cases to find the general solution from the particular solutions for a Riccati ODE depending on the number of particular solution(s) we have - 1, 2 or 3.
For more information, see Section 6 of “Methods of Solution of the Riccati Differential Equation” by D. R. Haaheim and F. M. Stein
System of ODEs¶
These functions are intended for internal use by
dsolve()
for system of differential equations.
- sympy.solvers.ode.ode._linear_2eq_order1_type6(x, y, t, r, eq)[source]¶
The equations of this type of ode are .
\[x' = f(t) x + g(t) y\]\[y' = a [f(t) + a h(t)] x + a [g(t) - h(t)] y\]This is solved by first multiplying the first equation by \(-a\) and adding it to the second equation to obtain
\[y' - a x' = -a h(t) (y - a x)\]Setting \(U = y - ax\) and integrating the equation we arrive at
\[y - ax = C_1 e^{-a \int h(t) \,dt}\]and on substituting the value of y in first equation give rise to first order ODEs. After solving for \(x\), we can obtain \(y\) by substituting the value of \(x\) in second equation.
- sympy.solvers.ode.ode._linear_2eq_order1_type7(x, y, t, r, eq)[source]¶
The equations of this type of ode are .
\[x' = f(t) x + g(t) y\]\[y' = h(t) x + p(t) y\]Differentiating the first equation and substituting the value of \(y\) from second equation will give a second-order linear equation
\[g x'' - (fg + gp + g') x' + (fgp - g^{2} h + f g' - f' g) x = 0\]This above equation can be easily integrated if following conditions are satisfied.
\(fgp - g^{2} h + f g' - f' g = 0\)
\(fgp - g^{2} h + f g' - f' g = ag, fg + gp + g' = bg\)
If first condition is satisfied then it is solved by current dsolve solver and in second case it becomes a constant coefficient differential equation which is also solved by current solver.
Otherwise if the above condition fails then, a particular solution is assumed as \(x = x_0(t)\) and \(y = y_0(t)\) Then the general solution is expressed as
\[x = C_1 x_0(t) + C_2 x_0(t) \int \frac{g(t) F(t) P(t)}{x_0^{2}(t)} \,dt\]\[y = C_1 y_0(t) + C_2 [\frac{F(t) P(t)}{x_0(t)} + y_0(t) \int \frac{g(t) F(t) P(t)}{x_0^{2}(t)} \,dt]\]where C1 and C2 are arbitrary constants and
\[F(t) = e^{\int f(t) \,dt}, P(t) = e^{\int p(t) \,dt}\]
- sympy.solvers.ode.systems.linear_ode_to_matrix(eqs, funcs, t, order)[source]¶
Convert a linear system of ODEs to matrix form
- Parameters:
eqs : list of SymPy expressions or equalities
The equations as expressions (assumed equal to zero).
funcs : list of applied functions
The dependent variables of the system of ODEs.
t : symbol
The independent variable.
order : int
The order of the system of ODEs.
- Returns:
The tuple
(As, b)
whereAs
is a tuple of matrices andb
is thethe matrix representing the rhs of the matrix equation.
- Raises:
ODEOrderError
When the system of ODEs have an order greater than what was specified
ODENonlinearError
When the system of ODEs is nonlinear
Explanation
Express a system of linear ordinary differential equations as a single matrix differential equation [1]. For example the system \(x' = x + y + 1\) and \(y' = x - y\) can be represented as
\[A_1 X' = A_0 X + b\]where \(A_1\) and \(A_0\) are \(2 \times 2\) matrices and \(b\), \(X\) and \(X'\) are \(2 \times 1\) matrices with \(X = [x, y]^T\).
Higher-order systems are represented with additional matrices e.g. a second-order system would look like
\[A_2 X'' = A_1 X' + A_0 X + b\]Examples
>>> from sympy import Function, Symbol, Matrix, Eq >>> from sympy.solvers.ode.systems import linear_ode_to_matrix >>> t = Symbol('t') >>> x = Function('x') >>> y = Function('y')
We can create a system of linear ODEs like
>>> eqs = [ ... Eq(x(t).diff(t), x(t) + y(t) + 1), ... Eq(y(t).diff(t), x(t) - y(t)), ... ] >>> funcs = [x(t), y(t)] >>> order = 1 # 1st order system
Now
linear_ode_to_matrix
can represent this as a matrix differential equation.>>> (A1, A0), b = linear_ode_to_matrix(eqs, funcs, t, order) >>> A1 Matrix([ [1, 0], [0, 1]]) >>> A0 Matrix([ [1, 1], [1, -1]]) >>> b Matrix([ [1], [0]])
The original equations can be recovered from these matrices:
>>> eqs_mat = Matrix([eq.lhs - eq.rhs for eq in eqs]) >>> X = Matrix(funcs) >>> A1 * X.diff(t) - A0 * X - b == eqs_mat True
If the system of equations has a maximum order greater than the order of the system specified, a ODEOrderError exception is raised.
>>> eqs = [Eq(x(t).diff(t, 2), x(t).diff(t) + x(t)), Eq(y(t).diff(t), y(t) + x(t))] >>> linear_ode_to_matrix(eqs, funcs, t, 1) Traceback (most recent call last): ... ODEOrderError: Cannot represent system in 1-order form
If the system of equations is nonlinear, then ODENonlinearError is raised.
>>> eqs = [Eq(x(t).diff(t), x(t) + y(t)), Eq(y(t).diff(t), y(t)**2 + x(t))] >>> linear_ode_to_matrix(eqs, funcs, t, 1) Traceback (most recent call last): ... ODENonlinearError: The system of ODEs is nonlinear.
See also
linear_eq_to_matrix
for systems of linear algebraic equations.
References
- sympy.solvers.ode.systems.canonical_odes(eqs, funcs, t)[source]¶
Function that solves for highest order derivatives in a system
- Parameters:
eqs : List
List of the ODEs
funcs : List
List of dependent variables
t : Symbol
Independent variable
- Returns:
List
Explanation
This function inputs a system of ODEs and based on the system, the dependent variables and their highest order, returns the system in the following form:
\[X'(t) = A(t) X(t) + b(t)\]Here, \(X(t)\) is the vector of dependent variables of lower order, \(A(t)\) is the coefficient matrix, \(b(t)\) is the non-homogeneous term and \(X'(t)\) is the vector of dependent variables in their respective highest order. We use the term canonical form to imply the system of ODEs which is of the above form.
If the system passed has a non-linear term with multiple solutions, then a list of systems is returned in its canonical form.
Examples
>>> from sympy import symbols, Function, Eq, Derivative >>> from sympy.solvers.ode.systems import canonical_odes >>> f, g = symbols("f g", cls=Function) >>> x, y = symbols("x y") >>> funcs = [f(x), g(x)] >>> eqs = [Eq(f(x).diff(x) - 7*f(x), 12*g(x)), Eq(g(x).diff(x) + g(x), 20*f(x))]
>>> canonical_eqs = canonical_odes(eqs, funcs, x) >>> canonical_eqs [[Eq(Derivative(f(x), x), 7*f(x) + 12*g(x)), Eq(Derivative(g(x), x), 20*f(x) - g(x))]]
>>> system = [Eq(Derivative(f(x), x)**2 - 2*Derivative(f(x), x) + 1, 4), Eq(-y*f(x) + Derivative(g(x), x), 0)]
>>> canonical_system = canonical_odes(system, funcs, x) >>> canonical_system [[Eq(Derivative(f(x), x), -1), Eq(Derivative(g(x), x), y*f(x))], [Eq(Derivative(f(x), x), 3), Eq(Derivative(g(x), x), y*f(x))]]
- sympy.solvers.ode.systems.linodesolve_type(A, t, b=None)[source]¶
Helper function that determines the type of the system of ODEs for solving with
sympy.solvers.ode.systems.linodesolve()
- Parameters:
A : Matrix
Coefficient matrix of the system of ODEs
b : Matrix or None
Non-homogeneous term of the system. The default value is None. If this argument is None, then the system is assumed to be homogeneous.
- Returns:
Dict
- Raises:
NotImplementedError
When the coefficient matrix does not have a commutative antiderivative
Explanation
This function takes in the coefficient matrix and/or the non-homogeneous term and returns the type of the equation that can be solved by
sympy.solvers.ode.systems.linodesolve()
.If the system is constant coefficient homogeneous, then “type1” is returned
If the system is constant coefficient non-homogeneous, then “type2” is returned
If the system is non-constant coefficient homogeneous, then “type3” is returned
If the system is non-constant coefficient non-homogeneous, then “type4” is returned
If the system has a non-constant coefficient matrix which can be factorized into constant coefficient matrix, then “type5” or “type6” is returned for when the system is homogeneous or non-homogeneous respectively.
Note that, if the system of ODEs is of “type3” or “type4”, then along with the type, the commutative antiderivative of the coefficient matrix is also returned.
If the system cannot be solved by
sympy.solvers.ode.systems.linodesolve()
, then NotImplementedError is raised.Examples
>>> from sympy import symbols, Matrix >>> from sympy.solvers.ode.systems import linodesolve_type >>> t = symbols("t") >>> A = Matrix([[1, 1], [2, 3]]) >>> b = Matrix([t, 1])
>>> linodesolve_type(A, t) {'antiderivative': None, 'type_of_equation': 'type1'}
>>> linodesolve_type(A, t, b=b) {'antiderivative': None, 'type_of_equation': 'type2'}
>>> A_t = Matrix([[1, t], [-t, 1]])
>>> linodesolve_type(A_t, t) {'antiderivative': Matrix([ [ t, t**2/2], [-t**2/2, t]]), 'type_of_equation': 'type3'}
>>> linodesolve_type(A_t, t, b=b) {'antiderivative': Matrix([ [ t, t**2/2], [-t**2/2, t]]), 'type_of_equation': 'type4'}
>>> A_non_commutative = Matrix([[1, t], [t, -1]]) >>> linodesolve_type(A_non_commutative, t) Traceback (most recent call last): ... NotImplementedError: The system does not have a commutative antiderivative, it cannot be solved by linodesolve.
See also
linodesolve
Function for which linodesolve_type gets the information
- sympy.solvers.ode.systems.matrix_exp_jordan_form(A, t)[source]¶
Matrix exponential \(\exp(A*t)\) for the matrix A and scalar t.
- Parameters:
A : Matrix
The matrix \(A\) in the expression \(\exp(A*t)\)
t : Symbol
The independent variable
Explanation
Returns the Jordan form of the \(\exp(A*t)\) along with the matrix \(P\) such that:
\[\exp(A*t) = P * expJ * P^{-1}\]Examples
>>> from sympy import Matrix, Symbol >>> from sympy.solvers.ode.systems import matrix_exp, matrix_exp_jordan_form >>> t = Symbol('t')
We will consider a 2x2 defective matrix. This shows that our method works even for defective matrices.
>>> A = Matrix([[1, 1], [0, 1]])
It can be observed that this function gives us the Jordan normal form and the required invertible matrix P.
>>> P, expJ = matrix_exp_jordan_form(A, t)
Here, it is shown that P and expJ returned by this function is correct as they satisfy the formula: P * expJ * P_inverse = exp(A*t).
>>> P * expJ * P.inv() == matrix_exp(A, t) True
References
- sympy.solvers.ode.systems.matrix_exp(A, t)[source]¶
Matrix exponential \(\exp(A*t)\) for the matrix
A
and scalart
.- Parameters:
A : Matrix
The matrix \(A\) in the expression \(\exp(A*t)\)
t : Symbol
The independent variable
Explanation
This functions returns the \(\exp(A*t)\) by doing a simple matrix multiplication:
\[\exp(A*t) = P * expJ * P^{-1}\]where \(expJ\) is \(\exp(J*t)\). \(J\) is the Jordan normal form of \(A\) and \(P\) is matrix such that:
\[A = P * J * P^{-1}\]The matrix exponential \(\exp(A*t)\) appears in the solution of linear differential equations. For example if \(x\) is a vector and \(A\) is a matrix then the initial value problem
\[\frac{dx(t)}{dt} = A \times x(t), x(0) = x0\]has the unique solution
\[x(t) = \exp(A t) x0\]Examples
>>> from sympy import Symbol, Matrix, pprint >>> from sympy.solvers.ode.systems import matrix_exp >>> t = Symbol('t')
We will consider a 2x2 matrix for comupting the exponential
>>> A = Matrix([[2, -5], [2, -4]]) >>> pprint(A) [2 -5] [ ] [2 -4]
Now, exp(A*t) is given as follows:
>>> pprint(matrix_exp(A, t)) [ -t -t -t ] [3*e *sin(t) + e *cos(t) -5*e *sin(t) ] [ ] [ -t -t -t ] [ 2*e *sin(t) - 3*e *sin(t) + e *cos(t)]
See also
matrix_exp_jordan_form
For exponential of Jordan normal form
References
- sympy.solvers.ode.systems.linodesolve(
- A,
- t,
- b=None,
- B=None,
- type='auto',
- doit=False,
- tau=None,
System of n equations linear first-order differential equations
- Parameters:
A : Matrix
Coefficient matrix of the system of linear first order ODEs.
t : Symbol
Independent variable in the system of ODEs.
b : Matrix or None
Non-homogeneous term in the system of ODEs. If None is passed, a homogeneous system of ODEs is assumed.
B : Matrix or None
Antiderivative of the coefficient matrix. If the antiderivative is not passed and the solution requires the term, then the solver would compute it internally.
type : String
Type of the system of ODEs passed. Depending on the type, the solution is evaluated. The type values allowed and the corresponding system it solves are: “type1” for constant coefficient homogeneous “type2” for constant coefficient non-homogeneous, “type3” for non-constant coefficient homogeneous, “type4” for non-constant coefficient non-homogeneous, “type5” and “type6” for non-constant coefficient homogeneous and non-homogeneous systems respectively where the coefficient matrix can be factorized to a constant coefficient matrix. The default value is “auto” which will let the solver decide the correct type of the system passed.
doit : Boolean
Evaluate the solution if True, default value is False
tau: Expression
Used to substitute for the value of \(t\) after we get the solution of the system.
- Returns:
List
- Raises:
ValueError
This error is raised when the coefficient matrix, non-homogeneous term or the antiderivative, if passed, are not a matrix or do not have correct dimensions
NonSquareMatrixError
When the coefficient matrix or its antiderivative, if passed is not a square matrix
NotImplementedError
If the coefficient matrix does not have a commutative antiderivative
Explanation
This solver solves the system of ODEs of the following form:
\[X'(t) = A(t) X(t) + b(t)\]Here, \(A(t)\) is the coefficient matrix, \(X(t)\) is the vector of n independent variables, \(b(t)\) is the non-homogeneous term and \(X'(t)\) is the derivative of \(X(t)\)
Depending on the properties of \(A(t)\) and \(b(t)\), this solver evaluates the solution differently.
When \(A(t)\) is constant coefficient matrix and \(b(t)\) is zero vector i.e. system is homogeneous, the system is “type1”. The solution is:
\[X(t) = \exp(A t) C\]Here, \(C\) is a vector of constants and \(A\) is the constant coefficient matrix.
When \(A(t)\) is constant coefficient matrix and \(b(t)\) is non-zero i.e. system is non-homogeneous, the system is “type2”. The solution is:
\[X(t) = e^{A t} ( \int e^{- A t} b \,dt + C)\]When \(A(t)\) is coefficient matrix such that its commutative with its antiderivative \(B(t)\) and \(b(t)\) is a zero vector i.e. system is homogeneous, the system is “type3”. The solution is:
\[X(t) = \exp(B(t)) C\]When \(A(t)\) is commutative with its antiderivative \(B(t)\) and \(b(t)\) is non-zero i.e. system is non-homogeneous, the system is “type4”. The solution is:
\[X(t) = e^{B(t)} ( \int e^{-B(t)} b(t) \,dt + C)\]When \(A(t)\) is a coefficient matrix such that it can be factorized into a scalar and a constant coefficient matrix:
\[A(t) = f(t) * A\]Where \(f(t)\) is a scalar expression in the independent variable \(t\) and \(A\) is a constant matrix, then we can do the following substitutions:
\[tau = \int f(t) dt, X(t) = Y(tau), b(t) = b(f^{-1}(tau))\]Here, the substitution for the non-homogeneous term is done only when its non-zero. Using these substitutions, our original system becomes:
\[Y'(tau) = A * Y(tau) + b(tau)/f(tau)\]The above system can be easily solved using the solution for “type1” or “type2” depending on the homogeneity of the system. After we get the solution for \(Y(tau)\), we substitute the solution for \(tau\) as \(t\) to get back \(X(t)\)
\[X(t) = Y(tau)\]Systems of “type5” and “type6” have a commutative antiderivative but we use this solution because its faster to compute.
The final solution is the general solution for all the four equations since a constant coefficient matrix is always commutative with its antidervative.
An additional feature of this function is, if someone wants to substitute for value of the independent variable, they can pass the substitution \(tau\) and the solution will have the independent variable substituted with the passed expression(\(tau\)).
Examples
To solve the system of ODEs using this function directly, several things must be done in the right order. Wrong inputs to the function will lead to incorrect results.
>>> from sympy import symbols, Function, Eq >>> from sympy.solvers.ode.systems import canonical_odes, linear_ode_to_matrix, linodesolve, linodesolve_type >>> from sympy.solvers.ode.subscheck import checkodesol >>> f, g = symbols("f, g", cls=Function) >>> x, a = symbols("x, a") >>> funcs = [f(x), g(x)] >>> eqs = [Eq(f(x).diff(x) - f(x), a*g(x) + 1), Eq(g(x).diff(x) + g(x), a*f(x))]
Here, it is important to note that before we derive the coefficient matrix, it is important to get the system of ODEs into the desired form. For that we will use
sympy.solvers.ode.systems.canonical_odes()
.>>> eqs = canonical_odes(eqs, funcs, x) >>> eqs [[Eq(Derivative(f(x), x), a*g(x) + f(x) + 1), Eq(Derivative(g(x), x), a*f(x) - g(x))]]
Now, we will use
sympy.solvers.ode.systems.linear_ode_to_matrix()
to get the coefficient matrix and the non-homogeneous term if it is there.>>> eqs = eqs[0] >>> (A1, A0), b = linear_ode_to_matrix(eqs, funcs, x, 1) >>> A = A0
We have the coefficient matrices and the non-homogeneous term ready. Now, we can use
sympy.solvers.ode.systems.linodesolve_type()
to get the information for the system of ODEs to finally pass it to the solver.>>> system_info = linodesolve_type(A, x, b=b) >>> sol_vector = linodesolve(A, x, b=b, B=system_info['antiderivative'], type=system_info['type_of_equation'])
Now, we can prove if the solution is correct or not by using
sympy.solvers.ode.checkodesol()
>>> sol = [Eq(f, s) for f, s in zip(funcs, sol_vector)] >>> checkodesol(eqs, sol) (True, [0, 0])
We can also use the doit method to evaluate the solutions passed by the function.
>>> sol_vector_evaluated = linodesolve(A, x, b=b, type="type2", doit=True)
Now, we will look at a system of ODEs which is non-constant.
>>> eqs = [Eq(f(x).diff(x), f(x) + x*g(x)), Eq(g(x).diff(x), -x*f(x) + g(x))]
The system defined above is already in the desired form, so we do not have to convert it.
>>> (A1, A0), b = linear_ode_to_matrix(eqs, funcs, x, 1) >>> A = A0
A user can also pass the commutative antiderivative required for type3 and type4 system of ODEs. Passing an incorrect one will lead to incorrect results. If the coefficient matrix is not commutative with its antiderivative, then
sympy.solvers.ode.systems.linodesolve_type()
raises a NotImplementedError. If it does have a commutative antiderivative, then the function just returns the information about the system.>>> system_info = linodesolve_type(A, x, b=b)
Now, we can pass the antiderivative as an argument to get the solution. If the system information is not passed, then the solver will compute the required arguments internally.
>>> sol_vector = linodesolve(A, x, b=b)
Once again, we can verify the solution obtained.
>>> sol = [Eq(f, s) for f, s in zip(funcs, sol_vector)] >>> checkodesol(eqs, sol) (True, [0, 0])
See also
linear_ode_to_matrix
Coefficient matrix computation function
canonical_odes
System of ODEs representation change
linodesolve_type
Getting information about systems of ODEs to pass in this solver
- sympy.solvers.ode.ode._nonlinear_2eq_order1_type1(x, y, t, eq)[source]¶
Equations:
\[x' = x^n F(x,y)\]\[y' = g(y) F(x,y)\]Solution:
\[x = \varphi(y), \int \frac{1}{g(y) F(\varphi(y),y)} \,dy = t + C_2\]where
if \(n \neq 1\)
\[\varphi = [C_1 + (1-n) \int \frac{1}{g(y)} \,dy]^{\frac{1}{1-n}}\]if \(n = 1\)
\[\varphi = C_1 e^{\int \frac{1}{g(y)} \,dy}\]where \(C_1\) and \(C_2\) are arbitrary constants.
- sympy.solvers.ode.ode._nonlinear_2eq_order1_type2(x, y, t, eq)[source]¶
Equations:
\[x' = e^{\lambda x} F(x,y)\]\[y' = g(y) F(x,y)\]Solution:
\[x = \varphi(y), \int \frac{1}{g(y) F(\varphi(y),y)} \,dy = t + C_2\]where
if \(\lambda \neq 0\)
\[\varphi = -\frac{1}{\lambda} log(C_1 - \lambda \int \frac{1}{g(y)} \,dy)\]if \(\lambda = 0\)
\[\varphi = C_1 + \int \frac{1}{g(y)} \,dy\]where \(C_1\) and \(C_2\) are arbitrary constants.
- sympy.solvers.ode.ode._nonlinear_2eq_order1_type3(x, y, t, eq)[source]¶
Autonomous system of general form
\[x' = F(x,y)\]\[y' = G(x,y)\]Assuming \(y = y(x, C_1)\) where \(C_1\) is an arbitrary constant is the general solution of the first-order equation
\[F(x,y) y'_x = G(x,y)\]Then the general solution of the original system of equations has the form
\[\int \frac{1}{F(x,y(x,C_1))} \,dx = t + C_1\]
- sympy.solvers.ode.ode._nonlinear_2eq_order1_type4(x, y, t, eq)[source]¶
Equation:
\[x' = f_1(x) g_1(y) \phi(x,y,t)\]\[y' = f_2(x) g_2(y) \phi(x,y,t)\]First integral:
\[\int \frac{f_2(x)}{f_1(x)} \,dx - \int \frac{g_1(y)}{g_2(y)} \,dy = C\]where \(C\) is an arbitrary constant.
On solving the first integral for \(x\) (resp., \(y\) ) and on substituting the resulting expression into either equation of the original solution, one arrives at a first-order equation for determining \(y\) (resp., \(x\) ).
- sympy.solvers.ode.ode._nonlinear_2eq_order1_type5(func, t, eq)[source]¶
Clairaut system of ODEs
\[x = t x' + F(x',y')\]\[y = t y' + G(x',y')\]The following are solutions of the system
\((i)\) straight lines:
\[x = C_1 t + F(C_1, C_2), y = C_2 t + G(C_1, C_2)\]where \(C_1\) and \(C_2\) are arbitrary constants;
\((ii)\) envelopes of the above lines;
\((iii)\) continuously differentiable lines made up from segments of the lines \((i)\) and \((ii)\).
- sympy.solvers.ode.ode._nonlinear_3eq_order1_type1(x, y, z, t, eq)[source]¶
Equations:
\[a x' = (b - c) y z, \enspace b y' = (c - a) z x, \enspace c z' = (a - b) x y\]First Integrals:
\[a x^{2} + b y^{2} + c z^{2} = C_1\]\[a^{2} x^{2} + b^{2} y^{2} + c^{2} z^{2} = C_2\]where \(C_1\) and \(C_2\) are arbitrary constants. On solving the integrals for \(y\) and \(z\) and on substituting the resulting expressions into the first equation of the system, we arrives at a separable first-order equation on \(x\). Similarly doing that for other two equations, we will arrive at first order equation on \(y\) and \(z\) too.
References
- sympy.solvers.ode.ode._nonlinear_3eq_order1_type2(x, y, z, t, eq)[source]¶
Equations:
\[a x' = (b - c) y z f(x, y, z, t)\]\[b y' = (c - a) z x f(x, y, z, t)\]\[c z' = (a - b) x y f(x, y, z, t)\]First Integrals:
\[a x^{2} + b y^{2} + c z^{2} = C_1\]\[a^{2} x^{2} + b^{2} y^{2} + c^{2} z^{2} = C_2\]where \(C_1\) and \(C_2\) are arbitrary constants. On solving the integrals for \(y\) and \(z\) and on substituting the resulting expressions into the first equation of the system, we arrives at a first-order differential equations on \(x\). Similarly doing that for other two equations we will arrive at first order equation on \(y\) and \(z\).
References
- sympy.solvers.ode.ode._nonlinear_3eq_order1_type3(x, y, z, t, eq)[source]¶
Equations:
\[x' = c F_2 - b F_3, \enspace y' = a F_3 - c F_1, \enspace z' = b F_1 - a F_2\]where \(F_n = F_n(x, y, z, t)\).
First Integral:
\[a x + b y + c z = C_1,\]where C is an arbitrary constant.
2. If we assume function \(F_n\) to be independent of \(t\),i.e, \(F_n\) = \(F_n (x, y, z)\) Then, on eliminating \(t\) and \(z\) from the first two equation of the system, one arrives at the first-order equation
\[\frac{dy}{dx} = \frac{a F_3 (x, y, z) - c F_1 (x, y, z)}{c F_2 (x, y, z) - b F_3 (x, y, z)}\]where \(z = \frac{1}{c} (C_1 - a x - b y)\)
References
- sympy.solvers.ode.ode._nonlinear_3eq_order1_type4(x, y, z, t, eq)[source]¶
Equations:
\[x' = c z F_2 - b y F_3, \enspace y' = a x F_3 - c z F_1, \enspace z' = b y F_1 - a x F_2\]where \(F_n = F_n (x, y, z, t)\)
First integral:
\[a x^{2} + b y^{2} + c z^{2} = C_1\]where \(C\) is an arbitrary constant.
2. Assuming the function \(F_n\) is independent of \(t\): \(F_n = F_n (x, y, z)\). Then on eliminating \(t\) and \(z\) from the first two equations of the system, one arrives at the first-order equation
\[\frac{dy}{dx} = \frac{a x F_3 (x, y, z) - c z F_1 (x, y, z)} {c z F_2 (x, y, z) - b y F_3 (x, y, z)}\]where \(z = \pm \sqrt{\frac{1}{c} (C_1 - a x^{2} - b y^{2})}\)
References
- sympy.solvers.ode.ode._nonlinear_3eq_order1_type5(x, y, z, t, eq)[source]¶
- \[x' = x (c F_2 - b F_3), \enspace y' = y (a F_3 - c F_1), \enspace z' = z (b F_1 - a F_2)\]
where \(F_n = F_n (x, y, z, t)\) and are arbitrary functions.
First Integral:
\[\left|x\right|^{a} \left|y\right|^{b} \left|z\right|^{c} = C_1\]where \(C\) is an arbitrary constant. If the function \(F_n\) is independent of \(t\), then, by eliminating \(t\) and \(z\) from the first two equations of the system, one arrives at a first-order equation.
References
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 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 fordsolve()
.
checkodesol()
- Checks if an equation is the solution to an ODE.
homogeneous_order()
- Returns the homogeneous order of an expression.
infinitesimals()
- Returns the infinitesimals of the Lie group of point transformations of an ODE, such that it is invariant.
checkinfsol()
- Checks if the given infinitesimals are the actual infinitesimals of a first order ODE.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.
Currently implemented solver methods
The following methods are implemented for solving ordinary differential
equations. See the docstrings of the various 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.
Power series solutions for first order differential equations.
Lie Group method of solving first order differential equations.
2nd order Liouville differential equations.
Power series solutions for second order differential equations at ordinary and regular singular points.
\(n\)th order differential equation that can be solved with algebraic rearrangement and integration.
\(n\)th order linear homogeneous differential equation with constant coefficients.
\(n\)th order linear inhomogeneous differential equation with constant coefficients using the method of undetermined coefficients.
\(n\)th 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 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 \(n\)th 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
HomogeneousCoeffBest
, 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 does not 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. 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.rst
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 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 sympy.core.basic.Basic.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/does not work). Be sure to call your
hint specifically in dsolve()
, that way the test
will not 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
\(n\)th order solver, but if your solver is slow, such 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.
These functions are not intended for end-user use.