是否可以在任何Python求解器中运行for循环？

Question

我试图找到一个只有一个未知数的可靠性工程非线性问题的解决方案。这应该相对容易些。但是，我遇到的问题在等式的求和部分之内。到目前为止，我的研究使我相信，不可能在Python的求解器中运行for循环。

我要解决的方程是：Equation

等式考虑了简单可靠性测试的结果。某些方程变量包含在列表中，而另一些则是单标量。方程变量为：

d =“修复有效性因子”，将通过纠正措施减轻的故障模式故障率的比例

v =测试期间采取纠正措施的时间

n =在每种故障模式下测试期间观察到的故障数

m =出现的故障模式总数

T =总测试时间（在这种情况下，测试175小时的25个单元= 4375小时）。

我正在尝试求解beta_hat的方程式。我已经使用其他方法/软件来求解方程，并且我知道beta_hat = 0.0016，但是，我需要使用Python来求解该方程，因为这是我在所有其他代码中所使用的。

保存每个方程式元素的值的列表是：

d = [0.5, 0., 0., 0., 0.8, 0., 0.7, 0.8, 0.5, 0.5, 0.5, 0., 0.]
v = [4375., 0., 0., 0., 4375., 0., 4375., 2500., 4375., 4375., 4375., 0., 0.]
n = [1, 3, 16, 2, 1, 4, 1, 3, 1, 1, 1, 8, 1]
m = 13
T = 4375

我尝试使用scipy.optimize（fsolve，root，least_squares）失败，但是我开始用尽了所有想法，我认为问题可能是我无法在Python的求解器中运行for循环像这个例子：

def f(x):
    for i in range(m):
        ((1 + x * T) / (x**2 * T)) * math.log(1 + x * T) * np.sum(n[i] / (1 / x + v[i] + (1 - d[i]) * (T - v[i]))) - m
    return f

result = optimize.root(f, 0.01)
print(result)

关于如何解决此问题的任何建议/想法？我是否有办法在求解器之外运行for循环？

Answer 1

正如我在评论中提到的，您的f()不会计算供求解器使用的数字。它似乎也没有实现链接到的方程式。这是一个完整的程序，包括一个人造玩具求解器：

def bisect(f, lo, hi, numiters):
    assert f(lo) < 0.0
    assert f(hi) > 0.0
    for _ in range(numiters):
        mid = lo + (hi - lo) / 2.0
        y = f(mid)
        if y < 0.0:
            lo = mid
        else:
            hi = mid
    return lo

d = [0.5, 0., 0., 0., 0.8, 0., 0.7, 0.8, 0.5, 0.5, 0.5, 0., 0.]
v = [4375., 0., 0., 0., 4375., 0., 4375., 2500., 4375., 4375., 4375., 0., 0.]
n = [1, 3, 16, 2, 1, 4, 1, 3, 1, 1, 1, 8, 1]
m = 13
T = 4375

def f(x):
    from math import log
    s = 0.0
    for i in range(m):
        s += n[i] / (1.0 / x + v[i] + (1.0 - d[i]) * (T - v[i]))
    return s * ((1.0 + x * T) / (x**2 * T)) * log(1.0 + x * T) - m

鉴于：

>>> f(5)
-12.979641343531911
>>> f(0.001)
3.975686546173577
>>> bisect(f, 5, 0.001, 53)
0.0016334432920776716

似乎与您期望的答案是0.0016相匹配。

因此，请考虑所有这一切：问题几乎肯定不是您要使用的求解器，而是您要传递给他们的函数。

Answer 2

给出Tim的评论和反馈，我修改了代码，并使用Solve和Root检查了答案。两次努力都达到了相同的结果（如预期）。

from math import log
from scipy.optimize import solve, root

d = [0.5, 0., 0., 0., 0.8, 0., 0.7, 0.8, 0.5, 0.5, 0.5, 0., 0.]
v = [4375., 0., 0., 0., 4375., 0., 4375., 2500., 4375., 4375., 4375., 0., 0.]
n = [1, 3, 16, 2, 1, 4, 1, 3, 1, 1, 1, 8, 1]
m = 13
T = 4375

def f(x):
    c = 0.0
    for i in range(m):
        c += n[i] / (1.0 / x + v[i] + (1.0 - d[i]) * (T - v[i]))
    return ((1.0 + x * T) / (x**2 * T)) * log(1.0 + x * T) * c - m

导致]

>>> fsolve(f, 0.001)
array([0.00163344])
>>> root(f, 0.001)
x: array([0.00163344])