数值求解一个或多个方程

使用 SymPy 数值求解一个或多个方程组。例如，数值求解 \(\cos(x) = x \) 返回 \( x \approx 0.739085133215161\)。

数值求解在以下情况下是有用的：

solve() 和 solveset() 不会尝试找到数值解，只会寻找数学上精确的符号解。因此，如果你想要数值解，请使用 nsolve()。

SymPy 是为符号数学设计的。如果你不需要进行符号操作，那么对于数值操作，你可以使用另一个免费且开源的包，如 NumPy 或 SciPy，它们会更快，支持数组操作，并且实现了更多的算法。使用 SymPy（或其依赖 mpmath）进行数值计算的主要原因是：

  • 在符号计算的上下文中使用 SymPy 进行简单的数值计算

  • 如果你需要任意精度的能力来获得比从 float64 中得到的更多精度的位数。

考虑的替代方案

数值求解方程的示例

以下是一个数值求解一个方程的示例：

>>> from sympy import cos, nsolve, Symbol
>>> x = Symbol('x')
>>> nsolve(cos(x) - x, x, 1)
0.739085133215161

指导

支持过定方程组。

查找实函数的复根

要解实函数的复根，请指定一个非实数（纯虚数或复数）的初始点：

>>> from sympy import nsolve
>>> from sympy.abc import x
>>> nsolve(x**2 + 2, 1) # Real initial point returns no root
Traceback (most recent call last):
    ...
ValueError: Could not find root within given tolerance. (4.18466446988997098217 > 2.16840434497100886801e-19)
Try another starting point or tweak arguments.
>>> from sympy import I
>>> nsolve(x**2 + 2, I) # Imaginary initial point returns a complex root
1.4142135623731*I
>>> nsolve(x**2 + 2, 1 + I) # Complex initial point returns a complex root
1.4142135623731*I

确保根在给定区间内

不能保证 nsolve() 会找到最接近初始点的根。这里，尽管根 -1 更接近初始点 -0.1，但 nsolve() 找到了根 1：

>>> from sympy import nsolve
>>> from sympy.abc import x
>>> nsolve(x**2 - 1, -0.1)
1.00000000000000

你可以通过使用 solver='bisect' 并在元组中指定区间，来确保找到的根在给定的区间内（如果存在这样的根）。在这里，指定区间 (-10, 0) 确保找到根 -1：

>>> from sympy import nsolve
>>> from sympy.abc import x
>>> nsolve(x**2 - 1, (-10, 0), solver='bisect')
-1.00000000000000

数值求解方程组

要解决一个多维函数系统，请提供一个元组。

  • 函数 (f1, f2)

  • 待求解的变量 (x1, x2)

  • 起始值 (-1, 1)

>>> from sympy import Symbol, nsolve
>>> 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]])

提高解决方案的精度

你可以使用 prec 来提高解决方案的精度：

>>> from sympy import Symbol, nsolve
>>> 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), prec=25))
Matrix([[-1.192873099352460791205211], [1.278444111699106966687122]])

创建一个可以用 SciPy 解决的函数

如上所述，SymPy 专注于符号计算，并未针对数值计算进行优化。如果你需要多次调用数值求解器，使用针对数值计算优化的求解器（如 SciPy 的 root_scalar()）会快得多。推荐的流程是：

  1. 使用 SymPy 生成（通过符号化简化或求解方程）数学表达式

  2. 使用 lambdify() 将其转换为 lambda 函数

  3. 使用诸如 SciPy 的数值库来生成数值解

>>> from sympy import simplify, cos, sin, lambdify
>>> from sympy.abc import x, y
>>> from scipy.optimize import root_scalar
>>> expr = cos(x * (x + x**2)/(x*sin(y)**2 + x*cos(y)**2 + x))
>>> simplify(expr) # 1. symbolically simplify expression
cos(x*(x + 1)/2)
>>> lam_f = lambdify(x, cos(x*(x + 1)/2)) # 2. lambdify
>>> sol = root_scalar(lam_f, bracket=[0, 2]) # 3. numerically solve using SciPy
>>> sol.root
1.3416277185114782

使用解决方案结果

将结果代入表达式

最佳实践是使用 evalf() 将数值代入表达式。以下代码演示了数值并不是精确的根，因为将其代回表达式会产生略微不同于零的结果：

>>> from sympy import cos, nsolve, Symbol
>>> x = Symbol('x')
>>> f = cos(x) - x
>>> x_value = nsolve(f, x, 1); x_value
0.739085133215161
>>> f.evalf(subs={x: x_value})
-5.12757857962640e-17

使用 subs 可能会由于精度错误导致结果不正确，这里实际上将 -5.12757857962640e-17 四舍五入为零：

>>> f.subs(x, x_value)
0

在替换值时，您也可以将一些符号保留为变量：

>>> from sympy import cos, nsolve, Symbol
>>> x = Symbol('x')
>>> f = cos(x) - x
>>> x_value = nsolve(f, x, 1); x_value
0.739085133215161
>>> y = Symbol('y')
>>> z = Symbol('z')
>>> g = x * y**2
>>> values = {x: x_value, y: 1}
>>> (x + y - z).evalf(subs=values)
1.73908513321516 - z

并非所有方程都能被求解

nsolve() 是一个数值求解函数，因此它通常能为无法代数求解的方程提供解决方案。

无解的方程

有些方程没有解，在这种情况下，SymPy 可能会返回一个错误。例如，方程 \(e^x = 0\)（在 SymPy 中为 exp(x)）没有解：

>>> from sympy import nsolve, exp
>>> from sympy.abc import x
>>> nsolve(exp(x), x, 1, prec=20)
Traceback (most recent call last):
...
ValueError: Could not find root within given tolerance. (5.4877893607115270300540019e-18 > 1.6543612251060553497428174e-24)
Try another starting point or tweak arguments.

报告一个错误

如果你在使用 nsolve() 时发现了一个错误，请在 SymPy 邮件列表 上发布这个问题。在问题解决之前，你可以使用 Alternatives to Consider 中列出的其他方法。