我正在尝试用Python在数值上求解这个方程式(numpy / scipy,所有功能都可用)



在这个公式中,K是常数,f和g是两个项,它们取决于E计数器(这是整数的离散表示),其中x是我要查找的变量。

举例来说,假设E为3个术语,则为:



f(E)和g(E)也是已知的。

我从numpy那里读到有关使用“ fsolve”的信息,尽管我不明白如何自动生成一个总和的函数。我可能会手动执行,但是要花50个学期才能用完一段时间,所以我也想学习一些新知识。

最佳答案

您可以使用scipy.optimize.fsolve,其中使用numpy.sum构造函数:

import numpy as np
import scipy.optimize

np.random.seed(123)

f = np.random.uniform(size=50)
g = np.random.uniform(size=f.size)
K = np.sum(f * np.exp(-g*np.pi))

def func(x, f, g, K):
    return np.sum(f * np.exp(-g*x), axis=0) - K

# The argument to the function is an array itself,
# so we need to introduce extra dimensions for f, g.
res = scipy.optimize.fsolve(func, x0=1, args=(f[:, None], g[:, None], K))


请注意,对于您的特定函数,您还可以通过提供函数的派生来辅助算法:

def derivative(x, f, g, K):
    return np.sum(-g*f * np.exp(-g*x), axis=0)

res = scipy.optimize.fsolve(func, fprime=derivative,
                            x0=1, args=(f[:, None], g[:, None], K))


寻找多个根

在某种意义上,该函数可以接受该过程,以使该函数接受N个输入(例如,对于每一行),并产生N个输出(对于每一行,仍为一个)。因此,输入和输出彼此独立,并且对应的雅可比矩阵是对角矩阵。这是一些示例代码:

import numpy as np
import scipy.optimize

np.random.seed(123)

image = np.random.uniform(size=(4000, 3000, 2))
f, g = image[:, :, 0], image[:, :, 1]
x_original = np.random.uniform(size=image.shape[0])  # Compute for each of the rows.
K = np.sum(f * np.exp(-g * x_original[:, None]), axis=1)

def func(x, f, g, K):
    return np.sum(f * np.exp(-g*x[:, None]), axis=1) - K

def derivative(x, f, g, K):
    return np.diag(np.sum(-g*f * np.exp(-g*x[:, None]), axis=1))

res = scipy.optimize.fsolve(func, fprime=derivative,
                            x0=0.5*np.ones(x_original.shape), args=(f, g, K))
assert np.allclose(res, x_original)

关于python - 寻找指数求和的解,我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/54258368/

10-10 20:21