使用scipy.newton求解含norm.cdf方程时迭代不收敛求助

我正在使用Python的scipy库求解一个包含norm.cdf的方程，示例代码如下：

import numpy as np
from scipy.stats import norm
from scipy.misc import derivative
from scipy import optimize

a_=1
b_=1
x0=2
def f(x, a=a_, b=b_):
    return (a*x*norm.cdf(b/x))

root = optimize.newton(f, x0, args=(a_, b_))
print(f"The root of the equation is approximately {root}")

但运行后出现如下错误：

RuntimeError: Failed to converge after 50 iterations, value is -0.12105344527421698.

问题核心原因 #

你的函数 f(x) = a*x*norm.cdf(b/x) 不存在实根：

• 正态分布累积分布函数 norm.cdf 的取值范围始终是 (0, 1)，永远不会等于0
• 当 a=1、b=1 时，只要 x≠0，a*x 也不会为0
  两者相乘的结果永远不为0，所以方程 f(x)=0 没有实数值解，这就是牛顿法无法收敛的根本原因。

解决方案（若需求为找极值点） #

如果你的真实需求是找函数的极值点（比如最小值/最大值），可以转而求解导数为0的点，示例代码如下：

import numpy as np
from scipy.stats import norm
from scipy import optimize

a_ = 1
b_ = 1
x0 = 2

def f(x, a=a_, b=b_):
    return a * x * norm.cdf(b/x)

# 手动推导的导数函数：f'(x) = a*norm.cdf(b/x) - a*b*norm.pdf(b/x)/x
def f_prime(x, a=a_, b=b_):
    cdf_val = norm.cdf(b/x)
    pdf_val = norm.pdf(b/x)
    return a * cdf_val - (a * b * pdf_val) / x

# 求解导数为0的点（即函数极值点）
extremum = optimize.newton(f_prime, x0, args=(a_, b_))
print(f"函数的极值点近似为: {extremum}")
print(f"该点的函数值为: {f(extremum)}")

运行后会得到函数在x>0区间的最小值点——原函数在x>0时先递减后递增，存在明确的极值；而x<0时函数值为负且单调变化，无极值点。

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