SymPy-1-13-中文文档-九-
SymPy 1.13 中文文档(九)
不等式求解器
对于一般情况应使用reduce_inequalities()
。其他函数是专用操作的子类别,将根据需要由reduce_inequalities
内部调用。
注意
对于一个以解代数方式减少一个或多个不等式的单变量友好指南,请参阅减少不等式的解数学含义。
注意
以下一些示例使用poly()
,它只是将表达式转换为多项式;它不改变表达式的数学含义。
sympy.solvers.inequalities.solve_rational_inequalities(eqs)
用有理系数解有理不等式系统。
示例
>>> from sympy.abc import x
>>> from sympy import solve_rational_inequalities, Poly
>>> solve_rational_inequalities([[
... ((Poly(-x + 1), Poly(1, x)), '>='),
... ((Poly(-x + 1), Poly(1, x)), '<=')]])
{1}
>>> solve_rational_inequalities([[
... ((Poly(x), Poly(1, x)), '!='),
... ((Poly(-x + 1), Poly(1, x)), '>=')]])
Union(Interval.open(-oo, 0), Interval.Lopen(0, 1))
参见
solve_poly_inequality
sympy.solvers.inequalities.solve_poly_inequality(poly, rel)
用有理系数解多项式不等式。
示例
>>> from sympy import solve_poly_inequality, Poly
>>> from sympy.abc import x
>>> solve_poly_inequality(Poly(x, x, domain='ZZ'), '==')
[{0}]
>>> solve_poly_inequality(Poly(x**2 - 1, x, domain='ZZ'), '!=')
[Interval.open(-oo, -1), Interval.open(-1, 1), Interval.open(1, oo)]
>>> solve_poly_inequality(Poly(x**2 - 1, x, domain='ZZ'), '==')
[{-1}, {1}]
参见
solve_poly_inequalities
sympy.solvers.inequalities.solve_poly_inequalities(polys)
用有理系数解多项式不等式。
示例
>>> from sympy import Poly
>>> from sympy.solvers.inequalities import solve_poly_inequalities
>>> from sympy.abc import x
>>> solve_poly_inequalities(((
... Poly(x**2 - 3), ">"), (
... Poly(-x**2 + 1), ">")))
Union(Interval.open(-oo, -sqrt(3)), Interval.open(-1, 1), Interval.open(sqrt(3), oo))
sympy.solvers.inequalities.reduce_rational_inequalities(exprs, gen, relational=True)
用有理系数解有理不等式系统。
示例
>>> from sympy import Symbol
>>> from sympy.solvers.inequalities import reduce_rational_inequalities
>>> x = Symbol('x', real=True)
>>> reduce_rational_inequalities([[x**2 <= 0]], x)
Eq(x, 0)
>>> reduce_rational_inequalities([[x + 2 > 0]], x)
-2 < x
>>> reduce_rational_inequalities([[(x + 2, ">")]], x)
-2 < x
>>> reduce_rational_inequalities([[x + 2]], x)
Eq(x, -2)
此函数找到非无限解集,因此如果未知符号声明为扩展实数而不是实数,则结果可能包括有限性条件:
>>> y = Symbol('y', extended_real=True)
>>> reduce_rational_inequalities([[y + 2 > 0]], y)
(-2 < y) & (y < oo)
sympy.solvers.inequalities.reduce_abs_inequality(expr, rel, gen)
解嵌套绝对值不等式。
示例
>>> from sympy import reduce_abs_inequality, Abs, Symbol
>>> x = Symbol('x', real=True)
>>> reduce_abs_inequality(Abs(x - 5) - 3, '<', x)
(2 < x) & (x < 8)
>>> reduce_abs_inequality(Abs(x + 2)*3 - 13, '<', x)
(-19/3 < x) & (x < 7/3)
参见
reduce_abs_inequalities
sympy.solvers.inequalities.reduce_abs_inequalities(exprs, gen)
解嵌套绝对值不等式系统。
示例
>>> from sympy import reduce_abs_inequalities, Abs, Symbol
>>> x = Symbol('x', extended_real=True)
>>> reduce_abs_inequalities([(Abs(3*x - 5) - 7, '<'),
... (Abs(x + 25) - 13, '>')], x)
(-2/3 < x) & (x < 4) & (((-oo < x) & (x < -38)) | ((-12 < x) & (x < oo)))
>>> reduce_abs_inequalities([(Abs(x - 4) + Abs(3*x - 5) - 7, '<')], x)
(1/2 < x) & (x < 4)
参见
reduce_abs_inequality
sympy.solvers.inequalities.reduce_inequalities(inequalities, symbols=[])
用有理系数解不等式系统。
示例
>>> from sympy.abc import x, y
>>> from sympy import reduce_inequalities
>>> reduce_inequalities(0 <= x + 3, [])
(-3 <= x) & (x < oo)
>>> reduce_inequalities(0 <= x + y*2 - 1, [x])
(x < oo) & (x >= 1 - 2*y)
sympy.solvers.inequalities.solve_univariate_inequality(expr, gen, relational=True, domain=Reals, continuous=False)
解实数单变量不等式。
参数:
表达式:关系型
目标不等式
gen:符号
解不等式的变量
关系型:布尔值
预期输出为关系类型或否
定义域:集合
解方程的定义域
连续: 布尔值
如果已知表达式在给定域上连续(因此不需要调用 continuous_domain()),则返回 True。
异常:
未实现错误
由于
sympy.solvers.solveset.solvify()
的限制,无法确定不等式的解。
注释
目前,由于 sympy.solvers.solveset.solvify()
的限制,我们无法解决所有不等式。此外,对于三角不等式返回的解受其周期间隔的限制。
示例
>>> from sympy import solve_univariate_inequality, Symbol, sin, Interval, S
>>> x = Symbol('x')
>>> solve_univariate_inequality(x**2 >= 4, x)
((2 <= x) & (x < oo)) | ((-oo < x) & (x <= -2))
>>> solve_univariate_inequality(x**2 >= 4, x, relational=False)
Union(Interval(-oo, -2), Interval(2, oo))
>>> domain = Interval(0, S.Infinity)
>>> solve_univariate_inequality(x**2 >= 4, x, False, domain)
Interval(2, oo)
>>> solve_univariate_inequality(sin(x) > 0, x, relational=False)
Interval.open(0, pi)
参见
sympy.solvers.solveset.solvify
solver 返回 solveset 解决方案与 solve 的输出 API
ODE
注意
针对解决 ODE 的初学者友好指南,请参阅代数解法解普通微分方程(ODE)。
用户函数
这些函数被导入到全局命名空间中,使用 from sympy import *
。这些函数(与下面的提示函数不同)旨在供 SymPy 的普通用户使用。
sympy.solvers.ode.dsolve(eq, func=None, hint='default', simplify=True, ics=None, xi=None, eta=None, x0=0, n=6, **kwargs)
解决任何(支持的)普通微分方程和普通微分方程组。
用于单个普通微分方程
当 eq
中方程的数量为一时,它被归类为此类。用法
dsolve(eq, f(x), hint)
-> 使用方法hint
解决普通微分方程eq
关于函数f(x)
。
详细信息
eq
可以是任何支持的普通微分方程(请参见
ode
文档字符串中支持的方法)。这可以是一个Equality
,或者一个假定等于0
的表达式。
f(x)
是一个单变量函数,其在那个变量组成普通微分方程
eq
。在许多情况下,提供这些变量并非必要;它将被自动检测(如果无法检测到,则会引发错误)。
hint
是您希望dsolve
使用的解法。使用
classify_ode(eq, f(x))
可以获取 ODE 的所有可能的提示。默认提示default
将使用classify_ode()
返回的第一个提示。有关更多选项,请参见下面的提示。
simplify
通过
odesimp()
实现简化。查看其文档字符串以获取更多信息。例如,可以关闭此选项,以禁用对func
的解或任意常数的简化。即使启用此选项,解可能包含比 ODE 阶数更多的任意常数。
xi
和eta
是普通微分方程。它们是使微分方程不变的点变换李群的无穷小量。用户可以指定无穷小量的值。如果未指定任何内容,则使用各种启发式方法使用
infinitesimals()
计算xi
和eta
。
ics
是微分方程的初值/边界条件集合。它应该以
{f(x0): x1, f(x).diff(x).subs(x, x2): x3}
的形式给出,依此类推。对于未指定初始条件的幂级数解,假定f(0)
为C0
并计算关于 0 的幂级数解。
x0
是微分方程的幂级数解的点。方程的解需要进行评估。
n
给出了依赖变量的指数,该指数是幂级数的上限。微分方程的解要进行评估。
提示
除了各种求解方法外,您还可以将一些元提示传递给
dsolve()
:
default
:这使用由
classify_ode()
返回的第一个提示。这是dsolve()
的默认参数。
all
:要使
dsolve()
应用所有相关的分类提示,请使用dsolve(ODE, func, hint="all")
。这将返回一个hint:solution
术语的字典。如果提示导致dsolve
抛出NotImplementedError
,则该提示键的值将是引发的异常对象。该字典还将包括一些特殊键:
order
: ODE 的阶数。参见ode_order()
在deutils.py
中的说明。best
: 最简单的提示;best
下面将会返回。best_hint
: 给出由best
给出的解决方案的提示。如果有多个提示产生最佳解决方案,则选择classify_ode()
返回的元组中的第一个提示。default
: 默认情况下会返回的解决方案。这是由classify_ode()
返回的元组中首个提示生成的解决方案。
all_Integral
:这与
all
相同,只是如果提示还有对应的_Integral
提示,它只返回_Integral
提示。如果由于困难或不可能的积分而导致dsolve()
无法继续,这会非常有用。这种元提示还比all
要快得多,因为integrate()
是一个昂贵的例程。
best
:要让
dsolve()
尝试所有方法并返回最简单的方法。这考虑到解决方案是否可在函数中解决,是否包含任何积分类(即无法评估的积分),以及哪个解决方案尺寸最小。另请参阅
classify_ode()
的文档字符串,以获取有关提示的更多信息,以及ode
的文档字符串,以获取所有支持的提示列表。
提示
你可以通过以下方式声明未知函数的导数:
>>> 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)
参见
test_ode.py
中的多个测试,这些测试也用作如何使用dsolve()
的示例集。
dsolve()
总是返回一个Equality
类(除非提示为all
或all_Integral
)。如果可能,它会为正在解决的函数明确解出解。否则,它会返回一个隐式解。任意常数是以
C1
,C2
等命名的符号。因为所有解应在数学上等价,某些提示可能会为 ODE 返回相同的结果。尽管如此,两种不同的提示通常会以不同的格式返回相同的解。还要注意,有时两个不同解中的任意常数值可能不相同,因为一个常数可能已经“吸收”了其他常数。
执行
help(ode.ode_<hintname>)
以获取有关特定提示的更多信息,其中<hintname>
是没有_Integral
的提示名称。
用于常规微分方程系统
用法
dsolve(eq, func)
-> 解决常规微分方程系统 eq
,其中 func
是包括 (x(t)), (y(t)), (z(t)) 在内的函数列表,列表中的函数数量取决于 eq
中提供的方程数量。
详细信息
eq
可以是任何支持的常规微分方程系统,可以是一个Equality
或一个假定为等于0
的表达式。
func
包含x(t)
和y(t)
作为一个变量函数,它们连同其一些导数组成常规微分方程系统eq
。不必提供这个;它将被自动检测(如果无法检测到,则引发错误)。
提示
提示由 classify_sysode 返回的参数组合形成,后续用于形成方法名称。
示例
>>> 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)
解决任何支持的常规微分方程系统
参数:
eqs:列表
要解决的常规微分方程系统
funcs:列表或无
组成常规微分方程系统的依赖变量列表
t:符号或无
常规微分方程系统中的自变量
ics:字典或无
系统的初始边界/条件集
doit:布尔值
如果为 True,则评估解。默认值为 True。如果积分评估花费太多时间和/或不是必需的,则可以设置为 false。
simplify: 布尔值
简化系统的解。默认值为 True。如果简化花费太多时间和/或不是必需的,则可以设置为 false。
返回:
方程的列表列表
引发:
NotImplementedError
当无法通过此函数解决常规微分方程系统时。
ValueError
当传递的参数不符合所需格式时。
解释
此函数接受常规微分方程系统作为输入,确定是否可以通过此函数解决,并在找到解决方案时返回。
此函数可以处理以下情况:1. 线性、一阶、恒定系数齐次常微分方程组 2. 线性、一阶、恒定系数非齐次常微分方程组 3. 线性、一阶、非恒定系数齐次常微分方程组 4. 线性、一阶、非恒定系数非齐次常微分方程组 5. 可以分解成上述四种形式的常微分方程系统的任意隐式系统 6. 可以简化为上述五种系统形式之一的任意高阶线性常微分方程系统。
上述描述的系统类型不受方程数量限制,即此函数可以解决上述类型的系统,而不管系统中方程的数量如何。但是,系统越大,解决系统所需的时间就越长。
此函数返回一个解的列表。每个解都是一个方程列表,其中 LHS 是依赖变量,RHS 是独立变量的表达式。
在非恒定系数类型中,并非所有系统都可以通过此函数求解。只有具有具有可交换的反导数的系数矩阵或者可以进一步分解以便分割系统可能具有具有可交换的反导数的系数矩阵的系统才能被解决。
示例
>>> 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))]]
您还可以传递常微分方程系统的初始条件:
>>> 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)]]
可选地,您可以传递依赖变量和要解决系统的自变量:
>>> 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))]]
让我们来看一个隐式的常微分方程系统:
>>> 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)
返回一个可能的dsolve()
常微分方程的分类元组。
元组是有序的,以便第一项是dsolve()
默认用于解决 ODE 的分类。一般来说,列表中靠近开头的分类比靠近末尾的分类更快地生成更好的解决方案,尽管总会有例外情况。要使dsolve()
使用不同的分类,使用dsolve(ODE, func, hint=<classification>)
。另请参阅dsolve()
的文档字符串以了解可以使用的不同元提示。
如果dict
为真,classify_ode()
将返回hint:match
表达式术语的字典。这是dsolve()
内部使用的。请注意,由于字典的顺序是任意的,这可能与元组的顺序不同。
您可以通过执行help(ode.ode_hintname)
获取关于不同提示的帮助,其中hintname
是提示的名称,不包括_Integral
。
请参阅 allhints
或 ode
的文档字符串,了解从 classify_ode()
返回的所有支持的提示列表。
注意事项
这些是关于提示名称的备注。
_Integral
如果一个分类在末尾有
_Integral
,它将返回一个表达式,其中包含一个未评估的Integral
类。请注意,如果integrate()
无法进行积分,某些提示可能会这样做,尽管仅使用_Integral
会更快。事实上,_Integral
提示总是比其对应的没有_Integral
的提示更快,因为integrate()
是一个昂贵的例程。如果dsolve()
卡住了,这可能是因为integrate()
在一个棘手或不可能的积分上卡住了。尝试使用_Integral
提示或all_Integral
来使其返回结果。请注意,一些提示没有
_Integral
对应项。这是因为在解决这些方法的 ODE 时没有使用integrate()
。例如,具有恒定系数的 n 阶线性齐次 ODE 不需要积分来解决,因此没有nth_linear_homogeneous_constant_coeff_Integrate
提示。您可以通过执行expr.doit()
轻松评估表达式中的任何未评估的Integral
。
序数
一些提示包含序数,例如
1st_linear
。这是为了帮助区分它们与尚未实现的其他提示及其他方法。如果一个提示中包含nth
,例如nth_linear
提示,则这意味着使用的方法适用于任意阶的 ODE。
indep
和 dep
一些提示包含
indep
或dep
一词。这些是指独立变量和依赖函数,分别。例如,如果一个 ODE 是关于 (f(x)) 的,则indep
将指 (x),dep
将指 (f)。
subs
如果提示中有
subs
一词,这意味着通过将subs
词后给出的表达式替换为单个虚拟变量来解决 ODE。这通常是在上述的indep
和dep
的术语中。替换的表达式将仅使用允许作为 Python 对象名称的字符编写,这意味着操作符将被拼写出来。例如,indep
/dep
将被写成indep_div_dep
。
coeff
提示中的
coeff
一词指的是 ODE 中某物的系数,通常是导数项的系数。有关各个方法的更多信息,请参阅单个方法的文档字符串(help(ode)
)。与coefficients
相对应,如undetermined_coefficients
,后者指的是一种方法的常见名称。
_best
具有多种基本解法的方法将为每种子方法提供提示和
_best
元分类。这将评估所有提示并返回最佳结果,使用与正常best
元提示相同的考虑。
示例
>>> 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)
代入sol
到ode
并检查结果是否为0
。
当func
是一个函数,如(f(x)),或者函数列表,如([f(x), g(x)]),当ode
是 ODE 系统时,这将有效。sol
可以是单个解或解列表。每个解可能是一个Equality
,表示解满足的等式,例如Eq(f(x), C1), Eq(f(x) + C1, 0)
;或者简单地是Expr
,例如f(x) - C1
。在大多数情况下,不需要显式地标识函数,但如果原方程中无法推断出函数,则可以通过func
参数提供。
如果传递了一系列解的序列,将使用相同类型的容器为每个解返回结果。
它尝试以下方法,按顺序执行,直到找到零等效:
-
在原方程中用解代替(f)。这仅在为(f)解决
ode
时有效。它会首先尝试解决它,除非solve_for_func == False
。 -
对解的(n)次导数,其中(n)是
ode
的阶数,进行检查,看看是否等于解。这仅适用于精确 ODE。 -
对解的第 1、第 2、…、第(n)阶导数进行操作,每次解决该阶导数的(f)(这总是可能的,因为(f)是线性操作符)。然后将每个导数反向代入
ode
。
此函数返回一个元组。元组中的第一项为True
,如果代换结果为0
,则为False
。元组中的第二项为代换结果。如果第一项为True
,则它应始终为0
。有时,即使表达式恰好等于0
,该函数也会返回False
。这是因为simplify()
无法将表达式简化为0
。如果此函数返回的表达式在恒等时消失,则sol
确实是 ODE 的解。
如果这个函数似乎卡住了,很可能是因为难以简化。
使用此函数来测试,测试元组的第一项。
示例
>>> 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)
如果(g)是齐次的,则返回阶数(n),如果不是齐次的,则返回None
。
确定函数是否是齐次的,并确定其阶数。如果函数(f(x, y, \cdots))是阶数为(n)的齐次函数,那么(f(t x, t y, \cdots) = t^n f(x, y, \cdots))。
如果函数是两个变量的函数(F(x, y)),那么(f)是任何阶齐次的等价于能够将(F(x, y))重写为(G(x/y))或(H(y/x))。这个事实用于解决系数同阶齐次的一阶常微分方程(参见HomogeneousCoeffSubsDepDivIndep
和HomogeneousCoeffSubsIndepDivDep
的文档字符串)。
符号可以是函数,但函数的每个参数必须是符号,并且出现在表达式中的函数参数必须与符号列表中给出的参数匹配。如果声明的函数出现具有不同参数的情况,则返回None
。
示例
>>> 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)
普通微分方程的微小函数(\xi(x,y))和(\eta(x,y))是使得微分方程不变的点变换的 Lie 群的微小量。因此,ODE (y'=f(x,y))将接受一个 Lie 群(x*=X(x,y;\varepsilon)=x+\varepsilon\xi(x,y)),(y=Y(x,y;\varepsilon)=y+\varepsilon\eta(x,y)),使得((y*)'=f(x, y^))。可以对新系统的坐标曲线(r(x,y))和(s(x,y))进行坐标变换,使得此 Lie 群变成平移群,(r*=r)且(s=s+\varepsilon)。它们是新系统的坐标曲线的切线。
考虑变换((x, y) \to (X, Y)),使得微分方程保持不变。(\xi)和(\eta)是变换坐标(X)和(Y)在(\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)| = \eta,]
可以通过解以下 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
解决上述提到的偏微分方程并不是简单的事情,只能通过对(\xi)和(\eta)(启发式)做出智能假设来解决。一旦找到一个微小量,尝试找到更多启发式就会停止。这样做是为了优化求解微分方程的速度。如果需要所有微小量的列表,则应将hint
标记为all
,这将给出所有微小量的完整列表。如果需要找到特定启发式的微小量,则可以将其作为标志传递给hint
。
示例
>>> 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}]
参考文献
- 通过对称群解微分方程,John Starrett,第 1 页 - 第 14 页
sympy.solvers.ode.checkinfsol(eq, infinitesimals, func=None, order=None)
此函数用于检查给定无穷小数是否是给定一阶微分方程的实际无穷小数。 此方法特定于 ODE 的李群解算器。
到目前为止,它只是通过替换偏微分方程中的无穷小数来检查。
[\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]
其中 (\eta) 和 (\xi) 是无穷小数,而 (h(x,y) = \frac{dy}{dx})
无穷小数应以 [{xi(x, y): inf, eta(x, y): inf}]
的字典列表形式给出,对应于函数无穷小数的输出。 它返回形式为 [(True/False, sol)]
的值列表,其中 sol
是在偏微分方程中替换无穷小数后得到的值。 如果为 True
,则 sol
将为 0。
sympy.solvers.ode.constantsimp(expr, constants)
简化带有任意常数的表达式。
此函数专门用于dsolve()
,并不适用于一般用途。
简化是通过将任意常数“吸收”到其他任意常数、数字和它们不独立的符号中来完成的。
所有符号必须以相同的名称并带有数字,例如 C1
、C2
、C3
。 这里的 symbolname
将是 ‘C
’,startnumber
将是 1,endnumber
将是 3。 如果任意常数与变量 x
独立无关,则独立符号将是 x
。 无需指定依赖函数,如 f(x)
,因为它已经包含独立符号 x
。
因为项被“吸收”到任意常数中,并且因为简化后常数被重新编号,所以表达式中的任意常数未必等于返回结果中同名的常数。
如果两个或更多的任意常数相加、相乘或互相提升到幂,则它们首先被一起吸收到单个任意常数中。 然后,如果必要,新常数将与其他项结合。
常数的吸收是有限辅助完成的:
-
在
Add
项中收集了条款,以尝试合并常数,因此 (e^x (C_1 \cos(x) + C_2 \cos(x))) 将简化为 (e^x C_1 \cos(x)); -
具有指数为
Add
的幂会被展开,因此 (e^{C_1 + x}) 将简化为 (C_1 e^x)。
使用constant_renumber()
在简化后重新编号常数,否则常数上可能出现任意数字,例如 (C_1 + C_3 x)。
在罕见情况下,单个常量可以“简化”为两个常量。每个微分方程解应具有与微分方程阶数相同数量的任意常数。这里的结果在技术上是正确的,但例如在表达式中可能有(C_1)和(C_2),而实际上(C_1)等于(C_2)。在这种情况下,请慎重使用,并利用在dsolve()
中使用提示的能力。
示例
>>> 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
``` ## 提示函数
这些函数是由`dsolve()`和其他函数内部使用的。与用户函数不同,上述函数不适用于普通 SymPy 用户的日常使用。相反,应该使用`dsolve()`等函数。尽管如此,这些函数在其文档字符串中包含有关各种 ODE 求解方法的有用信息。因此,它们在此处进行了文档化。
```py
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')
内置不可变序列。
如果未给出参数,则构造函数返回一个空元组。如果指定了可迭代对象,则从可迭代对象的项目初始化元组。
如果参数是一个元组,则返回的值是相同的对象。
sympy.solvers.ode.ode.odesimp(ode, eq, func, hint)
简化 ODE 的解,包括尝试解析func
并运行constantsimp()
。
它可能利用提示返回的解的类型知识来应用额外的简化。
它还尝试在表达式中集成任何Integral
,如果提示不是_Integral
提示。
此函数不应对dsolve()
返回的表达式产生影响,因为dsolve()
已调用odesimp()
,但是单独的提示函数不调用odesimp()
(因为dsolve()
包装器不会)。因此,该函数主要设计用于内部使用。
示例
>>> 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)
将expr
中的任意常量重新编号,以使用newconstants
中给定的符号名称。在此过程中,这将按标准方式重新排序表达式项。
如果未提供newconstants
,则新常量名称将为C1
,C2
等。否则,newconstants
应为可迭代对象,按顺序给出用于常量的新符号。
variables
参数是一个非常量符号的列表。假定expr
中的所有其他自由符号都是常量,并且将进行重新编号。如果未提供variables
,则假定以C
开头的任何编号符号(例如C1
)都是常量。
根据.sort_key()
对符号进行重新编号,因此它们应该按最终打印表达式中出现的顺序编号。请注意,此排序部分基于哈希,因此在不同机器上可能会产生不同的结果。
此函数的结构与constantsimp()
非常相似。
示例
>>> 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
variables
参数指定哪些是常数,以便其他符号不会重新编号:
>>> constant_renumber(expr, [C1, x])
C1*x**2 + C2 + C3*x
newconstants
参数用于指定在替换常数时使用哪些符号:
>>> 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)
返回一个扩展整数,表示 ODE 的解的简单程度。
按从最简单到最不简单的顺序考虑以下事项:
-
sol
被解决为func
。 -
sol
未被解决为func
,但如果传递给 solve(例如,由dsolve(ode, func, simplify=False
返回的解决方案)则可以。 -
如果
sol
未被解决为func
,则基于sol
的长度进行结果计算,由len(str(sol))
计算。 -
如果
sol
有任何未评估的Integral
,这将自动被认为比上述任何一种情况都不简单。
此函数返回一个整数,使得如果解 A 按上述度量标准比解 B 简单,则ode_sol_simplicity(sola, func) < ode_sol_simplicity(solb, func)
。
目前,以下是返回的数字,但如果启发式方法有所改进,可能会改变。仅保证排序。
简单性 | 返回 |
---|---|
sol 解决为func |
-2 |
sol 未解决为func 但可解 |
-1 |
sol 未解决也不可解决为func |
len(str(sol)) |
sol 包含Integral |
oo |
oo
在这里表示 SymPy 的无穷大,应该比任何整数都要大。
如果您已经知道solve()
无法解决sol
,则可以使用trysolving=False
跳过该步骤,这是唯一可能比较慢的步骤。例如,带有simplicity=False
标志的dsolve()
应该会这样做。
如果sol
是解的列表,则如果列表中的最坏解返回oo
,则返回该值,否则返回len(str(sol))
,即整个列表的字符串表示的长度。
示例
此函数旨在作为min
的关键参数传递,例如min(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)
解决具有可解因子的方程。
此函数用于解决具有因子的方程。因子可能是代数或 ode 类型。它将尝试独立解决每个因子。通过调用 dsolve 解决因子。我们将返回解的列表。
示例
>>> 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)
解决一阶确切常微分方程。
一个一阶微分方程如果是某个函数的全微分,则称为确切微分方程。也就是说,微分方程
[P(x, y) ,\partial{}x + Q(x, y) ,\partial{}y = 0]
如果存在某个函数 (F(x, y)),使得 (P(x, y) = \partial{}F/\partial{}x) 且 (Q(x, y) = \partial{}F/\partial{}y),则此微分方程是确切的。可以证明,第一阶微分方程确切的必要和充分条件是 (\partial{}P/\partial{}y = \partial{}Q/\partial{}x)。此时,解将如下所示:
>>> 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
第一偏导数的 (P) 和 (Q) 在单连通区域内存在且连续时。
注:SymPy 目前无法表示对表达式的惰性替换,因此提示 1st_exact_Integral
将返回带有 (dy) 的积分。这应表示您正在解决的函数。
例子
>>> 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)
参考资料
-
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 73
间接测试
class sympy.solvers.ode.single.HomogeneousCoeffBest(ode_problem)
返回从两个提示 1st_homogeneous_coeff_subs_dep_div_indep
和 1st_homogeneous_coeff_subs_indep_div_dep
中得到的常微分方程的最佳解。
这是由 ode_sol_simplicity()
决定的。
查看 HomogeneousCoeffSubsIndepDivDep
和 HomogeneousCoeffSubsDepDivIndep
的文档字符串,了解有关这些提示的更多信息。请注意,没有 ode_1st_homogeneous_coeff_best_Integral
提示。
例子
>>> 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
参考资料
-
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 59
间接测试
class sympy.solvers.ode.single.HomogeneousCoeffSubsDepDivIndep(ode_problem)
使用代换 (u_1 = \frac{\text{
这是一个微分方程
[P(x, y) + Q(x, y) dy/dx = 0]
这样的 (P) 和 (Q) 是同次的,并且同一个顺序。如果函数 (F(x, y)) 是顺序为 (n) 的同次函数,那么 (F(x t, y t) = t^n F(x, y))。等效地,(F(x, y)) 可以重写为 (G(y/x)) 或者 (H(x/y))。详见 homogeneous_order()
的文档字符串。
如果上述微分方程中的系数 (P) 和 (Q) 是同阶次的齐次函数,则可以证明将变量 (x) 和 (u) 分离的代换 (y = u_1 x)(即 (u_1 = y/x))将把微分方程变为关于变量 (x) 和 (u) 的可分离方程。如果 (h(u_1)) 是通过在 (P(x, f(x))) 上进行 (u_1 = f(x)/x) 替换得到的函数,(g(u_2)) 是在微分方程 (P(x, f(x)) + Q(x, f(x)) f'(x) = 0) 中对 (Q(x, f(x))) 进行替换得到的函数,则一般解为:
>>> 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)
|
/
其中 (u_1 h(u_1) + g(u_1) \ne 0),且 (x \ne 0)。
参见 HomogeneousCoeffBest
和 HomogeneousCoeffSubsIndepDivDep
的文档字符串。
示例
>>> 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
参考文献
-
M. Tenenbaum & H. Pollard,《普通微分方程》,多佛 1963 年,第 59 页
间接 doctest
class sympy.solvers.ode.single.HomogeneousCoeffSubsIndepDivDep(ode_problem)
使用替换 (u_2 = \frac{\text{独立变量}}{\text{因变量}}) 解决具有齐次系数的一阶微分方程。
这是一个微分方程
[P(x, y) + Q(x, y) dy/dx = 0]
其中 (P) 和 (Q) 是同阶次的齐次函数。如果函数 (F(x, y)) 是阶次为 (n) 的齐次函数,则满足 (F(x t, y t) = t^n F(x, y))。等价地,(F(x, y)) 可以重写为 (G(y/x)) 或 (H(x/y))。参见 homogeneous_order()
的文档字符串。
如果上述微分方程中的系数 (P) 和 (Q) 是同阶次的齐次函数,则可以证明将变量 (y) 和 (u_2) 分离的代换 (x = u_2 y)(即 (u_2 = x/y))将把微分方程变为关于变量 (y) 和 (u_2) 的可分离方程。如果 (h(u_2)) 是通过在 (P(x, f(x))) 上进行 (u_2 = x/f(x)) 替换得到的函数,(g(u_2)) 是在微分方程 (P(x, f(x)) + Q(x, f(x)) f'(x) = 0) 中对 (Q(x, f(x))) 进行替换得到的函数,则一般解为:
>>> 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
其中 (u_1 g(u_1) + h(u_1) \ne 0),且 (f(x) \ne 0)。
参见 HomogeneousCoeffBest
和 HomogeneousCoeffSubsDepDivIndep
的文档字符串。
示例
>>> 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
参考文献
-
M. Tenenbaum & H. Pollard,《普通微分方程》,多佛 1963 年,第 59 页
间接 doctest
class sympy.solvers.ode.single.FirstLinear(ode_problem)
解一阶线性微分方程。
这些是形如
[dy/dx + P(x) y = Q(x)\text{。}]
这类微分方程可以一般地解出。积分因子 (e^{\int P(x) ,dx}) 将方程变成可分离方程。一般解为:
>>> 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
| | |
\ / /
例子
>>> 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))
参考文献
-
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 92
间接 doctest
class sympy.solvers.ode.single.RationalRiccati(ode_problem)
给出至少有一个有理特解的一阶里卡提微分方程的一般解。
[y' = b_0(x) + b_1(x) y + b_2(x) y²]
其中 (b_0)、(b_1) 和 (b_2) 是 (x) 的有理函数且 (b_2 \ne 0)((b_2 = 0) 会使其成为贝努利方程)。
例子
>>> 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)
参考文献
-
Riccati ODE: 里卡提方程。
-
N. Thieu Vo - Rational and Algebraic Solutions of First-Order Algebraic ODEs: Algorithm 11, pp. 78 -
www3.risc.jku.at/publications/download/risc_5387/PhDThesisThieu.pdf
class sympy.solvers.ode.single.SecondLinearAiry(ode_problem)
给出 Airy 微分方程的解
[\frac{d²y}{dx²} + (a + b x) y(x) = 0]
关于 Airy 特殊函数 airyai 和 airybi 的问题。
例子
>>> 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)
给出贝塞尔微分方程的解
[x² \frac{d²y}{dx²} + x \frac{dy}{dx} y(x) + (x²-n²) y(x)]
若 (n) 是整数,则解的形式为 Eq(f(x), C0 besselj(n,x) + C1 bessely(n,x))
,因为这两个解是线性无关的;否则,若 (n) 是分数,则解的形式为 Eq(f(x), C0 besselj(n,x) + C1 besselj(-n,x))
,这也可以变形为 Eq(f(x), C0 besselj(n,x) + C1 bessely(n,x))
。
例子
>>> 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))
参考文献
class sympy.solvers.ode.single.Bernoulli(ode_problem)
解贝努利微分方程。
这些是形如
[dy/dx + P(x) y = Q(x) y^n\text{,}n \ne 1`\text{。}]
替换 (w = 1/y^{1-n}) 将这种类型的方程转化为线性方程(见 FirstLinear
的文档字符串)。一般解为:
>>> 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 |
|| | | | |
\\ / / / /
注意,当 (n = 1) 时,方程是可分离的(见 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 |
| /
/
例子
>>> 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
参考文献
-
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 95
间接 doctest
class sympy.solvers.ode.single.Liouville(ode_problem)
解二阶 Liouville 微分方程。
里卡提 ODE 的一般形式是
[\frac{d² y}{dx²} + g(y) \left(! \frac{dy}{dx}!\right)² + h(x) \frac{dy}{dx}\text{.}]
一般解是:
>>> 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
| |
/ /
示例
>>> 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) ]
参考资料
-
Goldstein and Braun, “Advanced Methods for the Solution of Differential Equations”, pp. 98
间接的 doctest
class sympy.solvers.ode.single.RiccatiSpecial(ode_problem)
一般的 Riccati 方程形式为
[dy/dx = f(x) y² + g(x) y + h(x)\text{.}]
尽管它没有通解 [1],但“特殊”形式 (dy/dx = a y² - b x^c) 在许多情况下有解 [2]。此例程返回解 (a(dy/dx) = b y² + c y/x + d/x²),通过使用适当的变量变换将其化简为特殊形式,并在 (a) 和 (b) 都不为零且 (c) 或 (d) 其中之一为零时有效。
>>> 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
参考资料
class sympy.solvers.ode.single.NthLinearConstantCoeffHomogeneous(ode_problem)
解出具有常系数的第 (n) 阶线性齐次微分方程。
这是一种形式的方程。
[a_n f^{(n)}(x) + a_{n-1} f^{(n-1)}(x) + \cdots + a_1 f'(x) + a_0 f(x) = 0\text{.}]
这些方程可以通过取特征方程的根来以一般方式解决,然后解将是每个根 (r) 的 (C_n x^i e^{r x}) 项的和,其中 (C_n) 是任意常数,(r) 是特征方程的根,(i) 是从 0 到根的重数减 1 的每个值(例如,重数为 2 的根 3 将创建项 (C_1 e^{3 x} + C_2 x e^{3 x}))。指数通常使用欧拉方程 (e^{I x} = \cos(x) + I \sin(x)) 来展开复根。在具有实系数的多项式中,复根总是成对共轭出现,因此两个根(在简化常数后)将表示为 (e^{a x} \left(C_1 \cos(b x) + C_2 \sin(b x)\right))。
如果 SymPy 无法找到特征方程的确切根,将返回ComplexRootOf
实例。
>>> 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))))
请注意,由于此方法不涉及积分,因此没有 nth_linear_constant_coeff_homogeneous_Integral
提示。
示例
>>> 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
参考资料
-
线性微分方程部分:非齐次方程与恒定系数
-
M. Tenenbaum & H. Pollard, “Ordinary Differential Equations”, Dover 1963, pp. 211
间接的 doctest
class sympy.solvers.ode.single.NthLinearConstantCoeffUndeterminedCoefficients(ode_problem)
使用特解法解决具有常系数的( n )阶线性微分方程。
该方法适用于形如
[a_n f^{(n)}(x) + a_{n-1} f^{(n-1)}(x) + \cdots + a_1 f'(x) + a_0 f(x) = P(x)\text{,}]
其中( P(x) )是一个具有有限个线性独立导数的函数。
符合这一要求的函数是形如( a x^i e^{b x} \sin(c x + d) )或( a x^i e^{b x} \cos(c x + d) )的有限和函数,其中( i )是非负整数,( a )、( b )、( c )和( d )是常数。例如,任何( x )的多项式,像( x² e^{2 x} )、( x \sin(x) )和( e^x \cos(x) )等函数都可以使用。( \sin ) 和 ( \cos ) 的乘积具有有限数量的导数,因为它们可以展开为( \sin(a x) )和( \cos(b x) )项。然而,SymPy 目前无法执行该展开,因此您需要手动将表达式重写为以上形式才能使用此方法。因此,例如,您需要手动将( \sin²(x) )转换为( (1 + \cos(2 x))/2 )才能正确地应用特解法。
该方法通过从表达式及其所有线性独立导数创建试验函数,并将它们代入原始 ODE 中来工作。每个项的系数将成为一组线性方程,这些方程将被解出并代入,从而得到解。如果任何试验函数与齐次方程的解线性相关,则它们将乘以足够的( x )使它们线性无关。
示例
>>> 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
参考文献
-
M. Tenenbaum & H. Pollard,《Ordinary Differential Equations》,Dover 1963 年,第 221 页
间接文档测试
class sympy.solvers.ode.single.NthLinearConstantCoeffVariationOfParameters(ode_problem)
使用参数变化法解决具有常系数的( n )阶线性微分方程。
该方法适用于形如
[f^{(n)}(x) + a_{n-1} f^{(n-1)}(x) + \cdots + a_1 f'(x) + a_0 f(x) = P(x)\text{。}]
该方法通过假设特解的形式来工作
[\sum_{x=1}^{n} c_i(x) y_i(x)\text{,}]
其中( y_i )是齐次方程的第( i )个解。然后使用朗斯基行列式和克莱姆法则解决方程。特解由以下给出
[\sum_{x=1}^n \left( \int \frac{W_i(x)}{W(x)} ,dx \right) y_i(x) \text{,}]
其中( W(x) )是基础系统的朗斯基行列式(齐次方程的( n )个线性独立解的系统),( W_i(x) )是将第( i )列替换为([0, 0, \cdots, 0, P(x)])后的基础系统的朗斯基行列式。
此方法足够通用,可解决任何具有常系数的 (n) 阶非齐次线性微分方程,但有时 SymPy 无法将朗斯基行列式简化到足够程度以进行积分。如果此方法挂起,请尝试使用 nth_linear_constant_coeff_variation_of_parameters_Integral
提示并手动简化积分。此外,当适用时,最好使用 nth_linear_constant_coeff_undetermined_coefficients
,因为它不使用积分,速度更快且更可靠。
警告:在 dsolve()
中使用 simplify=False
和 ‘nth_linear_constant_coeff_variation_of_parameters’ 可能导致其挂起,因为它不会在积分之前尝试简化朗斯基行列式。建议您仅在该方法中使用 simplify=False
和 ‘nth_linear_constant_coeff_variation_of_parameters_Integral’,特别是当齐次方程的解中包含三角函数时。
示例
>>> 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 ///
参考文献
-
M. Tenenbaum & H. Pollard,“Ordinary Differential Equations”,Dover 1963,pp. 233
间接的 doctest
class sympy.solvers.ode.single.NthLinearEulerEqHomogeneous(ode_problem)
解决 (n) 阶线性齐次变系数欧拉-欧拉方程。
这是一个形如 (0 = a_0 f(x) + a_1 x f'(x) + a_2 x² f''(x) \cdots) 的方程。
这些方程可以通过用 (f(x) = x^r) 形式的解代入,并为 (r) 派生特征方程得到一般解。当存在重根时,我们包括额外的形式为 (C_{r k} \ln^k(x) x^r) 的项,其中 (C_{r k}) 是任意积分常数,(r) 是特征方程的根,并且 (k) 取值范围包括 (r) 的重数。在根为复数的情况下,根据欧拉公式的展开,返回形式为 (C_1 x^a \sin(b \log(x)) + C_2 x^a \cos(b \log(x))) 的解。一般解是找到的各项之和。如果 SymPy 无法找到特征方程的确切根,则会返回 ComplexRootOf
实例。
>>> 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)))
请注意,因为此方法不涉及积分,所以没有 nth_linear_euler_eq_homogeneous_Integral
提示。
以下内容仅供内部使用:
-
returns = 'sol'
返回常微分方程的解。 -
returns = 'list'
返回线性独立解的列表,对应于基本解集,用于非齐次解法如参数变化法和待定系数法。请注意,尽管解应该是线性独立的,但此函数并不显式检查。您可以使用assert simplify(wronskian(sollist)) != 0
来检查线性独立性。还需要通过assert len(sollist) == order
。 -
returns = 'both'
,返回一个字典{'sol': <ODE 的解>, 'list': <线性独立解的列表>}
。
示例
>>> 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)
参考资料
-
C. Bender & S. Orszag, “Advanced Mathematical Methods for Scientists and Engineers”, Springer 1999, pp. 12
间接的 doctest
class sympy.solvers.ode.single.NthLinearEulerEqNonhomogeneousVariationOfParameters(ode_problem)
使用参数变化法解(n)阶线性非齐次 Cauchy-Euler 等价普通微分方程。
这是一个形式为(g(x) = a_0 f(x) + a_1 x f'(x) + a_2 x² f''(x) \cdots)的方程。
该方法通过假设特解采用以下形式运作
[\sum_{x=1}^{n} c_i(x) y_i(x) {a_n} {x^n} \text{, }]
其中(y_i)是齐次方程的第(i)个解。然后使用 Wronskian 和 Cramer 法则解出解。特解由以下给定的方程乘以(a_n x^{n})得到。
[\sum_{x=1}^n \left( \int \frac{W_i(x)}{W(x)} , dx \right) y_i(x) \text{, }]
其中(W(x))是基本系统的 Wronskian((n)个线性独立解组成的系统的 Wronskian),(W_i(x))是基本系统的 Wronskian,其中第(i)列被替换为([0, 0, \cdots, 0, \frac{x^{- n}}{a_n} g{\left(x \right)}])。
这种方法足够一般,可以解任意(n)阶非齐次线性微分方程,但有时 SymPy 无法将 Wronskian 简化到足够简单以进行积分。如果这种方法卡住,请尝试使用nth_linear_constant_coeff_variation_of_parameters_Integral
提示,并手动简化积分。此外,当适用时优先使用nth_linear_constant_coeff_undetermined_coefficients
,因为它不使用积分,更快速和可靠。
警告:在dsolve()
中与nth_linear_constant_coeff_variation_of_parameters
一起使用simplify=False
可能会导致卡住,因为在积分之前它不会尝试简化 Wronskian。建议您仅在使用‘simplify=False’与‘nth_linear_constant_coeff_variation_of_parameters_Integral’时才使用此方法,特别是如果齐次方程的解中含有三角函数。
示例
>>> 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)
使用待定系数法解(n)阶线性非齐次 Cauchy-Euler 等价普通微分方程。
这是一个形如( g(x) = a_0 f(x) + a_1 x f'(x) + a_2 x² f''(x) \cdots )的方程。
这些方程可以通过用( x = exp(t) )的形式的解来解决,然后推导形式为( g(exp(t)) = b_0 f(t) + b_1 f'(t) + b_2 f''(t) \cdots )的特征方程。如果( g(exp(t)) )有有限数量的线性独立导数,那么可以通过 n 次线性常系数特解法来解决。
符合此要求的函数是形如( a x^i e^{b x} \sin(c x + d) )或( a x^i e^{b x} \cos(c x + d) )的有限和函数,其中( i )是非负整数,( a ),( b ),( c ),和( d )是常数。例如,任何关于( x )的多项式,如( x² e^{2 x} ),( x \sin(x) ),和( e^x \cos(x) )都可以使用。( \sin )和( \cos )的乘积具有有限数量的导数,因为它们可以展开为( \sin(a x) )和( \cos(b x) )的项。然而,SymPy 目前无法执行这种展开,因此您需要手动将表达式转换为上述形式才能使用此方法。例如,您需要手动将( \sin²(x) )正确转换为( (1 + \cos(2 x))/2 )以便适当地应用待定系数法。
将( x )替换为( exp(t) )后,此方法通过从表达式及其所有线性独立导数中创建一个试探函数,并将它们代入原始的 ODE 中来运行。每个项的系数将形成一组线性方程,可以解出并代入,从而得到解。如果任何试探函数线性依赖于齐次方程的解,则它们将乘以足够的( x )使它们线性独立。
示例
>>> 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)
使用代数和积分解决 n 阶常微分方程。
此类方程的一般形式不存在,这种方程通过将微分视为可逆的代数函数来代数求解。
示例
>>> 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)]
请注意,此求解器可以返回没有任何积分常数的代数解(例如在上面的例子中( f(x) = 0 ))。
class sympy.solvers.ode.single.NthOrderReducible(ode_problem)
解决仅涉及因变量导数的 ODE,使用形式为( f^n(x) = g(x) )的替换。
例如,任何形如( f''(x) = h(f'(x), x) )的二阶 ODE 可以转化为一对一阶 ODEs ( g'(x) = h(g(x), x) ) 和 ( f'(x) = g(x) )。通常,( g )的一阶 ODE 更容易解决。如果给出了( g )的显式解,则通过积分即可找到( f )。
示例
>>> 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)
解决可分离的一阶微分方程。
这是可以写成(P(y) \tfrac{dy}{dx} = Q(x))形式的任何微分方程。解决方案可以通过重新排列项并进行积分找到:(\int P(y) ,dy = \int Q(x) ,dx)。此提示使用sympy.simplify.simplify.separatevars()
作为其后端,因此如果分离方程不被该解算器捕捉到,很可能是该函数的错误。separatevars()
足够智能,能够执行大多数展开和因式分解操作,以将可分离方程(F(x, y))转换为适当的形式(P(x)\cdot{}Q(y))。一般解是:
>>> 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)
| |
/ /
示例
>>> 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
参考文献
- M. Tenenbaum & H. Pollard,《普通微分方程》,多佛出版社 1963 年,第 52 页
间接 doctest
class sympy.solvers.ode.single.AlmostLinear(ode_problem)
解决了一个几乎线性的微分方程。
几乎线性微分方程的一般形式是
[a(x) g'(f(x)) f'(x) + b(x) g(f(x)) + c(x)]
这里的(f(x))是要求解的函数(因变量)。代换(g(f(x)) = u(x))导致了(u(x))的线性微分方程形式(a(x) u' + b(x) u + c(x) = 0)。通过first_linear
提示可以解出(u(x)),然后通过解(g(f(x)) = u(x))找到(f(x))。
示例
>>> 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/]
参见
sympy.solvers.ode.single.FirstLinear
参考文献
- Joel Moses,《符号积分 - 动荡的十年》,ACM 通讯,第 14 卷,第 8 号,1971 年 8 月,第 558 页
class sympy.solvers.ode.single.LinearCoefficients(ode_problem)
解决具有线性系数的微分方程。
具有线性系数的微分方程的一般形式是
[y' + F\left(!\frac{a_1 x + b_1 y + c_1}{a_2 x + b_2 y + c_2}!\right) = 0\text{,}]
其中(a_1)、(b_1)、(c_1)、(a_2)、(b_2)、(c_2)是常数,且(a_1 b_2 - a_2 b_1 \ne 0)。
可通过以下代换解决:
[ \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} ]
此代换将方程简化为齐次微分方程。
示例
>>> 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]
参见
sympy.solvers.ode.single.HomogeneousCoeffBest
,sympy.solvers.ode.single.HomogeneousCoeffSubsIndepDivDep
,sympy.solvers.ode.single.HomogeneousCoeffSubsDepDivIndep
参考文献
- Joel Moses,《符号积分 - 动荡的十年》,ACM 通讯,第 14 卷,第 8 号,1971 年 8 月,第 558 页
class sympy.solvers.ode.single.SeparableReduced(ode_problem)
解决可以化为可分离形式的微分方程。
这个方程的一般形式是
[y' + (y/x) H(x^n y) = 0\text{}.]
这可以通过代入 (u(y) = x^n y) 来解决。然后方程被化简为可分离形式 (\frac{u'}{u (\mathrm{power} - H(u))} - \frac{1}{x} = 0)。
一般解是:
>>> 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))
|
/
示例
>>> 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
参见
sympy.solvers.ode.single.Separable
参考文献
- Joel Moses, “符号积分 - 动荡的十年”, ACM 通信, 卷 14, 第 8 期, 1971 年 8 月, 第 558 页
class sympy.solvers.ode.single.LieGroup(ode_problem)
这个提示实现了通过解第一阶微分方程的 Lie 群方法。其目的是将给定的微分方程从给定坐标系转换为另一个坐标系,在这个坐标系下它变成了一个参数 Lie 群的不变量。转换后的 ODE 可以很容易地通过积分求解。它利用了 sympy.solvers.ode.infinitesimals()
函数返回变换的无穷小量。
坐标 (r) 和 (s) 可以通过解以下偏微分方程来找到。
[\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]
在新的坐标系中,微分方程变为可分离的形式
[\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}}]
找到积分解后,再通过用 (r) 和 (s) 用 (x) 和 (y) 表示的原始坐标系转换回来。
示例
>>> 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 /
参考文献
- 通过对称群解微分方程,John Starrett, 第 1 页至第 14 页
class sympy.solvers.ode.single.SecondHypergeometric(ode_problem)
解二阶线性微分方程。
它计算可以用 2F1, 1F1 或 0F1 超几何函数表示的特殊函数解。
[y'' + A(x) y' + B(x) y = 0\text{,}]
其中 (A) 和 (B) 是有理函数。
这类微分方程的解是非李乌维尔形式的。
给定线性 ODE 可以从 2F1 给出
[(x² - x) y'' + ((a + b + 1) x - c) y' + b a y = 0\text{,}]
其中 {a, b, c} 是任意常数。
注记
算法应该找到形式为
[y = P(x) _pF_q(..; ..;\frac{\alpha x^k + \beta}{\gamma x^k + \delta})\text{,}]
其中 pFq 是 2F1, 1F1 或 0F1 中的任意一个,(P) 是“任意函数”。目前只有 SymPy 中实现了 2F1 的情况,但是其他情况在论文中有描述,可能在将来实现(欢迎贡献!)。
示例
>>> 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)
参考文献
- “非李乌维尔形式的二阶线性 ODE 解”由 L. Chan, E.S. Cheb-Terrab
sympy.solvers.ode.ode.ode_1st_power_series(eq, func, order, match)
幂级数解是一种方法,它给出了微分方程解的泰勒级数展开。
对于一阶微分方程 (\frac{dy}{dx} = h(x, y)),在点 (x = x_{0}) 处存在幂级数解,如果 (h(x, y)) 在 (x_{0}) 处解析。解由以下公式给出
[y(x) = y(x_{0}) + \sum_{n = 1}^{\infty} \frac{F_{n}(x_{0},b)(x - x_{0})^n}{n!},]
其中 (y(x_{0}) = b) 表示 (x_{0}) 初始值时的 (y) 值。为了计算 (F_{n}(x_{0},b)) 的值,遵循以下算法,直到生成所需的项数。
-
(F_1 = h(x_{0}, b))
-
(F_{n+1} = \frac{\partial F_{n}}{\partial x} + \frac{\partial F_{n}}{\partial y}F_{1})
示例
>>> 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
参考文献
- 特拉维斯·W·沃克,《解决一阶微分方程的解析幂级数技术》,第 17, 18 页
sympy.solvers.ode.ode.ode_2nd_power_series_ordinary(eq, func, order, match)
在普通点处,提供了具有多项式系数的二阶齐次微分方程的幂级数解。齐次微分方程的形式为
[P(x)\frac{d²y}{dx²} + Q(x)\frac{dy}{dx} + R(x) y(x) = 0]
假设 (P(x)), (Q(x)) 和 (R(x)) 都是多项式,只需在 (x_{0}) 处存在 (\frac{Q(x)}{P(x)}) 和 (\frac{R(x)}{P(x)}) 即可。通过将 (y) 替换为 (\sum_{n=0}^\infty a_{n}x^{n}),并在微分方程中等价于第 (n) 项,得到递推关系。利用此关系可以生成各种项。
示例
>>> 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 /
参考文献
-
乔治·E·西蒙斯,《应用与历史注解的微分方程》,第 176 - 184 页
sympy.solvers.ode.ode.ode_2nd_power_series_regular(eq, func, order, match)
在正则点处,提供了具有多项式系数的二阶齐次微分方程的幂级数解。二阶齐次微分方程的形式为
[P(x)\frac{d²y}{dx²} + Q(x)\frac{dy}{dx} + R(x) y(x) = 0]
如果 (x - x0\frac{Q(x)}{P(x)}) 和 ((x - x0)^{2}\frac{R(x)}{P(x)}) 在 (x0) 处解析,则称点在 (x0) 处正则奇异。为简单起见,假设 (P(x)), (Q(x)) 和 (R(x)) 是多项式。寻找幂级数解的算法如下:
-
尝试将 ((x - x0)P(x)) 和 (((x - x0)^{2})Q(x)) 表达为关于 (x0) 的幂级数解。找到幂级数展开的常数 (p0) 和 (q0)。
-
解指数方程 (f(m) = m(m - 1) + m*p0 + q0),以获得指数方程的根 (m1) 和 (m2)。
-
如果 (m1 - m2) 是非整数,则存在两个级数解。如果 (m1 = m2),则只存在一个解。如果 (m1 - m2) 是整数,则确认存在一个解。另一个解可能存在也可能不存在。
幂级数解的形式为(x{m}\sum_{n=0}\infty a_{n}x^{n})。 系数由以下递推关系确定。 (a_{n} = -\frac{\sum_{k=0}^{n-1} q_{n-k} + (m + k)p_{n-k}}{f(m + n)})。 如果(m1 - m2)是整数,则可以从递推关系中看出,对于较低的根(m),当(n)等于两个根的差异时,分母变为零。 因此,如果分子不等于零,则存在第二个级数解。
示例
>>> 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
参考文献
- George E. Simmons,“应用与历史注释的微分方程”,p.p 176 - 184
Lie 启发式
这些函数用于 Lie 群求解器的内部使用。 尽管如此,它们在算法文档字符串中包含了关于各种启发式的有用信息。
sympy.solvers.ode.lie_group.lie_heuristic_abaco1_simple(match, comp=False)
第一条启发式使用了对(\xi)和(\eta)的以下四组假设
[\xi = 0, \eta = f(x)][\xi = 0, \eta = f(y)][\xi = f(x), \eta = 0][\xi = f(y), \eta = 0]
这一启发式的成功取决于代数因式分解。 对于第一个假设(\xi = 0)和(\eta)是(x)的函数,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]
减少至[f'(x) - f\frac{\partial h}{\partial y} = 0] 如果(\frac{\partial h}{\partial y})是(x)的函数,则通常可以轻松地进行积分。 对其他 3 个假设也采用了类似的想法。
参考文献
- E.S Cheb-Terrab, L.G.S Duarte and L.A,C.P da Mota, 使用对称方法的计算机代数解决一阶 ODE,pp. 8
sympy.solvers.ode.lie_group.lie_heuristic_abaco1_product(match, comp=False)
第二条启发式使用了对(\xi)和(\eta)的以下两个假设
[\eta = 0, \xi = f(x)g(y)][\eta = f(x)g(y), \xi = 0]
如果(\frac{1}{h^{2}}\frac{\partial²}{\partial x \partial y}\log(h))在(x)和(y)中是可分的,则此启发式的第一个假设成立,其中包含(x)的分离因子是(f(x)),(g(y))由以下四个假设得到
[e^{\int f\frac{\partial}{\partial x}\left(\frac{1}{f*h}\right),dy}]
假设[f\frac{\partial}{\partial x}\left(\frac{1}{f*h}\right))仅是(y)的函数。
如果(\frac{dy}{dx} = h(x, y))重写为(\frac{dy}{dx} = \frac{1}{h(y, x)}),并且第一个假设的相同特性得到满足,则第二个假设成立。 获得(f(x))和(g(y))后,再次交换坐标,得到(\eta)为(f(x)*g(y))
参考文献
- 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)
第三个启发性假设假设 (\xi) 和 (\eta) 是 (x) 和 (y) 的双变量多项式。下面的逻辑假设是,为了使无穷小变量成为双变量多项式,(h) 是 (x) 和 (y) 的有理函数,尽管这对于使无穷小变量成为双变量多项式并不是必要的。通过将它们代入偏微分方程中并分组相似项,可以找到无穷小的系数,这些系数是多项式,并且由于它们形成线性系统,可以解出并检查非平凡解。假设双变量的阶数增加直到某个最大值。
参考文献
- Lie Groups and Differential Equations pp. 327 - pp. 329
sympy.solvers.ode.lie_group.lie_heuristic_chi(match, comp=False)
第四个启发性的目标是找到满足偏微分方程 (\frac{d\chi}{dx} + h\frac{d\chi}{dx} - \frac{\partial h}{\partial y}\chi = 0) 的函数 (\chi(x, y))。
这假设 (\chi) 是 (x) 和 (y) 的双变量多项式。根据直觉,(h) 应该是 (x) 和 (y) 的有理函数。这里使用的方法是将 (\chi) 替换为一个一般的二项式,直到达到某个最大阶数。通过收集相同阶数的 (x) 和 (y) 的项来计算多项式的系数。
找到 (\chi) 后,下一步是使用 (\eta = \xi*h + \chi) 来确定 (\xi) 和 (\eta)。这可以通过将 (\chi) 除以 (h) 来完成,商是 (-\xi),余数是 (\eta)。
参考文献
- 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)
此启发式方法对 (\xi) 和 (\eta) 使用以下两个假设
[\eta = g(x), \xi = f(x)][\eta = f(y), \xi = g(y)]
对于第一个假设,
-
首先计算 (\frac{\frac{\partial h}{\partial y}}{\frac{\partial^{2} h}{ \partial yy}})。假设这个值是 A
-
如果这是一个常数,那么 (h) 就匹配到形式 (A(x) + B(x)e^{ \frac{y}{C}}),然后,(\frac{e^{\int \frac{A(x)}{C} ,dx}}{B(x)}) 给出 (f(x)),(A(x)*f(x)) 给出 (g(x))。
-
否则计算 (\frac{\frac{\partial A}{\partial X}}{\frac{\partial A}{ \partial Y}} = \gamma)。如果
a] (\gamma) 是 (x) 的一个函数。
b] (\frac{\gamma\frac{\partial h}{\partial y} - \gamma'(x) - \frac{ \partial h}{\partial x}}{h + \gamma} = G) 是 (x) 的一个函数。那么,(e^{\int G ,dx}) 给出 (f(x)),(-\gamma*f(x)) 给出 (g(x))。
如果 (\frac{dy}{dx} = h(x, y)) 被重新写成 (\frac{dy}{dx} = \frac{1}{h(y, x)}) 并且第一个假设的相同属性被满足。在获得 (f(x)) 和 (g(x)) 后,再次交换坐标,得到 (\xi) 为 (f(x^)) 和 (\eta) 为 (g(y^))。
参考文献
- 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)
此启发式方法对 (\xi) 和 (\eta) 使用以下两个假设
[\eta = 0, \xi = f(x) + g(y)][\eta = f(x) + g(y), \xi = 0]
此启发式的第一个假设成立,如果
[\frac{\partial}{\partial y}[(h\frac{\partial^{2}}{ \partial x{2}}(h))^{-1}]]
在(x)和(y)中分离,
-
含有(y)的分离因子是(\frac{\partial g}{\partial y})。从中可以确定(g(y))。
-
含有(x)的分离因子是(f''(x))。
-
(h\frac{\partial^{2}}{\partial x{2}}(h))等于(\frac{f''(x)}{f(x) + g(y)})。从此可以确定(f(x))。
如果(\frac{dy}{dx} = h(x, y))被重写为(\frac{dy}{dx} = \frac{1}{h(y, x)})并且第一个假设的相同属性得到满足后,第二个假设成立。在获得(f(x))和(g(y))之后,坐标再次互换,得到(\eta)为(f(x) + g(y))。
对于两种假设,常数因子在(g(y))和(f''(x))之间分离,从而从第 3] 获得的(f''(x))与从第 2] 获得的相同。如果不可能,则此启发式失败。
参考文献
- E.S. Cheb-Terrab, A.D. Roche, 对称性和一阶常微分方程模式, 第 7 - 第 8 页
sympy.solvers.ode.lie_group.lie_heuristic_abaco2_unique_unknown(match, comp=False)
此启发式假设存在未知函数或具有非整数幂的已知函数。
-
所有包含(x)和(y)的函数和非整数幂的列表
-
在列表中循环每个元素(f),找到(\frac{\frac{\partial f}{\partial x}}{ \frac{\partial f}{\partial x}} = R)
如果它在(x)和(y)中是分离的,则让(X)是包含(x)的因子。然后
a] 检查是否(\xi = X)和(\eta = -\frac{X}{R})满足偏微分方程。如果是,则返回。
(\xi)和(\eta)
b] 检查是否(\xi = \frac{-R}{X})和(\eta = -\frac{1}{X})满足偏微分方程。
如果是,则返回(\xi)和(\eta)
如果不是,则检查是否
a] (\xi = -R,\eta = 1)
b] (\xi = 1, \eta = -\frac{1}{R})
是解决方案。
参考文献
- E.S. Cheb-Terrab, A.D. Roche, 对称性和一阶常微分方程模式, 第 10 - 第 12 页
sympy.solvers.ode.lie_group.lie_heuristic_abaco2_unique_general(match, comp=False)
此启发式用于查找形式为(\eta = f(x)),(\xi = g(y))的无穷小量,而不对(h)做任何假设。
给出了下面提到的论文的完整步骤序列。
参考文献
- E.S. Cheb-Terrab, A.D. Roche, 对称性和一阶常微分方程模式, 第 10 - 第 12 页
sympy.solvers.ode.lie_group.lie_heuristic_linear(match, comp=False)
此启发式假设
-
(\xi = ax + by + c) 和
-
(\eta = fx + gy + h)
在将以下假设代入确定性偏微分方程中之后,它化简为
[f + (g - a)h - bh^{2} - (ax + by + c)\frac{\partial h}{\partial x} - (fx + gy + c)\frac{\partial h}{\partial y}]
解决简化的偏微分方程,使用特征线法,变得不切实际。所采用的方法是分组相似项并解决得到的线性方程组。与双变量启发式的区别在于在这种情况下(h)不必是有理函数。
参考文献
- E.S. Cheb-Terrab, A.D. Roche, 对称性和一阶常微分方程模式, 第 10 - 第 12 页
有理里卡蒂求解器
这些函数旨在内部用于解决至少具有一个有理特解的一阶里卡蒂微分方程。
sympy.solvers.ode.riccati.riccati_normal(w, x, b1, b2)
给定方程的解 (w(x))
[w'(x) = b_0(x) + b_1(x)w(x) + b_2(x)w(x)²]
和有理函数系数 (b_1(x)) 和 (b_2(x)),此函数转换解以给出其对应的正常 Riccati ODE 的解 (y(x))
[y'(x) + y(x)² = a(x)]
使用变换
[y(x) = -b_2(x)w(x) - b'_2(x)/(2b_2(x)) - b_1(x)/2]
sympy.solvers.ode.riccati.riccati_inverse_normal(y, x, b1, b2, bp=None)
逆转换正常 Riccati ODE 的解以获取 Riccati ODE 的解。
sympy.solvers.ode.riccati.riccati_reduced(eq, f, x)
将 Riccati ODE 转换为其相应的正常 Riccati ODE。
sympy.solvers.ode.riccati.construct_c(num, den, x, poles, muls)
辅助函数,根据每个极点上的函数估值计算 c 向量中的系数。
sympy.solvers.ode.riccati.construct_d(num, den, x, val_inf)
辅助函数,根据在 oo 处的函数估值计算 d 向量中的系数。
sympy.solvers.ode.riccati.rational_laurent_series(num, den, x, r, m, n)
该函数计算有理函数的 Laurent 级数系数。
参数:
num:是 f(x)
的分子的 Poly 对象。
den:是 f(x)
的分母的 Poly 对象。
x:级数展开的变量。
r:级数展开的点。
m:如果 r 是 f(x)
的极点,则其重数应为零。
否则为零。
n:扩展的术语的顺序。
返回:
series:具有术语幂的字典作为键
并将该项的系数作为值。
下面是如何计算的基本轮廓
在 (x_0) 处计算有理函数 (f(x)) -
- 将 (x_0) 放在 (x) 的位置。如果 (x_0)
(f(x)) 的极点,将表达式乘以 (x^m)
其中 (m) 是 (x_0) 的重数。表示
所得表达式作为 g(x)。我们通过这种替换
这样我们现在可以找到关于 g(x) 的 Laurent 级数
(x = 0)。
- 然后我们可以假设 (g(x)) 的 Laurent 级数
采取以下形式 -
[g(x) = \frac{num(x)}{den(x)} = \sum_{m = 0}^{\infty} a_m x^m]
其中 (a_m) 表示 Laurent 级数系数。
- 将分母乘以方程的 RHS
并形成系数 (a_m) 的递推关系。
sympy.solvers.ode.riccati.compute_m_ybar(x, poles, choice, N)
辅助函数用于计算 -
1. m - 必须为辅助微分方程找到的多项式解的度数限制。
2. ybar - 可使用极点 c 和 d 向量计算的解的一部分。
sympy.solvers.ode.riccati.solve_aux_eq(numa, dena, numy, deny, x, m)
辅助函数,用于找到辅助微分方程的 m 阶多项式解。
sympy.solvers.ode.riccati.remove_redundant_sols(sol1, sol2, x)
辅助函数,用于移除微分方程的冗余解。
sympy.solvers.ode.riccati.get_gen_sol_from_part_sol(part_sols, a, x)
” 从其特解计算 Riccati ODE 的一般解的辅助函数。
根据我们具有的特解数(1、2 或 3),有三种情况可以从特解中找到 Riccati ODE 的一般解。
更多信息,请参见 D. R. Haaheim 和 F. M. Stein 的《Riccati 微分方程的解法》第六部分。
sympy.solvers.ode.riccati.solve_riccati(fx, x, b0, b1, b2, gensol=False)
给出至少有一个有理特解的 Riccati ODE 的特解/一般解。
ODE 系统
这些函数用于内部由dsolve()
处理的常微分方程系统。
sympy.solvers.ode.ode._linear_2eq_order1_type6(x, y, t, r, eq)
这类 ODE 的方程是。
[x' = f(t) x + g(t) y][y' = a [f(t) + a h(t)] x + a [g(t) - h(t)] y]
首先将第一个方程乘以(-a)并加到第二个方程中得到解
[y' - a x' = -a h(t) (y - a x)]
设置(U = y - ax)并积分方程,我们得到
[y - ax = C_1 e^{-a \int h(t) ,dt}]
并且在第一个方程中用 y 的值代替后产生了一阶常微分方程。解出(x)后,我们可以通过将(x)的值代入第二个方程中获得(y)。
sympy.solvers.ode.ode._linear_2eq_order1_type7(x, y, t, r, eq)
这类 ODE 的方程是。
[x' = f(t) x + g(t) y][y' = h(t) x + p(t) y]
对第一个方程进行微分并用第二个方程的值替换将得到一个二阶线性方程
[g x'' - (fg + gp + g') x' + (fgp - g^{2} h + f g' - f' g) x = 0]
如果满足上述条件,则上述方程可以轻松积分。
-
(fgp - g^{2} h + f g' - f' g = 0)
-
(fgp - g^{2} h + f g' - f' g = ag, fg + gp + g' = bg)
如果满足第一个条件,则通过当前的 dsolve 求解器解决,在第二种情况下它变为一个常系数微分方程,也由当前求解器解决。
否则,如果上述条件不满足,则假定特解为(x = x_0(t))和(y = y_0(t))然后一般解表达为
[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]]
其中 C1 和 C2 是任意常数和
[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)
将一个线性常微分方程系统转换为矩阵形式
参数:
eqs : SymPy 表达式或等式的列表
表达式形式的方程(假设为零)。
funcs : 应用函数的列表
方程系统的因变量。
t : 符号
自变量。
order : 整数
常微分方程系统的阶数。
返回:
元组(As, b)
,其中As
是一组矩阵,b
是
表示矩阵方程的 rhs 矩阵。
引发:
ODEOrderError
当常微分方程系统的阶数大于指定的阶数时
ODENonlinearError
当常微分方程系统是非线性时
解释
将一组线性常微分方程表示为单个矩阵微分方程[1]。例如系统(x' = x + y + 1)和(y' = x - y)可以表示为
[A_1 X' = A_0 X + b]
其中(A_1)和(A_0)是(2 \times 2)矩阵,(b),(X)和(X')是(2 \times 1)矩阵,其中(X = [x, y]^T)。
高阶系统用额外的矩阵表示,例如二阶系统看起来像
[A_2 X'' = A_1 X' + A_0 X + b]
示例
>>> 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')
我们可以创建像
>>> 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
现在linear_ode_to_matrix
可以表示这个矩阵微分方程。
>>> (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]])
可以从这些矩阵中恢复原始方程:
>>> eqs_mat = Matrix([eq.lhs - eq.rhs for eq in eqs])
>>> X = Matrix(funcs)
>>> A1 * X.diff(t) - A0 * X - b == eqs_mat
True
如果方程组的最大阶大于指定系统的阶数,则引发 ODEOrderError 异常。
>>> 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
如果方程组是非线性的,则引发 ODENonlinearError。
>>> 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.
参见
linear_eq_to_matrix
用于线性代数方程组。
参考资料
[R884]
en.wikipedia.org/wiki/Matrix_differential_equation
sympy.solvers.ode.systems.canonical_odes(eqs, funcs, t)
解决系统中最高阶导数的函数
参数:
eqs : 列表
ODE 列表
funcs : 列表
依赖变量列表
t : 符号
自变量
返回:
列表
解释
此函数输入一组 ODEs,并根据系统、依赖变量及其最高阶返回以下形式的系统:
[X'(t) = A(t) X(t) + b(t)]
在这里,(X(t)) 是低阶依赖变量的向量,(A(t)) 是系数矩阵,(b(t)) 是非齐次项,(X'(t)) 是其各自最高阶依赖变量的向量。我们使用“规范形式”一词来暗示上述 ODE 系统。
如果传递的系统具有带有多个解的非线性项,则以其规范形式返回系统列表。
示例
>>> 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)
确定用于与sympy.solvers.ode.systems.linodesolve()
解决的 ODE 系统类型的辅助函数
参数:
A : 矩阵
ODE 系统的系数矩阵
b : 矩阵或 None
系统的非齐次项。默认值为 None。如果该参数为 None,则假定系统为齐次的。
返回:
字典
引发:
NotImplementedError
当系数矩阵没有交换反导数时
解释
此函数接受系数矩阵和/或非齐次项,并返回可以由sympy.solvers.ode.systems.linodesolve()
解决的方程类型。
如果系统是常数系数齐次的,则返回“type1”
如果系统是常数系数非齐次的,则返回“type2”
如果系统是非常数系数的齐次方程组,则返回“type3”
如果系统是非常数系数非齐次的,则返回“type4”
如果系统具有可以因式分解为常数系数矩阵的非常数系数矩阵,则分别针对系统是齐次或非齐次返回“type5”或“type6”。
注意,如果常微分方程组是“type3”或“type4”,则还将返回系数矩阵的可交换反导数,以及类型。
如果系统不能通过 sympy.solvers.ode.systems.linodesolve()
解决,则会引发 NotImplementedError 异常。
例子
>>> 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.
另请参见
linodesolve
linodesolve_type 函数获取信息的函数
sympy.solvers.ode.systems.matrix_exp_jordan_form(A, t)
矩阵指数 (\exp(A*t)) 适用于矩阵 A 和标量 t。
参数:
A : 矩阵
表达式 (\exp(A*t)) 中的矩阵 (A)
t : 符号
自变量
说明
返回 (\exp(A*t)) 的约当形式,以及矩阵 (P),使得:
[\exp(A*t) = P * expJ * P^{-1}]
例子
>>> from sympy import Matrix, Symbol
>>> from sympy.solvers.ode.systems import matrix_exp, matrix_exp_jordan_form
>>> t = Symbol('t')
我们将考虑一个 2x2 不完全矩阵。这表明我们的方法即使对于不完全矩阵也有效。
>>> A = Matrix([[1, 1], [0, 1]])
可以观察到这个函数给出了我们所需的约当标准型以及必要的可逆矩阵 (P)。
>>> P, expJ = matrix_exp_jordan_form(A, t)
在这里,显示了此函数返回的 P 和 expJ 是正确的,因为它们满足公式:P * expJ * P_inverse = exp(A*t)。
>>> P * expJ * P.inv() == matrix_exp(A, t)
True
参考文献
[R885]
zh.wikipedia.org/wiki/Defective_matrix
[R886]
zh.wikipedia.org/wiki/Jordan_matrix
[R887]
zh.wikipedia.org/wiki/Jordan_normal_form
sympy.solvers.ode.systems.matrix_exp(A, t)
矩阵指数 (\exp(A*t)) 适用于矩阵 A
和标量 t
。
参数:
A : 矩阵
表达式 (\exp(A*t)) 中的矩阵 (A)
t : 符号
自变量
说明
此函数通过简单的矩阵乘法返回了 (\exp(A*t)):
[\exp(A*t) = P * expJ * P^{-1}]
其中 (expJ) 是 (\exp(J*t))。 (J) 是矩阵 (A) 的约当标准型,(P) 是使得以下等式成立的矩阵:
[A = P * J * P^{-1}]
矩阵指数 (\exp(A*t)) 出现在线性微分方程的解中。例如,如果 (x) 是向量,(A) 是矩阵,则初始值问题
[\frac{dx(t)}{dt} = A \times x(t), x(0) = x0]
有唯一解
[x(t) = \exp(A t) x0]
例子
>>> from sympy import Symbol, Matrix, pprint
>>> from sympy.solvers.ode.systems import matrix_exp
>>> t = Symbol('t')
我们将考虑一个 2x2 矩阵来计算指数
>>> A = Matrix([[2, -5], [2, -4]])
>>> pprint(A)
[2 -5]
[ ]
[2 -4]
现在,exp(A*t) 给出如下:
>>> 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)]
另请参见
matrix_exp_jordan_form
对于约当标准形的指数
参考文献
[R888]
zh.wikipedia.org/wiki/Jordan_normal_form
[R889]
zh.wikipedia.org/wiki/Matrix_exponential
sympy.solvers.ode.systems.linodesolve(A, t, b=None, B=None, type='auto', doit=False, tau=None)
线性一阶微分方程组的 n 个方程
参数:
A : 矩阵
线性一阶常微分方程组的系数矩阵。
t : 符号
常微分方程组中的自变量。
b:矩阵或 None
ODE 系统中的非齐次项。如果传递了 None,则假定是齐次的 ODE 系统。
B:矩阵或 None
系数矩阵的反导数。如果未传递反导数且解决方案需要该项,则求解器将在内部计算它。
type:字符串
传递的 ODE 系统的类型。根据类型评估解决方案。允许的类型值及其解决的对应系统是:“type1” 用于常系数齐次、“type2” 用于常系数非齐次、“type3” 用于非常系数齐次、“type4” 用于非常系数非齐次,“type5” 和 “type6” 用于系数矩阵可分解为常系数矩阵的非常系数齐次和非齐次系统,分别。“auto” 是默认值,让求解器决定传递系统的正确类型。
doit:布尔值
如果为 True,则评估解决方案,默认值为 False
tau: Expression
用于在得到系统解的后替换 (t) 的值。
返回:
列表
引发:
ValueError
当系数矩阵、非齐次项或传递的反导数不是矩阵或没有正确维度时引发此错误
NonSquareMatrixError
当系数矩阵或其反导数传递的不是方阵时
NotImplementedError
如果系数矩阵没有可交换的反导数
解释
此求解器解决以下形式的 ODE 系统:
[X'(t) = A(t) X(t) + b(t)]
这里,(A(t)) 是系数矩阵,(X(t)) 是 n 个独立变量的向量,(b(t)) 是非齐次项,(X'(t)) 是 (X(t)) 的导数
根据 (A(t)) 和 (b(t)) 的特性,此求解器以不同方式评估解决方案。
当 (A(t)) 是常系数矩阵且 (b(t)) 是零向量即系统是齐次的时候,系统是“type1”。解决方案是:
[X(t) = \exp(A t) C]
这里,(C) 是常数向量,(A) 是常系数矩阵。
当 (A(t)) 是常系数矩阵且 (b(t)) 是非零向量即系统是非齐次的时候,系统是“type2”。解决方案是:
[X(t) = e^{A t} ( \int e^{- A t} b ,dt + C)]
当 (A(t)) 是系数矩阵且其可与其反导数 (B(t)) 交换且 (b(t)) 是零向量即系统是齐次的时候,系统是“type3”。解决方案是:
[X(t) = \exp(B(t)) C]
当 (A(t)) 与其反导数 (B(t)) 可交换且 (b(t)) 是非零向量即系统是非齐次的时候,系统是“type4”。解决方案是:
[X(t) = e^{B(t)} ( \int e^{-B(t)} b(t) ,dt + C)]
当 (A(t)) 是一个系数矩阵,可以分解为标量和常系数矩阵:
[A(t) = f(t) * A]
其中(f(t))是独立变量(t)的标量表达式,(A)是常数矩阵,然后我们可以进行以下替换:
[tau = \int f(t) dt, X(t) = Y(\tau), b(t) = b(f^{-1}(\tau))]
在非齐次项为非零时,仅执行对非齐次项的替换。通过这些替换,我们的原始系统变为:
[Y'(\tau) = A * Y(\tau) + b(\tau)/f(\tau)]
上述系统可以根据系统的齐次性使用“type1”或“type2”的解法轻松解决。在得到(Y(\tau))的解后,我们将(tau)的解作为(t)进行替换,以得到(X(t))。
[X(t) = Y(\tau)]
“type5”和“type6”的系统具有可交换的反导数,但我们使用这个解决方案因为它计算速度更快。
最终的解是所有四个方程的一般解,因为常数系数矩阵始终与其反导数可交换。
此函数的另一个特点是,如果有人想要替换独立变量的值,他们可以传递替换(tau),解决方案将使用传递的表达式((tau))进行独立变量的替换。
示例
要直接使用此函数解决 ODE 系统,必须按正确顺序执行几件事情。对函数的错误输入将导致不正确的结果。
>>> 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))]
在我们推导系数矩阵之前,重要的是将 ODE 系统转换为所需的形式。为此,我们将使用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))]]
现在,我们将使用sympy.solvers.ode.systems.linear_ode_to_matrix()
来获取系数矩阵和非齐次项(如果有的话)。
>>> eqs = eqs[0]
>>> (A1, A0), b = linear_ode_to_matrix(eqs, funcs, x, 1)
>>> A = A0
我们已经准备好系数矩阵和非齐次项。现在,我们可以使用sympy.solvers.ode.systems.linodesolve_type()
来获取 ODE 系统的信息,最终将其传递给求解器。
>>> 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'])
现在,我们可以通过使用sympy.solvers.ode.checkodesol()
来证明解是否正确。
>>> sol = [Eq(f, s) for f, s in zip(funcs, sol_vector)]
>>> checkodesol(eqs, sol)
(True, [0, 0])
我们还可以使用 doit 方法来评估函数传递的解决方案。
>>> sol_vector_evaluated = linodesolve(A, x, b=b, type="type2", doit=True)
现在,我们将看看一个非恒定的 ODE 系统。
>>> eqs = [Eq(f(x).diff(x), f(x) + x*g(x)), Eq(g(x).diff(x), -x*f(x) + g(x))]
上面定义的系统已经处于所需形式,因此我们无需对其进行转换。
>>> (A1, A0), b = linear_ode_to_matrix(eqs, funcs, x, 1)
>>> A = A0
用户还可以传递需要 type3 和 type4 ODEs 系统的可交换反导数。如果传递错误的反导数,将导致不正确的结果。如果系数矩阵与其反导数不可交换,则 sympy.solvers.ode.systems.linodesolve_type()
会引发 NotImplementedError。如果具有可交换的反导数,则函数只返回关于系统的信息。
>>> system_info = linodesolve_type(A, x, b=b)
现在,我们可以将反导数作为参数传递以获得解。如果未传递系统信息,则求解器将在内部计算所需的参数。
>>> sol_vector = linodesolve(A, x, b=b)
再次,我们可以验证所得解。
>>> sol = [Eq(f, s) for f, s in zip(funcs, sol_vector)]
>>> checkodesol(eqs, sol)
(True, [0, 0])
另见
linear_ode_to_matrix
系数矩阵计算函数
canonical_odes
ODEs 系统表示更改
linodesolve_type
获取关于 ODEs 系统的信息以在此求解器中传递
sympy.solvers.ode.ode._nonlinear_2eq_order1_type1(x, y, t, eq)
方程:
[x' = x^n F(x,y)][y' = g(y) F(x,y)]
解:
[x = \varphi(y), \int \frac{1}{g(y) F(\varphi(y),y)} ,dy = t + C_2]
其中
如果 (n \neq 1)
[\varphi = [C_1 + (1-n) \int \frac{1}{g(y)} ,dy]^{\frac{1}{1-n}}]
如果 (n = 1)
[\varphi = C_1 e^{\int \frac{1}{g(y)} ,dy}]
其中 (C_1) 和 (C_2) 是任意常数。
sympy.solvers.ode.ode._nonlinear_2eq_order1_type2(x, y, t, eq)
方程:
[x' = e^{\lambda x} F(x,y)][y' = g(y) F(x,y)]
解:
[x = \varphi(y), \int \frac{1}{g(y) F(\varphi(y),y)} ,dy = t + C_2]
其中
如果 (\lambda \neq 0)
[\varphi = -\frac{1}{\lambda} \log(C_1 - \lambda \int \frac{1}{g(y)} ,dy)]
如果 (\lambda = 0)
[\varphi = C_1 + \int \frac{1}{g(y)} ,dy]
其中 (C_1) 和 (C_2) 是任意常数。
sympy.solvers.ode.ode._nonlinear_2eq_order1_type3(x, y, t, eq)
自主系统的一般形式
[x' = F(x,y)][y' = G(x,y)]
假设 (y = y(x, C_1)),其中 (C_1) 是任意常数是一阶方程的一般解
[F(x,y) y'_x = G(x,y)]
然后,原方程组的一般解形式为
[\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)
方程:
[x' = f_1(x) g_1(y) \phi(x,y,t)][y' = f_2(x) g_2(y) \phi(x,y,t)]
第一个积分:
[\int \frac{f_2(x)}{f_1(x)} ,dx - \int \frac{g_1(y)}{g_2(y)} ,dy = C]
其中 (C) 是任意常数。
在解 (x)(或 (y))的第一个积分和将得到的表达式代入原解的任一方程后,得到用于确定 (y)(或 (x))的一阶方程。
sympy.solvers.ode.ode._nonlinear_2eq_order1_type5(func, t, eq)
Clairaut ODEs 系统
[x = t x' + F(x',y')][y = t y' + G(x',y')]
系统的解如下:
((i)) 直线:
[x = C_1 t + F(C_1, C_2), y = C_2 t + G(C_1, C_2)]
其中 (C_1) 和 (C_2) 是任意常数;
((ii)) 上述线条的包络;
((iii)) 由直线 ((i)) 和 ((ii)) 的线段连续可微而成。
sympy.solvers.ode.ode._nonlinear_3eq_order1_type1(x, y, z, t, eq)
方程:
[a x' = (b - c) y z, \enspace b y' = (c - a) z x, \enspace c z' = (a - b) x y]
第一积分:
[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]
其中 (C_1) 和 (C_2) 是任意常数。通过解 (y) 和 (z) 的积分,并将得到的表达式代入系统的第一方程,我们得到关于 (x) 的分离的一阶方程。类似地,对其他两个方程进行操作,我们也会得到关于 (y) 和 (z) 的一阶方程。
参考资料
-eqworld.ipmnet.ru/en/solutions/sysode/sode0401.pdf
sympy.solvers.ode.ode._nonlinear_3eq_order1_type2(x, y, z, t, eq)
方程:
[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)]
第一积分:
[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]
其中 (C_1) 和 (C_2) 是任意常数。通过解 (y) 和 (z) 的积分,并将得到的表达式代入系统的第一方程,我们得到关于 (x) 的一阶微分方程。同样地,对其他两个方程进行类似操作,我们会得到关于 (y) 和 (z) 的一阶方程。
参考资料:
-eqworld.ipmnet.ru/en/solutions/sysode/sode0402.pdf
sympy.solvers.ode.ode._nonlinear_3eq_order1_type3(x, y, z, t, eq)
方程:
[x' = c F_2 - b F_3, \enspace y' = a F_3 - c F_1, \enspace z' = b F_1 - a F_2]
其中 (F_n = F_n(x, y, z, t))。
- 第一积分:
[a x + b y + c z = C_1,]
其中 (C) 是任意常数。
2. 如果我们假设函数 (F_n) 独立于 (t), 即 (F_n = F_n (x, y, z)),那么在从系统的前两个方程中消除 (t) 和 (z) 后,我们得到关于 (x) 的一阶方程。
[\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)}]
其中 (z = \frac{1}{c} (C_1 - a x - b y))
参考资料:
-eqworld.ipmnet.ru/en/solutions/sysode/sode0404.pdf
sympy.solvers.ode.ode._nonlinear_3eq_order1_type4(x, y, z, t, eq)
方程:
[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]
其中 (F_n = F_n (x, y, z, t))
- 第一积分:
[a x^{2} + b y^{2} + c z^{2} = C_1]
其中 (C) 是任意常数。
2. 假设函数 (F_n) 独立于 (t):(F_n = F_n (x, y, z))。然后,在从系统的前两个方程中消除 (t) 和 (z) 后,我们得到一阶方程
[\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)}]
其中 (z = \pm \sqrt{\frac{1}{c} (C_1 - a x^{2} - b y^{2})})
参考资料:
-eqworld.ipmnet.ru/en/solutions/sysode/sode0405.pdf
sympy.solvers.ode.ode._nonlinear_3eq_order1_type5(x, y, z, t, eq)
[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)]
其中 (F_n = F_n (x, y, z, t)) 且为任意函数。
第一积分:
[\left|x\right|^{a} \left|y\right|^{b} \left|z\right|^{c} = C_1]
其中 (C) 是任意常数。如果函数 (F_n) 不依赖于 (t),则通过从系统的前两个方程中消除 (t) 和 (z),可以得到一个一阶方程。
参考文献
-eqworld.ipmnet.ru/en/solutions/sysode/sode0406.pdf
ode 模块信息
此模块包含 dsolve()
和它使用的不同辅助函数。
dsolve()
解决普通微分方程。查看各个函数的文档字符串以了解它们的用途。注意,偏微分方程的支持在 pde.py
中。注意,提示函数有描述其各种方法的文档字符串,但它们仅供内部使用。使用 dsolve(ode, func, hint=hint)
使用特定提示来解决 ODE。另请参阅 dsolve()
的文档字符串。
此模块中的函数
这些是此模块中的用户函数:
dsolve()
- 解决 ODE。classify_ode()
- 将 ODE 分类为可能的提示供dsolve()
使用。checkodesol()
- 检查方程是否是 ODE 的解。homogeneous_order()
- 返回表达式的齐次阶数。infinitesimals()
- 返回 ODE 的点变换 Lie 群的无穷小,使其保持不变。checkinfsol()
- 检查给定的无穷小是否是一阶 ODE 的实际无穷小。这些是用于内部使用的非解算器辅助函数。用户应使用
dsolve()
的各种选项来获得这些函数提供的功能:
odesimp()
- 进行所有形式的 ODE 简化。ode_sol_simplicity()
- 用于比较解的简易性的关键函数。constantsimp()
- 简化任意常数。constant_renumber()
- 重新编号任意常数。_handle_Integral()
- 评估未评估的积分。参见这些函数的文档字符串。
当前实现的求解器方法
下列方法用于解决普通微分方程。更多信息请查看各种提示函数的文档字符串(运行 help(ode)
):
- 一阶可分离微分方程。
- 系数或 \(dx\) 和 \(dy\) 为同阶函数的一阶微分方程。
- 一阶精确微分方程。
- 一阶线性微分方程。
- 一阶 Bernoulli 微分方程。
- 一阶微分方程的幂级数解法。
- 用 Lie Group 方法解决一阶微分方程。
- 二阶 Liouville 微分方程。
- 二阶微分方程的幂级数解法,在普通和正则奇点处。
- 可以通过代数重排和积分解决的 \(n\) 阶微分方程。
- 具有常系数的 \(n\) 阶线性齐次微分方程。
- 具有常系数的 \(n\) 阶线性非齐次微分方程,使用特解系数方法。
- 具有常系数的 \(n\) 阶线性非齐次微分方程,使用参数变化法。
该模块的理念
该模块设计用于简化添加新的 ODE 求解方法,无需干扰其他方法的求解代码。其思想在于有一个 classify_ode()
函数,输入一个 ODE 并告诉您什么提示(如果有)可以解决该 ODE。它在不尝试解决 ODE 的情况下执行此操作,因此速度很快。每个求解方法都是一个提示,并且有自己的函数,命名为 ode_<hint>
。该函数接受 ODE 和由 classify_ode()
收集的任何匹配表达式,并返回解决的结果。如果此结果中有任何积分,提示函数将返回一个未评估的 Integral
类。dsolve()
是用户包装器函数,将对结果调用 odesimp()
,它将尝试解决依赖变量(我们正在解决的函数)的方程,简化表达式中的任意常数,并评估任何积分,如果提示允许的话。
如何添加新的解法方法
如果你有一个希望dsolve()
能够解决的 ODE,请尽量避免在这里添加特殊情况的代码。相反,尝试找到一种通用的方法来解决你的 ODE,以及其他的 ODE。这样,ode
模块将变得更加健壮,而且不受特殊情况的限制。WolphramAlpha 和 Maple 的 DETools[odeadvisor] 函数是你可以用来分类特定 ODE 的两个资源。如果可能的话,一个方法最好是能处理一个(n)阶的 ODE,而不仅仅是特定阶数的 ODE。
要添加一个新的方法,你需要做几件事情。首先,你需要为你的方法取一个提示名字。尽量命名你的提示,以便与所有其他方法都不会产生歧义,包括可能尚未实现的方法。如果你的方法使用积分,还要包括一个 hint_Integral
提示。如果有多种方法可以解决 ODE,那么每种方法都应包括一个提示,以及一个 <hint>_best
提示。你的 ode_<hint>_best()
函数应该使用 ode_sol_simplicity
作为关键参数选择最佳方法。例如参见HomogeneousCoeffBest
。使用你的方法的函数将被称为 ode_<hint>()
,因此提示必须仅使用 Python 函数名允许的字符(字母数字字符和下划线‘_
’字符)。除了 _Integral
提示(dsolve()
会自动处理这些),每个提示都应包括一个函数。提示名应全部小写,除非是一个常见的大写单词(如 Integral 或 Bernoulli)。如果你有一个不希望与 all_Integral
一起运行的提示,没有 _Integral
对应项(例如会破坏 all_Integral
目的的最佳提示),你需要在 dsolve()
代码中手动移除它。另请参阅 classify_ode()
的文档字符串,了解编写提示名称的指导方针。
总的来说,要确定您的方法返回的解决方案与其他可能解决相同 ODE 的方法相比如何。然后,将您的提示放入allhints
元组中,以确定调用它们的顺序。此元组的排序决定了默认的提示顺序。请注意,异常是可以接受的,因为用户很容易使用dsolve()
选择单个提示。一般来说,_Integral
变体应放在列表的末尾,而_best
变体应放在其适用的各种提示之前。例如,undetermined_coefficients
提示在variation_of_parameters
提示之前,因为虽然参数变化比不定系数更一般,但不定系数通常为其能解决的 ODE 生成更干净的结果,并且不需要积分,因此速度更快。
接下来,你需要有一个匹配表达式或函数,匹配 ODE 的类型,你应该将其放在classify_ode()
中(如果匹配函数不止几行。应尽可能匹配 ODE 而不解决它,以便classify_ode()
保持快速,并且不受解决代码中错误的影响。务必考虑边界情况。例如,如果你的解决方法涉及除以某些东西,请确保排除除数为 0 的情况。
在大多数情况下,ODE 的匹配还将为您提供解决它所需的各个部分。您应该将其放在一个字典中(.match()
将为您执行此操作),并在classify_ode()
的相关部分中添加matching_hints['hint'] = matchdict
。classify_ode()
然后会将此发送到dsolve()
,后者将其作为match
参数发送到您的函数。您的函数应命名为pyode_<hint>(eq, func, order, match)`. If you need to send more information, put it in the ``match
字典。例如,如果您必须在classify_ode()
中替换一个虚拟变量以匹配 ODE,则需要使用(match)字典将其传递给您的函数。您可以使用func.args[0]
访问自变量,使用func.func
访问因变量(您试图解决的函数)。如果在试图解决 ODE 时发现无法解决,则引发NotImplementedError
。dsolve()
将使用all
元提示捕获此错误,而不会导致整个例程失败。
为你的函数添加一个描述所采用方法的文档字符串。与 SymPy 中的其他任何内容一样,你需要在文档字符串中添加一个 doctest,并在 test_ode.py
中添加真实的测试。尽量保持与其他提示函数的文档字符串一致性。将你的方法添加到此文档字符串顶部的列表中。同时,将你的方法添加到 docs/src
目录中的 ode.rst
中,以便 Sphinx 文档将其文档字符串引入主 SymPy 文档中。确保通过在文档目录中运行 make html
来生成 Sphinx 文档,以验证文档字符串的格式是否正确。
如果你的解决方法涉及积分,请使用Integral
而不是integrate()
。这使用户可以通过使用你提示的 _Integral
变体来避开复杂/慢速的积分。在大多数情况下,调用sympy.core.basic.Basic.doit()
将对你的解决方案进行积分。如果不是这种情况,你需要在_handle_Integral()
中编写特殊代码。任意常数应命名为C1
、C2
等符号。所有解决方法都应返回一个等式实例。如果需要任意数量的任意常数,可以使用constants = numbered_symbols(prefix='C', cls=Symbol, start=1)
。如果可能以一般方式解决依赖函数,请这样做。否则,请尽力而为,但不要在你的 ode_<hint>()
函数中调用 solve
。odesimp()
将尝试为你解决方案解决问题,因此你不需要这样做。最后,如果你的 ODE 有一个可以应用于解决方案的常见简化,可以在odesimp()
中添加一个特殊情况。例如,从 1st_homogeneous_coeff
提示返回的解经常有许多log
项,因此 odesimp()
在它们上调用 logcombine()
(在这种情况下,将任意常数写为 log(C1)
而不是 C1
也很有帮助)。还要考虑常见的解决方案重新排列的方法,以便 constantsimp()
更好地利用它。最好将简化放在 odesimp()
中而不是在你的方法中,因为可以在 dsolve()
的 simplify
标志中关闭它。如果在你的函数中有任何多余的简化,请确保只在 if match.get('simplify', True):
下运行它,特别是如果它可能很慢或者可能减少解决方案的域。
最后,与对 SymPy 的每一次贡献一样,您的方法将需要进行测试。在test_ode.py
中为每种方法添加一个测试。遵循那里的惯例,即使用dsolve(eq, f(x), hint=your_hint)
测试求解器,并使用checkodesol()
测试解决方案(如果运行过慢或者不起作用,可以将这些放在单独的测试中并跳过/XFAIL)。确保在dsolve()
中特别调用您的提示,这样测试不会仅仅因为引入另一个匹配提示而失败。如果您的方法适用于高阶(>1)ODEs,您将需要对每个解运行sol = constant_renumber(sol, 'C', 1, order)
,其中order
是 ODE 的阶数。这是因为constant_renumber
通过打印顺序重新编号任意常数,这取决于平台。尽量测试您求解器的每个边界情况,包括一系列阶数,如果它是一个\(n\)阶求解器的话,但如果您的求解器运行速度较慢,比如涉及到复杂的积分,尽量保持测试运行时间短。
请随意重构现有提示,以避免重复代码或创建不一致。如果您能证明您的方法完全复制了现有方法,包括获取解决方案的简易性和速度,则可以删除旧的、不太普遍的方法。现有代码在test_ode.py
中经过了广泛测试,因此如果有任何问题,其中的某个测试肯定会失败。
这些函数不适用于最终用户使用。
sympy.solvers.ode.ode._handle_Integral(expr, func, hint)
将具有积分的解转换为实际解决方案。
对于大多数提示,这只是运行expr.doit()
。
偏微分方程
用户函数
这些是导入全局命名空间的函数 from sympy import *
。它们供用户使用。
sympy.solvers.pde.pde_separate(eq, fun, sep, strategy='mul')
通过加法或乘法分离方法分离偏微分方程中的变量。它尝试重写方程,使指定的变量之一出现在与其他变量不同的方程一侧。
参数:
-
eq – 偏微分方程
-
fun – 原始函数 F(x, y, z)
-
分离 – 分离函数列表 [X(x), u(y, z)]
-
策略 – 分离策略。您可以选择加法分离(‘add’)和乘法分离(‘mul’),默认为乘法分离。
示例
>>> from sympy import E, Eq, Function, pde_separate, Derivative as D
>>> from sympy.abc import x, t
>>> u, X, T = map(Function, 'uXT')
>>> eq = Eq(D(u(x, t), x), E**(u(x, t))*D(u(x, t), t))
>>> pde_separate(eq, u(x, t), [X(x), T(t)], strategy='add')
[exp(-X(x))*Derivative(X(x), x), exp(T(t))*Derivative(T(t), t)]
>>> eq = Eq(D(u(x, t), x, 2), D(u(x, t), t, 2))
>>> pde_separate(eq, u(x, t), [X(x), T(t)], strategy='mul')
[Derivative(X(x), (x, 2))/X(x), Derivative(T(t), (t, 2))/T(t)]
另请参阅
pde_separate_add
, pde_separate_mul
sympy.solvers.pde.pde_separate_add(eq, fun, sep)
用于搜索加法可分离解的辅助函数。
考虑具有两个独立变量 x、y 和一个依赖变量 w 的方程式,我们寻找依赖于不同参数的两个函数的乘积:
(w(x, y, z) = X(x) + y(y, z))
示例
>>> from sympy import E, Eq, Function, pde_separate_add, Derivative as D
>>> from sympy.abc import x, t
>>> u, X, T = map(Function, 'uXT')
>>> eq = Eq(D(u(x, t), x), E**(u(x, t))*D(u(x, t), t))
>>> pde_separate_add(eq, u(x, t), [X(x), T(t)])
[exp(-X(x))*Derivative(X(x), x), exp(T(t))*Derivative(T(t), t)]
sympy.solvers.pde.pde_separate_mul(eq, fun, sep)
用于搜索乘法可分离解的辅助函数。
考虑具有两个独立变量 x、y 和一个依赖变量 w 的方程式,我们寻找依赖于不同参数的两个函数的乘积:
(w(x, y, z) = X(x)*u(y, z))
示例
>>> from sympy import Function, Eq, pde_separate_mul, Derivative as D
>>> from sympy.abc import x, y
>>> u, X, Y = map(Function, 'uXY')
>>> eq = Eq(D(u(x, y), x, 2), D(u(x, y), y, 2))
>>> pde_separate_mul(eq, u(x, y), [X(x), Y(y)])
[Derivative(X(x), (x, 2))/X(x), Derivative(Y(y), (y, 2))/Y(y)]
sympy.solvers.pde.pdsolve(eq, func=None, hint='default', dict=False, solvefun=None, **kwargs)
解决任何(支持的)类型的偏微分方程。
用法
`pdsolve(eq, f(x,y), hint) -> 解决偏微分方程 eq,得到函数 f(x,y),使用方法 hint。
详细信息
eq
可以是任何支持的偏微分方程(请参见支持的方法的偏微分方程文档字符串)。这可以是一个等式,也可以是一个假设等于 0 的表达式。
f(x,y)
是一个具有该变量导数的两个变量函数。变量构成偏微分方程。在许多情况下,无需提供它;它将被自动检测到(如果无法检测到,则会引发错误)。
hint
是您希望 pdsolve 使用的解决方法。使用使用 classify_pde(eq, f(x,y)) 来获取偏微分方程的所有可能提示。默认提示‘default’将使用 classify_pde() 返回的第一个提示。有关可以用作提示的更多选项,请参见下面的提示。
solvefun
是返回任意函数的约定。由 PDE 求解器。如果用户未设置,它默认为 F。
提示
除了各种解决方法外,还有一些元提示,您可以传递给 pdsolve():
“默认”:
这使用由 classify_pde() 返回的第一个提示。这是 pdsolve() 的默认参数。
“全部”:
要使 pdsolve 应用所有相关的分类提示,请使用 pdsolve(PDE, func, hint=”all”)。这将返回一个提示:解决方案字典条目。如果提示导致 pdsolve 抛出 NotImplementedError,则该提示键的值将是引发的异常对象。该字典还将包括一些特殊键:
- order: PDE 的阶数。另请参阅 deutils.py 中的 ode_order()。
- 默认:默认情况下将返回的解决方案。这是由 classify_pde() 返回的元组中首次出现的提示产生的解决方案。
“all_Integral”:
这与“all”相同,只是如果提示还有相应的“_Integral”提示,则仅返回“_Integral”提示。如果“all”由于难以或不可能的积分而导致 pdsolve() 挂起,这将非常有用。这个元提示也比“all”快得多,因为 integrate() 是一个昂贵的例程。
另请参阅 classify_pde() 的文档字符串以获取有关提示的更多信息,以及 pde 的文档字符串以获取所有支持的提示列表。
提示
-
您可以这样声明未知函数的导数:
>>> from sympy import Function, Derivative >>> from sympy.abc import x, y # x and y are the independent variables >>> f = Function("f")(x, y) # f is a function of x and y >>> # fx will be the partial derivative of f with respect to x >>> fx = Derivative(f, x) >>> # fy will be the partial derivative of f with respect to y >>> fy = Derivative(f, y)
-
请查看 test_pde.py 进行许多测试,这也可以作为如何使用 pdsolve() 的示例集。
-
pdsolve 总是返回一个 Equality 类(除非提示为“all”或“all_Integral”)。请注意,无法像 ODE 的情况那样获得 f(x, y) 的显式解。
-
执行 help(pde.pde_hintname) 以获取有关特定提示的更多信息
示例
>>> from sympy.solvers.pde import pdsolve
>>> from sympy import Function, Eq
>>> from sympy.abc import x, y
>>> f = Function('f')
>>> u = f(x, y)
>>> ux = u.diff(x)
>>> uy = u.diff(y)
>>> eq = Eq(1 + (2*(ux/u)) + (3*(uy/u)), 0)
>>> pdsolve(eq)
Eq(f(x, y), F(3*x - 2*y)*exp(-2*x/13 - 3*y/13))
sympy.solvers.pde.classify_pde(eq, func=None, dict=False, *, prep=True, **kwargs)
返回一个 PDE 的可能 pdsolve() 分类的元组。
元组被排序,以便第一项是 pdsolve() 默认用于解决 PDE 的分类。一般来说,列表开头附近的分类比列表末尾附近的分类更快地产生更好的解决方案,尽管总会有例外。要使 pdsolve 使用不同的分类,请使用 pdsolve(PDE, func, hint=
如果 dict
为真,则 classify_pde() 将返回一个提示:匹配表达式字典条目。这是为 pdsolve() 的内部使用而设计的。请注意,由于字典的顺序是任意的,因此这很可能不会与元组的顺序相同。
您可以通过执行 help(pde.pde_hintname) 来获取有关不同提示的帮助,其中 hintname 是不带“_Integral”的提示名称。
请查看 sympy.pde.allhints 或 sympy.pde 文档字符串,以获取从 classify_pde 返回的所有支持的提示列表。
示例
>>> from sympy.solvers.pde import classify_pde
>>> from sympy import Function, Eq
>>> from sympy.abc import x, y
>>> f = Function('f')
>>> u = f(x, y)
>>> ux = u.diff(x)
>>> uy = u.diff(y)
>>> eq = Eq(1 + (2*(ux/u)) + (3*(uy/u)), 0)
>>> classify_pde(eq)
('1st_linear_constant_coeff_homogeneous',)
sympy.solvers.pde.checkpdesol(pde, sol, func=None, solve_for_func=True)
检查给定解是否满足偏微分方程。
pde 是可以以方程或表达式形式给出的偏微分方程。sol 是要检查 pde 是否满足的解。这也可以以方程或表达式形式给出。如果未提供函数,则将使用 deutils 中的辅助函数 _preprocess 来识别函数。
如果传递了一系列解决方案,则将使用相同类型的容器来为每个解决方案返回结果。
目前正在实施以下方法来检查解是否满足 PDE:
- 直接将解代入偏微分方程并检查。如果尚未解出(f),则会解出(f),前提是未将 solve_for_func 设置为 False。
如果解满足 PDE,则返回元组(True, 0)。否则返回元组(False, expr),其中 expr 是将解代入 PDE 后获得的值。但如果已知解返回 False,则可能是因为 doit()无法将其简化为零。
举例
>>> from sympy import Function, symbols
>>> from sympy.solvers.pde import checkpdesol, pdsolve
>>> x, y = symbols('x y')
>>> f = Function('f')
>>> eq = 2*f(x,y) + 3*f(x,y).diff(x) + 4*f(x,y).diff(y)
>>> sol = pdsolve(eq)
>>> assert checkpdesol(eq, sol)[0]
>>> eq = x*f(x,y) + f(x,y).diff(x)
>>> checkpdesol(eq, sol)
(False, (x*F(4*x - 3*y) - 6*F(4*x - 3*y)/25 + 4*Subs(Derivative(F(_xi_1), _xi_1), _xi_1, 4*x - 3*y))*exp(-6*x/25 - 8*y/25))
提示方法
这些函数用于内部使用。但它们包含有关各种求解方法的有用信息。
sympy.solvers.pde.pde_1st_linear_constant_coeff_homogeneous(eq, func, order, match, solvefun)
解决具有恒定系数的一阶线性齐次偏微分方程。
这个偏微分方程的一般形式是
[a \frac{\partial f(x,y)}{\partial x} + b \frac{\partial f(x,y)}{\partial y} + c f(x,y) = 0]
其中(a)、(b)和(c)是常数。
一般解的形式为:
[f(x, y) = F(- a y + b x ) e^{- \frac{c (a x + b y)}{a² + b²}}]
并且可以在 SymPy 中使用pdsolve
找到:
>>> from sympy.solvers import pdsolve
>>> from sympy.abc import x, y, a, b, c
>>> from sympy import Function, pprint
>>> f = Function('f')
>>> u = f(x,y)
>>> ux = u.diff(x)
>>> uy = u.diff(y)
>>> genform = a*ux + b*uy + c*u
>>> pprint(genform)
d d
a*--(f(x, y)) + b*--(f(x, y)) + c*f(x, y)
dx dy
>>> pprint(pdsolve(genform))
-c*(a*x + b*y)
---------------
2 2
a + b
f(x, y) = F(-a*y + b*x)*e
举例
>>> from sympy import pdsolve
>>> from sympy import Function, pprint
>>> from sympy.abc import x,y
>>> f = Function('f')
>>> pdsolve(f(x,y) + f(x,y).diff(x) + f(x,y).diff(y))
Eq(f(x, y), F(x - y)*exp(-x/2 - y/2))
>>> pprint(pdsolve(f(x,y) + f(x,y).diff(x) + f(x,y).diff(y)))
x y
- - - -
2 2
f(x, y) = F(x - y)*e
参考文献
- Viktor Grigoryan, “Partial Differential Equations” Math 124A - Fall 2010, pp.7
sympy.solvers.pde.pde_1st_linear_constant_coeff(eq, func, order, match, solvefun)
解决具有恒定系数的一阶线性偏微分方程。
这个偏微分方程的一般形式是
[a \frac{\partial f(x,y)}{\partial x} + b \frac{\partial f(x,y)}{\partial y} + c f(x,y) = G(x,y)]
其中(a)、(b)和(c)是常数,而(G(x, y))可以是(x)和(y)的任意函数。
偏微分方程的一般解为:
[\begin{split}f(x, y) = \left. \left[F(\eta) + \frac{1}{a² + b²} \int\limits^{a x + b y} G\left(\frac{a \xi + b \eta}{a² + b²}, \frac{- a \eta + b \xi}{a² + b²} \right) e^{\frac{c \xi}{a² + b²}}, d\xi\right] e^{- \frac{c \xi}{a² + b²}} \right|_{\substack{\eta=- a y + b x\ \xi=a x + b y }}, ,\end{split}]
其中(F(\eta))是任意单值函数。可以在 SymPy 中使用pdsolve
找到解:
>>> from sympy.solvers import pdsolve
>>> from sympy.abc import x, y, a, b, c
>>> from sympy import Function, pprint
>>> f = Function('f')
>>> G = Function('G')
>>> u = f(x, y)
>>> ux = u.diff(x)
>>> uy = u.diff(y)
>>> genform = a*ux + b*uy + c*u - G(x,y)
>>> pprint(genform)
d d
a*--(f(x, y)) + b*--(f(x, y)) + c*f(x, y) - G(x, y)
dx dy
>>> pprint(pdsolve(genform, hint='1st_linear_constant_coeff_Integral'))
// a*x + b*y \ \|
|| / | ||
|| | | ||
|| | c*xi | ||
|| | ------- | ||
|| | 2 2 | ||
|| | /a*xi + b*eta -a*eta + b*xi\ a + b | ||
|| | G|------------, -------------|*e d(xi)| ||
|| | | 2 2 2 2 | | ||
|| | \ a + b a + b / | -c*xi ||
|| | | -------||
|| / | 2 2||
|| | a + b ||
f(x, y) = ||F(eta) + -------------------------------------------------------|*e ||
|| 2 2 | ||
\\ a + b / /|eta=-a*y + b*x, xi=a*x + b*y
举例
>>> from sympy.solvers.pde import pdsolve
>>> from sympy import Function, pprint, exp
>>> from sympy.abc import x,y
>>> f = Function('f')
>>> eq = -2*f(x,y).diff(x) + 4*f(x,y).diff(y) + 5*f(x,y) - exp(x + 3*y)
>>> pdsolve(eq)
Eq(f(x, y), (F(4*x + 2*y)*exp(x/2) + exp(x + 4*y)/15)*exp(-y))
参考文献
- Viktor Grigoryan, “Partial Differential Equations” Math 124A - Fall 2010, pp.7
sympy.solvers.pde.pde_1st_linear_variable_coeff(eq, func, order, match, solvefun)
解决具有可变系数的一阶线性偏微分方程。这个偏微分方程的一般形式是
[a(x, y) \frac{\partial f(x, y)}{\partial x} + b(x, y) \frac{\partial f(x, y)}{\partial y} + c(x, y) f(x, y) = G(x, y)]
其中(a(x, y))、(b(x, y))、(c(x, y))和(G(x, y))是(x)和(y)的任意函数。通过以下变换将这个偏微分方程转换为 ODE:
-
(\xi)作为(x)
-
(\eta)作为解的常数,满足偏微分方程(\frac{dy}{dx} = -\frac{b}{a})
进行前述替换后,将其简化为线性 ODE
[a(\xi, \eta)\frac{du}{d\xi} + c(\xi, \eta)u - G(\xi, \eta) = 0]
可以使用dsolve
解决。
>>> from sympy.abc import x, y
>>> from sympy import Function, pprint
>>> a, b, c, G, f= [Function(i) for i in ['a', 'b', 'c', 'G', 'f']]
>>> u = f(x,y)
>>> ux = u.diff(x)
>>> uy = u.diff(y)
>>> genform = a(x, y)*u + b(x, y)*ux + c(x, y)*uy - G(x,y)
>>> pprint(genform)
d d
-G(x, y) + a(x, y)*f(x, y) + b(x, y)*--(f(x, y)) + c(x, y)*--(f(x, y))
dx dy
举例
>>> from sympy.solvers.pde import pdsolve
>>> from sympy import Function, pprint
>>> from sympy.abc import x,y
>>> f = Function('f')
>>> eq = x*(u.diff(x)) - y*(u.diff(y)) + y**2*u - y**2
>>> pdsolve(eq)
Eq(f(x, y), F(x*y)*exp(y**2/2) + 1)
参考文献
- Viktor Grigoryan, “Partial Differential Equations” Math 124A - Fall 2010, pp.7
pde 模块信息
此模块包含 pdsolve()及其使用的不同辅助函数。它受 ode 模块的启发,因此基本架构保持不变。
此模块中的函数
这些是此模块中的用户函数:
- pdsolve() - 解决 PDE’s
- classify_pde() - 将偏微分方程分类为 dsolve()可能的提示。
- pde_separate() - 通过变量分离偏微分方程。
- 加法或乘法分离方法。
这些是此模块中的辅助函数:
- pde_separate_add() - 用于搜索加法可分解解的辅助函数。
- pde_separate_mul() - 用于搜索乘法可分解解的辅助函数。
- 可分离解。
当前实现的解算器方法
实现了以下方法以解决偏微分方程。有关每个函数的详细信息,请参阅各种 pde_hint()函数的文档字符串(运行 help(pde)):
- 常系数的一阶线性齐次偏微分方程。
- 常系数的一阶线性常规偏微分方程。
- 变系数的一阶线性偏微分方程。
解算器
SymPy 中的solvers模块实现了解方程的方法。
注意
对于专注于解决常见类型方程的初学者指南,请参阅 Solve Equations。
注意
solve()
是一个更老、更成熟的通用函数,用于解决许多类型的方程。 solve()
具有许多选项,并且在内部使用不同的方法来确定您传递的方程类型,因此如果您知道正在处理的方程类型,可能希望使用更新的solveset()
来解决一元方程,linsolve()
来解决线性方程组,以及nonlinsolve()
来解决非线性方程组。
代数方程
使用solve()
来解代数方程。我们假设所有方程都等于 0,因此解决 x**2 == 1 可以转换为以下代码:
>>> from sympy.solvers import solve
>>> from sympy import Symbol
>>> x = Symbol('x')
>>> solve(x**2 - 1, x)
[-1, 1]
solve()
的第一个参数是一个等于零的方程,第二个参数是我们想要解的符号。
sympy.solvers.solvers.solve(f, *symbols, **flags)
代数地解方程和方程组。
参数:
f:
- 一个单一的 Expr 或 Poly 必须为零
- 一个等式
- 一个关系表达式
- 一个布尔值
- 以上的一个或多个可迭代对象
符号:(要解决的对象)指定为
- 未指定(将使用其他非数字对象)
- 单一符号
- 符号的去嵌套列表(例如,
solve(f, x, y)
)- 符号的有序迭代对象(例如,
solve(f, [x, y])
)
标志:
dict=True(默认为 False)
返回映射解列表(可能为空)。
set=True(默认为 False)
返回符号列表和解的元组集合。
exclude=[](默认)
不要尝试解决排除列表中的任何自由符号;如果给出表达式,其中的自由符号将自动提取。
check=True(默认)
如果为 False,则不进行任何解的测试。如果要包含使任何分母为零的解,这可能很有用。
numerical=True(默认)
如果f仅有一个符号,进行快速数值检查。
minimal=True(默认为 False)
一种非常快速、最小的测试。
warn=True(默认为 False)
如果
checksol()
无法得出结论,则显示警告。simplify=True(默认)
在返回它们之前简化所有三阶及以上的多项式,并且(如果检查不为 False)在将它们代入应为零的函数时使用通用简化函数对解决方案和表达式进行简化。
force=True(默认为 False)
使所有没有关于符号的假设的符号为正。
rational=True(默认)
重置浮点数为有理数;如果不使用此选项,包含浮点数的系统可能因多项式问题而无法解决。如果 rational=None,则浮点数将被重置为有理数,但答案将被重置为浮点数。如果标志为 False,则不会对浮点数进行任何处理。
manual=True(默认为 False)
不要使用多项式/矩阵方法来解决方程组,按顺序解决它们,就像“手动”一样。
implicit=True(默认为 False)
允许
solve
返回模式的解,以其他包含该模式的函数来表示;仅在该模式位于某些可逆函数内部(如 cos,exp 等)时才需要。particular=True(默认为 False)
指示
solve
尝试找到线性系统的特定解,其中尽可能多的解为零;这是非常昂贵的。quick=True(默认为 False;
particular
必须为 True)选择一个快速的启发式方法来找到具有许多零的解,而值为 False 则使用保证找到尽可能多零的非常缓慢的方法。
cubics=True(默认)
遇到立方表达式时返回显式解。当为 False 时,四次方和五次方也被禁用。
quartics=True(默认)
遇到四次表达式时返回显式解。当为 False 时,五次方也被禁用。
quintics=True(默认)
遇到五次表达式时返回显式解(如果可能)。
说明
当前支持:
-
多项式
-
无穷的
-
上述内容的分段组合
-
线性和多项式方程组
-
包含关系表达式的系统
-
由未定系数暗示的系统
示例
默认输出根据输入而变化,可能是一个列表(可能为空)、一个字典、一个列表字典或元组,或者是涉及关系的表达式。具体关于可能出现的不同形式的输出,请参见 Solve Output by Type。在这里,可以说为了从solve
获取统一的输出,使用dict=True
或set=True
(见下文)足够了。
>>> from sympy import solve, Poly, Eq, Matrix, Symbol
>>> from sympy.abc import x, y, z, a, b
传递的表达式可以是 Expr、Equality 或 Poly 类(或相同类型的列表);矩阵被视为包含矩阵所有元素的列表:
>>> solve(x - 3, x)
[3]
>>> solve(Eq(x, 3), x)
[3]
>>> solve(Poly(x - 3), x)
[3]
>>> solve(Matrix([[x, x + y]]), x, y) == solve([x, x + y], x, y)
True
如果没有指示感兴趣的符号,并且方程式是单变量的,则返回值列表;否则,字典中的键将指示在表达式中找到的哪些(所有变量中使用的变量和解决方案)变量:
>>> solve(x**2 - 4)
[-2, 2]
>>> solve((x - a)*(y - b))
[{a: x}, {b: y}]
>>> solve([x - 3, y - 1])
{x: 3, y: 1}
>>> solve([x - 3, y**2 - 1])
[{x: 3, y: -1}, {x: 3, y: 1}]
如果你传递寻找解的符号,输出将根据你传递的符号数量、是否传递表达式列表以及是否解决了线性系统而变化。通过使用dict=True
或set=True
来获得统一的输出。
>>> #### *** feel free to skip to the stars below *** #### >>> from sympy import TableForm >>> h = [None, ';|;'.join(['e', 's', 'solve(e, s)', 'solve(e, s, dict=True)', ... 'solve(e, s, set=True)']).split(';')] >>> t = [] >>> for e, s in [ ... (x - y, y), ... (x - y, [x, y]), ... (x**2 - y, [x, y]), ... ([x - 3, y -1], [x, y]), ... ]: ... how = [{}, dict(dict=True), dict(set=True)] ... res = [solve(e, s, **f) for f in how] ... t.append([e, '|', s, '|'] + [res[0], '|', res[1], '|', res[2]]) ... >>> # ******************************************************* # >>> TableForm(t, headings=h, alignments="<") e | s | solve(e, s) | solve(e, s, dict=True) | solve(e, s, set=True) --------------------------------------------------------------------------------------- x - y | y | [x] | [{y: x}] | ([y], {(x,)}) x - y | [x, y] | [(y, y)] | [{x: y}] | ([x, y], {(y, y)}) x**2 - y | [x, y] | [(x, x**2)] | [{y: x**2}] | ([x, y], {(x, x**2)}) [x - 3, y - 1] | [x, y] | {x: 3, y: 1} | [{x: 3, y: 1}] | ([x, y], {(3, 1)})
- 如果任何方程不依赖于给定的符号,则将从方程组中排除它,并且可能以对不感兴趣的变量隐式给出答案:
>>> solve([x - y, y - 3], x) {x: y}
当你传递除一个自由符号外的所有自由符号时,将尝试基于不定系数法找到单个解。如果成功,将返回值的字典。如果要对一个或多个符号的表达式进行代数解,请将其传递给要解决的列表:
>>> e = a*x + b - 2*x - 3
>>> solve(e, [a, b])
{a: 2, b: 3}
>>> solve([e], [a, b])
{a: -b/x + (2*x + 3)/x}
当没有任何符号的解使得所有表达式都为零时,将返回空列表(或在set=True
时返回空集):
>>> from sympy import sqrt
>>> solve(3, x)
[]
>>> solve(x - 3, y)
[]
>>> solve(sqrt(x) + 1, x, set=True)
([x], set())
如果作为符号给出的对象不是符号,则会进行代数隔离,并可能获得隐式解。这主要是为了方便你节省将对象替换为符号并解决该符号的时间。它只在指定对象可以使用subs
方法替换为符号时才起作用:
>>> from sympy import exp, Function >>> f = Function('f')
>>> solve(f(x) - x, f(x)) [x] >>> solve(f(x).diff(x) - f(x) - x, f(x).diff(x)) [x + f(x)] >>> solve(f(x).diff(x) - f(x) - x, f(x)) [-x + Derivative(f(x), x)] >>> solve(x + exp(x)**2, exp(x), set=True) ([exp(x)], {(-sqrt(-x),), (sqrt(-x),)})
>>> from sympy import Indexed, IndexedBase, Tuple >>> A = IndexedBase('A') >>> eqs = Tuple(A[1] + A[2] - 3, A[1] - A[2] + 1) >>> solve(eqs, eqs.atoms(Indexed)) {A[1]: 1, A[2]: 2}
- 要解决导数中的函数,请使用
dsolve()
。
要隐式解符号,请使用 implicit=True:
>>> solve(x + exp(x), x)
[-LambertW(1)]
>>> solve(x + exp(x), x, implicit=True)
[-exp(x)]
可以解决表达式中任何可以使用subs
替换为符号的内容。
>>> solve(x + 2 + sqrt(3), x + 2) [-sqrt(3)] >>> solve((x + 2 + sqrt(3), x + 4 + y), y, x + 2) {y: -2 + sqrt(3), x + 2: -sqrt(3)}
在这种隐式求解中没有采取任何英雄主义行为,因此可能最终解中仍然会出现符号:
>>> eqs = (x*y + 3*y + sqrt(3), x + 4 + y) >>> solve(eqs, y, x + 2) {y: -sqrt(3)/(x + 3), x + 2: -2*x/(x + 3) - 6/(x + 3) + sqrt(3)/(x + 3)} >>> solve(eqs, y*x, x) {x: -y - 4, x*y: -3*y - sqrt(3)}
如果尝试解决一个数字,请记住,你获得的数字并不一定意味着该值等同于所获得的表达式:
>>> solve(sqrt(2) - 1, 1) [sqrt(2)] >>> solve(x - y + 1, 1) # /!\ -1 is targeted, too [x/(y - 1)] >>> [_.subs(z, -1) for _ in solve((x - y + 1).subs(-1, z), 1)] [-x + y]
附加示例
使用check=True
(默认情况下)的solve()
将通过符号标签来消除不需要的解。如果不包括任何假设,则将返回所有可能的解:
>>> x = Symbol("x")
>>> solve(x**2 - 1)
[-1, 1]
通过设置positive
标志,将只返回一个解:
>>> pos = Symbol("pos", positive=True)
>>> solve(pos**2 - 1)
[1]
当检查解时,使得任何分母为零的解将自动被排除。如果不想排除这样的解,请使用check=False
选项:
>>> from sympy import sin, limit
>>> solve(sin(x)/x) # 0 is excluded
[pi]
如果check=False
,则会找到使得分子为零的解,但是当(x = 0)时,(\sin(x)/x)具有众所周知的极限(无间断),为 1:
>>> solve(sin(x)/x, check=False)
[0, pi]
在下面的情况中,限制存在并等于在check=True
时被排除的(x = 0)的值:
>>> eq = x**2*(1/x - z**2/x)
>>> solve(eq, x)
[]
>>> solve(eq, x, check=False)
[0]
>>> limit(eq, x, 0, '-')
0
>>> limit(eq, x, 0, '+')
0
解决关系
当传递给solve
的一个或多个表达式为关系式时,将返回一个关系结果(并且dict
和set
标志将被忽略):
>>> solve(x < 3)
(-oo < x) & (x < 3)
>>> solve([x < 3, x**2 > 4], x)
((-oo < x) & (x < -2)) | ((2 < x) & (x < 3))
>>> solve([x + y - 3, x > 3], x)
(3 < x) & (x < oo) & Eq(x, 3 - y)
虽然未对关系中的符号进行检查假设,但设置假设将影响某些关系可能自动简化的方式:
>>> solve(x**2 > 4)
((-oo < x) & (x < -2)) | ((2 < x) & (x < oo))
>>> r = Symbol('r', real=True)
>>> solve(r**2 > 4)
(2 < r) | (r < -2)
目前在 SymPy 中没有算法允许您使用关系解析多个变量。因此,以下内容不能确定q < 0
(尝试解析r
和q
会引发错误):
>>> from sympy import symbols
>>> r, q = symbols('r, q', real=True)
>>> solve([r + q - 3, r > 3], r)
(3 < r) & Eq(r, 3 - q)
您可以直接调用solve
遇到关系时调用的例程:reduce_inequalities()
。它将 Expr 视为 Equality。
>>> from sympy import reduce_inequalities
>>> reduce_inequalities([x**2 - 4])
Eq(x, -2) | Eq(x, 2)
如果每个关系只包含一个感兴趣的符号,则可以为多个符号处理表达式:
>>> reduce_inequalities([0 <= x - 1, y < 3], [x, y])
(-oo < y) & (1 <= x) & (x < oo) & (y < 3)
但是,如果任何关系具有超过一个感兴趣的符号,则会引发错误:
>>> reduce_inequalities([0 <= x*y - 1, y < 3], [x, y])
Traceback (most recent call last):
...
NotImplementedError:
inequality has more than one symbol of interest.
禁用高阶显式解决方案
解决多项式表达式时,您可能不希望获得显式解决方案(这可能会很长)。如果表达式是单变量的,则将返回CRootOf
实例:
>>> solve(x**3 - x + 1)
[-1/((-1/2 - sqrt(3)*I/2)*(3*sqrt(69)/2 + 27/2)**(1/3)) -
(-1/2 - sqrt(3)*I/2)*(3*sqrt(69)/2 + 27/2)**(1/3)/3,
-(-1/2 + sqrt(3)*I/2)*(3*sqrt(69)/2 + 27/2)**(1/3)/3 -
1/((-1/2 + sqrt(3)*I/2)*(3*sqrt(69)/2 + 27/2)**(1/3)),
-(3*sqrt(69)/2 + 27/2)**(1/3)/3 -
1/(3*sqrt(69)/2 + 27/2)**(1/3)]
>>> solve(x**3 - x + 1, cubics=False)
[CRootOf(x**3 - x + 1, 0),
CRootOf(x**3 - x + 1, 1),
CRootOf(x**3 - x + 1, 2)]
如果表达式是多变量的,则可能不会返回任何解决方案:
>>> solve(x**3 - x + a, x, cubics=False)
[]
有时即使标志为 False,也会获得解决方案,因为可能对表达式进行了因式分解。在以下示例中,方程可以因式分解为线性因子和二次因子的乘积,因此获得了显式解决方案(这不需要解决三次表达式):
>>> eq = x**3 + 3*x**2 + x - 1
>>> solve(eq, cubics=False)
[-1, -1 + sqrt(2), -sqrt(2) - 1]
解决涉及根式的方程
由于 SymPy 使用主根的原则,某些根式方程的解将被错过,除非 check=False:
>>> from sympy import root
>>> eq = root(x**3 - 3*x**2, 3) + 1 - x
>>> solve(eq)
[]
>>> solve(eq, check=False)
[1/3]
在上述示例中,方程只有一个解。其他表达式将产生虚假根,必须手动检查;给奇次幂根提供负参数的根也需要特别检查:
>>> from sympy import real_root, S
>>> eq = root(x, 3) - root(x, 5) + S(1)/7
>>> solve(eq) # this gives 2 solutions but misses a 3rd
[CRootOf(7*x**5 - 7*x**3 + 1, 1)**15,
CRootOf(7*x**5 - 7*x**3 + 1, 2)**15]
>>> sol = solve(eq, check=False)
>>> [abs(eq.subs(x,i).n(2)) for i in sol]
[0.48, 0.e-110, 0.e-110, 0.052, 0.052]
第一个解为负数,因此必须使用real_root
来查看它是否满足表达式:
>>> abs(real_root(eq.subs(x, sol[0])).n(2))
0.e-110
如果方程的根不是实数,则需要更多的注意找到根,特别是对于高阶方程。考虑以下表达式:
>>> expr = root(x, 3) - root(x, 5)
我们将通过选择每个根式的第 1 个根在 x = 3 处构造此表达式的已知值:
>>> expr1 = root(x, 3, 1) - root(x, 5, 1)
>>> v = expr1.subs(x, -3)
solve
函数无法找到此方程的任何精确根:
>>> eq = Eq(expr, v); eq1 = Eq(expr1, v)
>>> solve(eq, check=False), solve(eq1, check=False)
([], [])
函数unrad
可以用于获得方程的一种形式,从而可以找到数值根:
>>> from sympy.solvers.solvers import unrad
>>> from sympy import nroots
>>> e, (p, cov) = unrad(eq)
>>> pvals = nroots(e)
>>> inversion = solve(cov, x)[0]
>>> xvals = [inversion.subs(p, i) for i in pvals]
尽管eq
或eq1
可以用于查找xvals
,但只有通过expr1
才能验证解决方案:
>>> z = expr - v
>>> [xi.n(chop=1e-9) for xi in xvals if abs(z.subs(x, xi).n()) < 1e-9]
[]
>>> z1 = expr1 - v
>>> [xi.n(chop=1e-9) for xi in xvals if abs(z1.subs(x, xi).n()) < 1e-9]
[-3.0]
另请参阅
rsolve
用于解决递归关系
dsolve
用于解决微分方程
sympy.solvers.solvers.solve_linear(lhs, rhs=0, symbols=[], exclude=[])
返回一个从f = lhs - rhs
派生的元组,其中一个是以下之一:(0, 1)
、(0, 0)
、(symbol, solution)
、(n, d)
。
解释
(0, 1)
表示f
与symbols中不在exclude中的符号无关。
(0, 0)
表示在给定的符号中没有方程的解。如果元组的第一个元素不为零,则函数保证依赖于symbols中的符号。
(symbol, solution)
其中符号在 f
的分子中线性出现,在 symbols 中(如果给定),且不在 exclude 中(如果给定)。除了 mul=True
扩展之外,f
不对 f
进行任何简化,因此解决方案将严格对应于唯一解决方案。
(n, d)
其中 n
和 d
是 f
的分子和分母,当分子不是任何感兴趣符号的线性时;除非找到该符号的解决方案(在这种情况下,第二个元素是解决方案,而不是分母)。
示例
>>> from sympy import cancel, Pow
f
与 symbols 中不在 exclude 中的符号无关:
>>> from sympy import cos, sin, solve_linear
>>> from sympy.abc import x, y, z
>>> eq = y*cos(x)**2 + y*sin(x)**2 - y # = y*(1 - 1) = 0
>>> solve_linear(eq)
(0, 1)
>>> eq = cos(x)**2 + sin(x)**2 # = 1
>>> solve_linear(eq)
(0, 1)
>>> solve_linear(x, exclude=[x])
(0, 1)
变量 x
在以下每个中作为线性变量出现:
>>> solve_linear(x + y**2)
(x, -y**2)
>>> solve_linear(1/x - y**2)
(x, y**(-2))
当 x
或 y
中不是线性时,然后返回分子和分母:
>>> solve_linear(x**2/y**2 - 3)
(x**2 - 3*y**2, y**2)
如果表达式的分子是一个符号,则如果该符号的解决方案将任何分母设为 0,则返回 (0, 0)
:
>>> eq = 1/(1/x - 2)
>>> eq.as_numer_denom()
(x, 1 - 2*x)
>>> solve_linear(eq)
(0, 0)
但是自动重写可能会导致分母中的符号出现在分子中,因此将返回一个解决方案:
>>> (1/x)**-1
x
>>> solve_linear((1/x)**-1)
(x, 0)
使用未评估的表达式来避免这种情况:
>>> solve_linear(Pow(1/x, -1, evaluate=False))
(0, 0)
如果允许 x
在以下表达式中取消,则似乎在 x
中是线性的,但 solve_linear
不会执行此类取消,因此解决方案将始终满足原始表达式,而不会引发除零错误。
>>> eq = x**2*(1/x - z**2/x)
>>> solve_linear(cancel(eq))
(x, 0)
>>> solve_linear(eq)
(x**2*(1 - z**2), x)
可以给出希望解决方案的符号列表:
>>> solve_linear(x + y + z, symbols=[y])
(y, -x - z)
还可以给出一个要忽略的符号列表:
>>> solve_linear(x + y + z, exclude=[x])
(y, -x - z)
(因为它是从符号的规范排序列表中第一个具有线性解的变量,所以得到了 y
的解。)
sympy.solvers.solvers.solve_linear_system(system, *symbols, **flags)
解决具有 (N) 个线性方程和 (M) 个变量的系统,这意味着支持欠定和过度定系统。
解释
可能的解数是零、一或无限。相应地,此过程将返回 None 或一个具有解的字典。在欠定系统的情况下,所有任意参数都将被跳过。这可能导致返回一个空字典的情况。在这种情况下,可以为所有符号分配任意值。
此函数的输入是一个 (N\times M + 1) 矩阵,这意味着它必须是增广形式。如果您喜欢输入 (N) 个方程和 (M) 个未知数,那么请使用 solve(Neqs, *Msymbols)
。注意:此例程会制作矩阵的本地副本,因此传递的矩阵不会被修改。
这里使用的算法是无分数高斯消元法,在消元后得到一个上三角矩阵。然后使用回代法找到解。这种方法比高斯-约旦方法更有效和更紧凑。
示例
>>> from sympy import Matrix, solve_linear_system
>>> from sympy.abc import x, y
解决以下系统:
x + 4 y == 2
-2 x + y == 14
>>> system = Matrix(( (1, 4, 2), (-2, 1, 14)))
>>> solve_linear_system(system, x, y)
{x: -6, y: 2}
一个退化系统返回一个空字典:
>>> system = Matrix(( (0,0,0), (0,0,0) ))
>>> solve_linear_system(system, x, y)
{}
sympy.solvers.solvers.solve_linear_system_LU(matrix, syms)
使用 LUsolve
解决增广矩阵系统,并返回一个字典,其中解按顺序键入为 syms 的符号。
解释
矩阵必须可逆。
示例
>>> from sympy import Matrix, solve_linear_system_LU
>>> from sympy.abc import x, y, z
>>> solve_linear_system_LU(Matrix([
... [1, 2, 0, 1],
... [3, 2, 2, 1],
... [2, 0, 0, 1]]), [x, y, z])
{x: 1/2, y: 1/4, z: -1/2}
另请参阅
LUsolve
sympy.solvers.solvers.solve_undetermined_coeffs(equ, coeffs, *syms, **flags)
解决一个由匹配coeffs
变量中的系数形成的包含(k)个参数的方程组,这些系数依赖于其余变量(或由syms
显式给出)。
解释
此函数的结果是一个字典,其中包含关于(q)中系数的符号值 - 如果没有解或系数不出现在方程中,则为空 - 否则为 None(如果系统未被识别)。如果存在多个解,解将作为列表传递。输出可以使用与solve
相同的语义进行修改,因为传递的标志直接发送到solve
,例如标志dict=True
将始终返回作为字典的解的列表。
此函数接受 Equality 和 Expr 类实例。当指定符号以确定要确定的参数时,解决过程最有效,但将尝试确定它们(如果不存在)。如果未获得预期的解决方案(并且未指定符号),请尝试指定它们。
示例
>>> from sympy import Eq, solve_undetermined_coeffs
>>> from sympy.abc import a, b, c, h, p, k, x, y
>>> solve_undetermined_coeffs(Eq(a*x + a + b, x/2), [a, b], x)
{a: 1/2, b: -1/2}
>>> solve_undetermined_coeffs(a - 2, [a])
{a: 2}
方程可以在符号中是非线性的:
>>> X, Y, Z = y, x**y, y*x**y
>>> eq = a*X + b*Y + c*Z - X - 2*Y - 3*Z
>>> coeffs = a, b, c
>>> syms = x, y
>>> solve_undetermined_coeffs(eq, coeffs, syms)
{a: 1, b: 2, c: 3}
系数在非线性的情况下也可以是非线性的,但如果只有一个解,则将其作为字典返回:
>>> eq = a*x**2 + b*x + c - ((x - h)**2 + 4*p*k)/4/p
>>> solve_undetermined_coeffs(eq, (h, p, k), x)
{h: -b/(2*a), k: (4*a*c - b**2)/(4*a), p: 1/(4*a)}
多个解始终作为列表返回:
>>> solve_undetermined_coeffs(a**2*x + b - x, [a, b], x)
[{a: -1, b: 0}, {a: 1, b: 0}]
使用标志dict=True
(符合solve()
中的语义)将强制结果始终为包含任何解的列表。
>>> solve_undetermined_coeffs(a*x - 2*x, [a], dict=True)
[{a: 2}]
sympy.solvers.solvers.nsolve(*args, dict=False, **kwargs)
数值解非线性方程组:nsolve(f, [args,] x0, modules=['mpmath'], **kwargs)
。
解释
f
是表示系统的符号表达式的向量函数。args是变量。如果只有一个变量,则可以省略此参数。x0
是接近解的起始向量。
使用 modules 关键字指定应用于评估函数和雅可比矩阵的模块。确保使用支持矩阵的模块。有关语法的更多信息,请参阅lambdify
的文档字符串。
如果关键字参数包含dict=True
(默认为 False),nsolve
将返回一个解映射的列表(可能为空)。如果想要使用nsolve
作为solve
的替代解决方案,这可能特别有用,因为对于两种方法都使用 dict 参数会产生一致类型结构的返回值。请注意:为了保持与solve
的一致性,即使nsolve
(目前至少)一次只找到一个解,解也将作为列表返回。
支持超定系统。
示例
>>> from sympy import Symbol, nsolve
>>> import mpmath
>>> mpmath.mp.dps = 15
>>> x1 = Symbol('x1')
>>> x2 = Symbol('x2')
>>> f1 = 3 * x1**2 - 2 * x2**2 - 1
>>> f2 = x1**2 - 2 * x1 + x2**2 + 2 * x2 - 8
>>> print(nsolve((f1, f2), (x1, x2), (-1, 1)))
Matrix([[-1.19287309935246], [1.27844411169911]])
对于一维函数,语法更简单:
>>> from sympy import sin, nsolve
>>> from sympy.abc import x
>>> nsolve(sin(x), x, 2)
3.14159265358979
>>> nsolve(sin(x), 2)
3.14159265358979
要以高于默认精度解决,请使用 prec 参数:
>>> from sympy import cos
>>> nsolve(cos(x) - x, 1)
0.739085133215161
>>> nsolve(cos(x) - x, 1, prec=50)
0.73908513321516064165531208767387340401341175890076
>>> cos(_)
0.73908513321516064165531208767387340401341175890076
要解决实函数的复根,必须指定非实数的初始点:
>>> from sympy import I
>>> nsolve(x**2 + 2, I)
1.4142135623731*I
使用mpmath.findroot
,您可以找到更详细的文档,特别是关于关键字参数和可用求解器的部分。但请注意,对于在根附近非常陡峭的函数,验证解可能会失败。在这种情况下,您应该使用标志verify=False
并独立验证解决方案。
>>> from sympy import cos, cosh
>>> f = cos(x)*cosh(x) - 1
>>> nsolve(f, 3.14*100)
Traceback (most recent call last):
...
ValueError: Could not find root within given tolerance. (1.39267e+230 > 2.1684e-19)
>>> ans = nsolve(f, 3.14*100, verify=False); ans
312.588469032184
>>> f.subs(x, ans).n(2)
2.1e+121
>>> (f/f.diff(x)).subs(x, ans).n(2)
7.4e-15
如果根的边界已知且使用二分法,可以安全地跳过验证:
>>> bounds = lambda i: (3.14*i, 3.14*(i + 1))
>>> nsolve(f, bounds(100), solver='bisect', verify=False)
315.730061685774
或者,当忽略分母时,函数可能行为更佳。然而,并非总是如此,因此使用哪个函数的决定留给用户决定。
>>> eq = x**2/(1 - x)/(1 - 2*x)**2 - 100
>>> nsolve(eq, 0.46)
Traceback (most recent call last):
...
ValueError: Could not find root within given tolerance. (10000 > 2.1684e-19)
Try another starting point or tweak arguments.
>>> nsolve(eq.as_numer_denom()[0], 0.46)
0.46792545969349058
sympy.solvers.solvers.checksol(f, symbol, sol=None, **flags)
检查 sol 是否为方程f == 0的解。
解释
输入可以是单个符号及其对应值,也可以是符号和值的字典。当作为字典给出且标志simplify=True
时,字典中的值将被简化。 f 可以是单个方程或方程的可迭代对象。解必须满足f中的所有方程才被认为是有效的;如果一个解不满足任何方程,则返回 False;如果一个或多个检查无法得出结论(且没有 False),则返回 None。
示例
>>> from sympy import checksol, symbols
>>> x, y = symbols('x,y')
>>> checksol(x**4 - 1, x, 1)
True
>>> checksol(x**4 - 1, x, 0)
False
>>> checksol(x**2 + y**2 - 5**2, {x: 3, y: 4})
True
要使用checksol()
检查表达式是否为零,将其作为f传递,并为symbol发送空字典:
>>> checksol(x**2 + x - x*(x + 1), {})
True
如果checksol()
无法得出结论,则返回 None。
标志:
‘numerical=True(默认)’
如果f
只有一个符号,则进行快速数值检查。
‘minimal=True(默认为 False)’
一个非常快速的、最小化的测试。
‘warn=True(默认为 False)’
如果checksol()
无法得出结论,则显示警告。
‘simplify=True(默认)’
在将解代入函数之前简化解,以及在尝试特定简化之前简化函数
‘force=True(默认为 False)’
使所有符号为正,不假设符号的符号。
sympy.solvers.solvers.unrad(eq, *syms, **flags)
去除具有符号参数的根并返回(eq, cov)
,None 或引发错误。
解释
如果没有根可去除,则返回 None。
如果有根并且无法去除或者无法解决重写系统为多项式所需的原始符号与变量变化关系,则引发 NotImplementedError。
否则返回元组(eq, cov)
,其中:
eq, cov
eq 是一个没有根的方程(在感兴趣的符号中),其解集是原始表达式的超集。 eq 可以用新变量重新表达;与原始变量的关系由cov
给出,其中包含v
和v**p - b
的列表,其中p
是清除根所需的幂,b
是现在用感兴趣的符号表示的根的多项式。例如,对于 sqrt(2 - x),元组将是(c, c**2 - 2 + x)
。 eq 的解将包含原方程的解(如果有的话)。
syms
如果提供了一个符号的可迭代对象,则限制根式的集中焦点:只有带有一个或多个感兴趣符号的根式将被清除。如果未设置syms,则使用所有自由符号。
flags在递归调用期间用于内部通信。还识别两个选项:
take
,如果定义了,则被解释为一个单参数函数,如果给定的 Pow 应该处理,则返回 True。
如果表达式中有根式,可以将其去除:
- 所有根的底数都相同;在这种情况下进行变量更改。
- 如果表达式中所有根都出现在一个项中。
- 只有四个带有 sqrt()因子的项,或者少于四个具有 sqrt()因子的项。
- 只有两个带有根式的项。
示例
>>> from sympy.solvers.solvers import unrad
>>> from sympy.abc import x
>>> from sympy import sqrt, Rational, root
>>> unrad(sqrt(x)*x**Rational(1, 3) + 2)
(x**5 - 64, [])
>>> unrad(sqrt(x) + root(x + 1, 3))
(-x**3 + x**2 + 2*x + 1, [])
>>> eq = sqrt(x) + root(x, 3) - 2
>>> unrad(eq)
(_p**3 + _p**2 - 2, [_p, _p**6 - x])
普通微分方程(ODEs)
参见 ODE。
偏微分方程(PDEs)
参见 PDE。
Deutils(用于解决 ODE 和 PDE 的实用程序)
sympy.solvers.deutils.ode_order(expr, func)
返回给定微分方程关于 func 的阶数。
此函数采用递归实现。
示例
>>> from sympy import Function
>>> from sympy.solvers.deutils import ode_order
>>> from sympy.abc import x
>>> f, g = map(Function, ['f', 'g'])
>>> ode_order(f(x).diff(x, 2) + f(x).diff(x)**2 +
... f(x).diff(x), f(x))
2
>>> ode_order(f(x).diff(x, 2) + g(x).diff(x, 3), f(x))
2
>>> ode_order(f(x).diff(x, 2) + g(x).diff(x, 3), g(x))
3
递归方程
sympy.solvers.recurr.rsolve(f, y, init=None)
解决具有有理系数的一元递归。
给定(k)阶线性递归(\operatorname{L} y = f),或等效地:
[a_{k}(n) y(n+k) + a_{k-1}(n) y(n+k-1) + \cdots + a_{0}(n) y(n) = f(n)]
其中(a_{i}(n)),对于(i=0, \ldots, k),是关于(n)的多项式或有理函数,而(f)是超几何函数或在(n)上具有固定数量不同超几何项的和,找到所有解或返回None
,如果没有找到。
初始条件可以作为字典的两种形式给出:
{ n_0 : v_0, n_1 : v_1, ..., n_m : v_m}
{y(n_0) : v_0, y(n_1) : v_1, ..., y(n_m) : v_m}
或作为值列表L
:
L = [v_0, v_1, ..., v_m]
其中L[i] = v_i
,对于(i=0, \ldots, m),映射到(y(n_i))。
示例
让我们考虑以下递归:
[(n - 1) y(n + 2) - (n² + 3 n - 2) y(n + 1) + 2 n (n + 1) y(n) = 0]
>>> from sympy import Function, rsolve
>>> from sympy.abc import n
>>> y = Function('y')
>>> f = (n - 1)*y(n + 2) - (n**2 + 3*n - 2)*y(n + 1) + 2*n*(n + 1)*y(n)
>>> rsolve(f, y(n))
2**n*C0 + C1*factorial(n)
>>> rsolve(f, y(n), {y(0):0, y(1):3})
3*2**n - 3*factorial(n)
参见
rsolve_poly
, rsolve_ratio
, rsolve_hyper
sympy.solvers.recurr.rsolve_poly(coeffs, f, n, shift=0, **hints)
给定具有多项式系数和非齐次方程(\operatorname{L} y = f)的(k)阶线性递归运算符(\operatorname{L}),其中(f)是一个多项式,在特征为零的域(K)上寻找所有多项式解。
算法执行两个基本步骤:
- 计算一般多项式解的度(N)。
- 找到所有(N)次或更低次的多项式(\operatorname{L} y = f)。
有两种方法来计算多项式解。如果度约束相对较小,即小于或等于递归的阶数,则使用未知系数法。这将得到一个带有(N+1)未知数的代数方程组。
另一种情况下,算法将初始方程转化为等价的方程,使得代数方程组仅有(r)个不定元。这种方法相对于朴素方法更为复杂,并由 Abramov、Bronstein 和 Petkovsek 共同发明。
可以将此处实现的算法推广到线性(q)-差分和微分方程的情况。
假设我们想计算到常数的第(m)个伯努利多项式。为此,我们可以使用(b(n+1) - b(n) = m n^{m-1})的递推关系,其解为(b(n) = B_m + C)。例如:
>>> from sympy import Symbol, rsolve_poly
>>> n = Symbol('n', integer=True)
>>> rsolve_poly([-1, 1], 4*n**3, n)
C0 + n**4 - 2*n**3 + n**2
参考文献
[R890]
S. A. Abramov, M. Bronstein 和 M. Petkovsek,关于线性算子方程的多项式解,见:T. Levelt,编,Proc. ISSAC ‘95,ACM Press,New York,1995,290-296。
[R891]
M. Petkovsek,具有多项式系数的线性递推的超几何解,J. Symbolic Computation,14 (1992),243-264。
[R892]
- Petkovsek, H. S. Wilf, D. Zeilberger,A = B,1996。
sympy.solvers.recurr.rsolve_ratio(coeffs, f, n, **hints)
给定具有多项式系数的阶为(k)的线性递推算子(\operatorname{L})和非齐次方程(\operatorname{L} y = f),其中(f)是多项式,我们寻求特征零域(K)上所有有理解。
此过程仅接受多项式,但如果您有兴趣解决有理系数的递推,则使用rsolve
,它将预处理给定的方程并使用多项式参数运行此过程。
此算法执行两个基本步骤:
- 计算可以作为方程(\operatorname{L} y = f)任何有理解的通用分母的多项式(v(n))。
- 通过替换(y(n) = u(n)/v(n))构造新的线性差分方程,并解出(u(n))找到其所有的多项式解。如果找不到任何解,则返回
None
。
此处实现的算法是原始 Abramov 算法的修订版,于 1989 年开发。新方法实现更简单,整体效率更佳。此方法可以轻松地适应(q)-差分方程的情况。
除了单独找到有理解外,此函数还是 Hyper 算法的重要部分,用于寻找递推的非齐次部分的特解。
示例
>>> from sympy.abc import x
>>> from sympy.solvers.recurr import rsolve_ratio
>>> rsolve_ratio([-2*x**3 + x**2 + 2*x - 1, 2*x**3 + x**2 - 6*x,
... - 2*x**3 - 11*x**2 - 18*x - 9, 2*x**3 + 13*x**2 + 22*x + 8], 0, x)
C0*(2*x - 3)/(2*(x**2 - 1))
参见
rsolve_hyper
参考文献
[R893]
S. A. Abramov,具有多项式系数的线性差分和(q)-差分方程的有理解,见:T. Levelt,编,Proc. ISSAC ‘95,ACM Press,New York,1995,285-289。
sympy.solvers.recurr.rsolve_hyper(coeffs, f, n, **hints)
给定具有多项式系数的阶为(k)的线性递推算子(\operatorname{L})和非齐次方程(\operatorname{L} y = f),我们寻求所有特征零域(K)上的超几何解。
不齐次部分可以是超几何的,也可以是一组固定数量的两两不相似超几何项的总和。
该算法执行三个基本步骤:
- 将 (\operatorname{L} y = f) 不齐次部分中类似的超几何项分组,并使用 Abramov 算法找到特解。
- 计算 (\operatorname{L}) 的生成集并找到其基础,以便所有解线性无关。
- 用与 (\operatorname{L}) 基础维数相等的任意常数形成最终解。
如果 a(n) 由具有多项式系数的一阶线性差分方程消灭,或者更简单地说,如果连续项比是有理函数,则其为超几何项。
此过程的输出是固定数量的超几何项的线性组合。但是,底层方法可以生成更大类别的解 - D’Alembert 项。
注意,此方法不仅计算不齐次方程的核,还将其缩减为基础,以便通过此过程生成的解决方案是线性无关的。
示例
>>> from sympy.solvers import rsolve_hyper
>>> from sympy.abc import x
>>> rsolve_hyper([-1, -1, 1], 0, x)
C0*(1/2 - sqrt(5)/2)**x + C1*(1/2 + sqrt(5)/2)**x
>>> rsolve_hyper([-1, 1], 1 + x, x)
C0 + x*(x + 1)/2
参考文献
[R894]
M. Petkovsek,带有多项式系数的线性递推超几何解,J. Symbolic Computation,14 (1992),243-264。
[R895]
- Petkovsek,H. S. Wilf,D. Zeilberger,A = B,1996。
多项式方程组
sympy.solvers.polysys.solve_poly_system(seq, *gens, strict=False, **args)
返回多项式方程组的解列表,否则返回 None。
参数:
seq: 列表/元组/集合
列出需要解决的所有方程
gens: 生成器
seq 的方程的生成器,我们希望得到解的生成器
strict: 布尔值(默认为 False)
如果 strict 为 True,则在已知解决方案可能不完整时会引发 NotImplementedError(这可能发生在不能用根式表示所有解的情况下)
args: 关键字参数
解方程的特殊选项。
返回:
List[Tuple]
元组列表,其中元素为按 gens 传递顺序解决方案的符号
无
当计算的基础仅包含地面时返回 None。
示例
>>> from sympy import solve_poly_system
>>> from sympy.abc import x, y
>>> solve_poly_system([x*y - 2*y, 2*y**2 - x**2], x, y)
[(0, 0), (2, -sqrt(2)), (2, sqrt(2))]
>>> solve_poly_system([x**5 - x + y**3, y**2 - 1], x, y, strict=True)
Traceback (most recent call last):
...
UnsolvableFactorError
sympy.solvers.polysys.solve_triangulated(polys, *gens, **args)
使用 Gianni-Kalkbrenner 算法解多项式系统。
该算法通过在地面域中计算一个 Groebner 基础,然后通过在地面域的适当构造的代数扩展中迭代地计算多项式因式来进行。
参数:
polys: 列表/元组/集合
列出需要解决的所有方程
gens: 生成器
polys 中我们想要解的方程的生成器
args: 关键字参数
解方程的特殊选项
返回:
List[Tuple]
元组列表。满足 polys 中列出方程的符号的解决方案
示例
>>> from sympy import solve_triangulated
>>> from sympy.abc import x, y, z
>>> F = [x**2 + y + z - 1, x + y**2 + z - 1, x + y + z**2 - 1]
>>> solve_triangulated(F, x, y, z)
[(0, 0, 1), (0, 1, 0), (1, 0, 0)]
参考文献
1. Patrizia Gianni,Teo Mora,使用 Groebner 基础代数解多项式方程组,应用代数,代数算法和纠错编码的 AAECC-5,LNCS 356 247–257,1989
丢番图方程(DEs)
参见丢番图方程
不等式
参见不等式求解器
线性规划(优化)
sympy.solvers.simplex.lpmax(f, constr)
返回线性方程f
在使用 Ge、Le 或 Eq 表达的线性约束下的最大值。
所有变量都是未约束的,除非受到约束。
示例
>>> from sympy.solvers.simplex import lpmax
>>> from sympy import Eq
>>> from sympy.abc import x, y
>>> lpmax(x, [2*x - 3*y >= -1, Eq(x+ 3*y,2), x <= 2*y])
(4/5, {x: 4/5, y: 2/5})
变量的负值是允许的,除非明确排除:
>>> lpmax(x, [x <= -1])
(-1, {x: -1})
如果为x
添加了非负约束,则没有可能的解决方案:
>>> lpmax(x, [x <= -1, x >= 0])
Traceback (most recent call last):
...
sympy.solvers.simplex.InfeasibleLPError: inconsistent/False constraint
参见
linprog
, lpmin
sympy.solvers.simplex.lpmin(f, constr)
返回在使用 Ge、Le 或 Eq 表达的线性约束下的线性方程f
的最小值。
所有变量都是未约束的,除非受到约束。
示例
>>> from sympy.solvers.simplex import lpmin
>>> from sympy import Eq
>>> from sympy.abc import x, y
>>> lpmin(x, [2*x - 3*y >= -1, Eq(x + 3*y, 2), x <= 2*y])
(1/3, {x: 1/3, y: 5/9})
变量的负值是允许的,除非明确排除,因此最小化x
对于x <= 3
是一个无约束问题,而以下问题有一个有界解:
>>> lpmin(x, [x >= 0, x <= 3])
(0, {x: 0})
如果没有指明x
是非负的,这个目标没有最小值:
>>> lpmin(x, [x <= 3])
Traceback (most recent call last):
...
sympy.solvers.simplex.UnboundedLPError:
Objective function can assume arbitrarily large values!
参见
linprog
, lpmax
sympy.solvers.simplex.linprog(c, A=None, b=None, A_eq=None, b_eq=None, bounds=None)
返回在给定约束A*x <= b
和A_eq*x = b_eq
下c*x
的最小化值。除非给出边界,否则变量在解中将具有非负值。
如果没有给出A
,那么系统的维度将由C
的长度确定。
默认情况下,所有变量都将是非负的。如果bounds
作为单个元组(lo, hi)
给出,则所有变量将被限制在lo
和hi
之间。使用 None 表示lo
或hi
在负或正方向上没有约束,例如(None, 0)
表示非正值。要设置单个范围,传递一个长度等于A
列数的列表,每个元素都是一个元组;如果只有少数变量取非默认值,则可以作为字典传递,键给出相应分配变量的列,例如bounds={2: (1, 4)}
将限制第三个变量的值在[1, 4]
范围内。
示例
>>> from sympy.solvers.simplex import linprog
>>> from sympy import symbols, Eq, linear_eq_to_matrix as M, Matrix
>>> x = x1, x2, x3, x4 = symbols('x1:5')
>>> X = Matrix(x)
>>> c, d = M(5*x2 + x3 + 4*x4 - x1, x)
>>> a, b = M([5*x2 + 2*x3 + 5*x4 - (x1 + 5)], x)
>>> aeq, beq = M([Eq(3*x2 + x4, 2), Eq(-x1 + x3 + 2*x4, 1)], x)
>>> constr = [i <= j for i,j in zip(a*X, b)]
>>> constr += [Eq(i, j) for i,j in zip(aeq*X, beq)]
>>> linprog(c, a, b, aeq, beq)
(9/2, [0, 1/2, 0, 1/2])
>>> assert all(i.subs(dict(zip(x, _[1]))) for i in constr)
参见
lpmin
, lpmax