python 使用fsolve查找解

33qvvth1  于 2023-01-01  发布在  Python
关注(0)|答案(2)|浏览(198)
import numpy as np
from scipy.optimize import fsolve

musun = 132712000000
T = 365.25 * 86400 * 2 / 3
e = 581.2392124070273

def f(x):
    return ((T * musun ** 2 / (2 * np.pi)) ** (1 / 3) * np.sqrt(1 - x ** 2)
        - np.sqrt(.5 * musun ** 2 / e * (1 - x ** 2)))

x = fsolve(f, 0.01)
f(x)

print x

这个代码有什么问题吗?好像不起作用。

vhmi4jdf

vhmi4jdf1#

因为sqrt返回负参数NaN,所以函数f(x)对于所有真实的x都是不可计算的。我将函数改为numpy.emath.sqrt(),它可以在参数〈0时输出复数值,并返回表达式的绝对值。

import numpy as np
from scipy.optimize import fsolve
sqrt = np.emath.sqrt

musun = 132712000000
T = 365.25 * 86400 * 2 / 3
e = 581.2392124070273

def f(x):
    return np.abs((T * musun ** 2 / (2 * np.pi)) ** (1 / 3) * sqrt(1 - x ** 2)
        - sqrt(.5 * musun ** 2 / e * (1 - x ** 2)))

x = fsolve(f, 0.01)
x, f(x)

这样就可以得到正确的结果:

(array([ 1.]), array([ 121341.22302275]))

解非常接近真根,但是f(x)仍然非常大,因为f(x)具有非常大的因子:木孙。

8nuwlpux

8nuwlpux2#

fsolve()返回f(x) = 0的根(请参见here)。
当我绘制xf(x)值(范围-1到1)时,我发现在x = -1x = 1处有根。但是,如果是x > 1x < -1,则两个sqrt()函数都将传递一个负参数,这将导致错误invalid value encountered in sqrt
fsolve()找不到位于函数有效范围两端的根,这并不令我感到惊讶。
我发现,在试图找到一个函数的根之前,先画出它的图形总是一个好主意,因为这可以表明它被任何求根算法找到根的可能性有多大(或者在这种情况下,不太可能)。

相关问题