Python使用sympy.nonlinsolve求解三元非线性方程组无收敛解报错如何解决

你遇到的收敛报错来源于两方面：

• 输入参数都是浮点数，sympy.nonlinsolve默认的符号求解逻辑对浮点精度的适配性差
• 三个点几乎在同一水平线上，外接圆圆心的z坐标绝对值很大，默认求解步长和迭代次数不足以收敛到正确根

方案1：先消元转线性方程组求解（最稳定，无收敛问题） #

三个方程都是圆的标准方程，两两相减即可消去二次项，得到两个关于x_m、z_m的线性方程，求解后直接计算r即可：

import sympy as sp

x_1 = 0.0
z_1 = 1.0
x_2 = 15.81
z_2 = 0.99
x_3 = 23.8
z_3 = 0.98

x_m, z_m = sp.symbols('x_m, z_m', real=True)
# 两两方程相减消去r²和二次项
Eq_1 = sp.Eq((x_1 - x_m) ** 2 + (z_1 - z_m) ** 2, (x_2 - x_m) ** 2 + (z_2 - z_m) ** 2)
Eq_2 = sp.Eq((x_2 - x_m) ** 2 + (z_2 - z_m) ** 2, (x_3 - x_m) ** 2 + (z_3 - z_m) ** 2)
# 求解线性方程组
sol = sp.solve([Eq_1, Eq_2], [x_m, z_m])
# 计算半径r
r = sp.sqrt((x_1 - sol[x_m])**2 + (z_1 - sol[z_m])**2)
print(f"圆心坐标: ({sol[x_m]}, {sol[z_m]})")
print(f"半径: {r}")

方案2：使用数值求解器nsolve指定初始猜测值 #

如果要直接求解原始非线性方程组，用nsolve并给出合理的初始值即可，该场景下三个点近似水平，圆心z坐标很大，初始值可以设为r=1000, x_m=12, z_m=1000：

import sympy as sp

x_1 = 0.0
z_1 = 1.0
x_2 = 15.81
z_2 = 0.99
x_3 = 23.8
z_3 = 0.98

r, x_m, z_m = sp.symbols('r, x_m, z_m', real=True)
Eq_1 = sp.Eq((x_1 - x_m) ** 2 + (z_1 - z_m) ** 2 - r ** 2, 0)
Eq_2 = sp.Eq((x_2 - x_m) ** 2 + (z_2 - z_m) ** 2 - r ** 2, 0)
Eq_3 = sp.Eq((x_3 - x_m) ** 2 + (z_3 - z_m) ** 2 - r ** 2, 0)
# 指定初始猜测值求解
ans = sp.nsolve([Eq_1, Eq_2, Eq_3], [r, x_m, z_m], [1000, 12, 1000])
print(ans)

方案3：按报错提示调整nonlinsolve参数 #

将输入的浮点数转成sympy有理数避免精度干扰，同时调大迭代次数maxsteps：

import sympy as sp

# 转成sympy有理数
x_1 = sp.Rational(0)
z_1 = sp.Rational(1,1)
x_2 = sp.Rational('15.81')
z_2 = sp.Rational('0.99')
x_3 = sp.Rational('23.8')
z_3 = sp.Rational('0.98')

r, x_m, z_m = sp.symbols('r, x_m, z_m', real=True)
Eq_1 = sp.Eq((x_1 - x_m) ** 2 + (z_1 - z_m) ** 2 - r ** 2, 0)
Eq_2 = sp.Eq((x_2 - x_m) ** 2 + (z_2 - z_m) ** 2 - r ** 2, 0)
Eq_3 = sp.Eq((x_3 - x_m) ** 2 + (z_3 - z_m) ** 2 - r ** 2, 0)
ans = sp.nonlinsolve([Eq_1, Eq_2, Eq_3], [r, x_m, z_m], dict=True, maxsteps=200)
print(ans)

内容的提问来源于stack exchange，提问作者Digital_Citizen