Python牛顿方法-爆破

Question

我遇到了牛顿法的一些问题。我试着调试，但不知道根本原因是什么。我想我会把它发到这里看看是否有人能帮忙。

我有一个简单的函数来计算油密度，rho_p0，它先进行初始猜测，然后计算rho_a，然后迭代生成rho_p0的最终值。

from scipy.optimize import newton
def mccain_hill_rhop0(rho_p0, Rs,Sg,So):
        try:
            a0 = -49.8930
            a1 = 85.0149
            a2 = - 3.70373
            a3 = 0.0479818
            a4 = 2.98914
            a5 = - 0.0356888
            print(rho_p0)
            rho_a=(a0+a1*Sg+a2*Sg*rho_p0+a3*Sg*rho_p0**2+a4*rho_p0+a5*rho_p0**2)
            rho_p0=(Rs*Sg+4600*So)/(73.71+Rs*Sg/rho_a)
            return rho_p0
        except RuntimeError:
            return None
Rs=1000
Sg=0.55
So=0.8
rho_p01=52.8-0.01*Rs
rho_p0= newton(mccain_hill_rhop0, rho_p01, args=(Rs,Sg,So,), tol=10**(-3) , maxiter=100)
print(rho_p0)

这是错误消息。我不知道为什么会在第二次迭代之后爆炸。有什么想法吗？谢谢你的帮忙

RuntimeError: Tolerance of 2.5470966514398777e+33 reached. Failed to converge after 7 iterations, value is 2.5470966514398786e+33.
42.8
42.80438
-452.7973894722485
1393.5747708318825
672642.4222763113
1505488387.2917664
8.483536060177454e+17
2.5470966514398786e+33

Answer 1

我重新安排了返回函数，使我正在求解的变量在分子中工作。

from scipy.optimize import newton
def mccain_hill_rhop0(rho_p0, Rs,Sg,So):
        try:
            a0 = -49.8930
            a1 = 85.0149
            a2 = - 3.70373
            a3 = 0.0479818
            a4 = 2.98914
            a5 = - 0.0356888
            rho_a=(a0+a1*Sg+a2*Sg*rho_p0+a3*Sg*rho_p0**2+a4*rho_p0+a5*rho_p0**2)
            print((rho_p0, rho_a))
            rho_p0=rho_a*(Rs*Sg+4600*So-73.71*rho_p0)/(Rs*Sg)
            return rho_p0
        except RuntimeError:
            return None

输出

(42.8, 20.580650689599977)
(42.80438, 20.581334274602213)
(59.200734064755764, 20.63974342124446)
(57.3824271698717, 20.879752959108018)
(57.3871110909038, 20.879213691644452)
57.38705738842851