来自 scipy.optimize 的 fsolve 失败

Question

我正在尝试计算向量值函数的根，但我从 fsolve 得到的解决方案是完全错误的。我确实收到了警告：175: RuntimeWarning: 迭代没有取得很好的进展，正如最近十次迭代的改进所衡量的那样。 warnings.warn（消息，运行时警告）

这是工作示例：


import numpy as np 
from scipy.optimize import fsolve


def AMOC(amoc_state,
         gamma= 1/0.34,
         theta = 1,
         mu = 7.5,
         sigma = 0.85):
    T = amoc_state[0]
    S = amoc_state[1]
    dT = -gamma * (T-theta) - T * (1+ mu*np.abs(T-S)) 
    dS = sigma-S*(1+mu*np.abs(T-S))
    return (dT, dS) 

test = fsolve(AMOC, (0.8,0.3), 2.3,xtol = 1e-7, factor = 0.1, maxfev = 1000)

这给了我 test = array([0.57372733, 0.47660354]) 这显然不是 AMOC 函数的根，因为

AMOC((0.57372733, 0.47660354), gamma = 2.3) returns (-0.01121948095040759, 0.026224895476309684)

非常感谢任何帮助！

Answer 1

正如评论中已经提到的，abs 函数不可微分，您的函数 AMOC 也是不可微分的。这很重要，因为 fsolve 的基础算法使用导数信息来求解 F(x) = 0，因此假设函数 F 是可微分的。

但是，作为解决方法，您可以将 abs(x) 替换为平滑近似值 sqrt(x**2 + 1.0e-4):

def AMOC(amoc_state, gamma= 1/0.34, theta = 1, mu = 7.5, sigma = 0.85):
    T, S = amoc_state
    smooth_abs = np.sqrt((T-S)**2 + 1.0e-4)
    dT = -gamma * (T-theta) - T * (1+ mu*smooth_abs) 
    dS = sigma-S*(1+mu*smooth_abs)
    return np.array((dT, dS))

接下来，您可以尝试沃伦已经提到的另一个初步猜测。否则，您可以通过scipy.optimize.minimize将问题作为等效最小化问题来解决。根据我的经验，求解器更强大：

from scipy.optimize import minimize

res = minimize(lambda x: np.sum(AMOC(x, 2.3)**2), x0=(0.8,0.3))

产生解决方案 array([2.25836418e-01, 1.38576534e-06])。

来自 scipy.optimize 的 fsolve 失败

fsolve from scipy.optimize fails

python

solver

scipy