Scipy.optimize.newton方法求解多变量方程组报错问题咨询及替代方案询问

嘿，我来帮你捋清楚这个问题哈！首先得明确一个关键点：scipy.optimize.newton本质上是为单变量优化设计的函数！你不带fprime的时候能跑通多变量，其实是因为它退化成了多变量版的割线法（本质是逐维度单独处理），但这并不是它的设计目标，所以当你传入针对多变量的雅可比矩阵时，函数内部的维度处理逻辑就会混乱，直接导致了你遇到的索引错误。

那该怎么正确用牛顿法求解多变量方程组呢？给你几个靠谱的方案：

方案一：用scipy专门的多变量方程组求解函数——scipy.optimize.root #

这个函数就是为多变量非线性方程组量身打造的，完美支持牛顿法，只需要指定method='newton'并传入你的雅可比矩阵函数就行。

把你的代码修改一下就能用：

import numpy as np
from scipy.optimize import root

def f(x):
    F = [0, 0]
    F[0] = 5 * x[0] + 7 * x[1] - 6
    F[1] = 10 * x[0] - 3 * x[1] - 46
    return F

def jacobian(x):
    return np.array([[5, 7], [10, -3]])  # 直接返回二维数组，写法更简洁

q = [10, -10]
result = root(f, q, jac=jacobian, method='newton', options={'maxiter': 10000})

print("求解结果：", result.x)
print("是否收敛：", result.success)
print("收敛信息：", result.message)

运行这个代码就能得到正确结果，而且确实是用牛顿法求解的，完全适配多变量场景。

方案二：手动实现牛顿-拉夫逊迭代 #

如果想自己完全控制迭代过程（比如后续处理更复杂的问题时需要自定义逻辑），也可以手动写牛顿法的核心逻辑，步骤很清晰：

1. 初始化迭代起点
2. 每次迭代计算当前点的函数值F(x)和雅可比矩阵J(x)
3. 解线性方程组J(x) * delta_x = -F(x)
4. 更新迭代点x = x + delta_x
5. 直到满足收敛条件（比如函数值的模足够小，或者迭代步长达标）

示例代码：

import numpy as np

def f(x):
    return np.array([5*x[0]+7*x[1]-6, 10*x[0]-3*x[1]-46])

def jacobian(x):
    return np.array([[5,7],[10,-3]])

# 初始化参数
x = np.array([10, -10], dtype=np.float64)
tolerance = 1e-8
max_iter = 10000
iter_count = 0

while iter_count < max_iter:
    F = f(x)
    # 判断是否收敛
    if np.linalg.norm(F) < tolerance:
        break
    J = jacobian(x)
    # 解线性方程组 J * delta_x = -F
    delta_x = np.linalg.solve(J, -F)
    x += delta_x
    iter_count += 1

print("手动牛顿法求解结果：", x)
print("迭代次数：", iter_count)

这个方法自由度很高，适合你后续拓展复杂问题时灵活调整。

其他替代方案 #

• 如果你需要处理的是最小二乘形式的问题，也可以尝试scipy.optimize.least_squares，它也支持传入雅可比矩阵，但对于你这种严格的方程组求解，root函数会更直接。
• 另外，scipy.optimize.fsolve也是常用的多变量方程组求解函数，默认用hybrd方法，如果你对求解方法没有严格的牛顿法要求，也可以试试这个。

备注：内容来源于stack exchange，提问作者Nisunt