超越方程

Question

我遇到了一个需要解决的超越方程的问题：

K = K0*(exp(-t*B)/1+L*B)

变量'B'是未知的。作为第一步，我需要用下面的表达式来计算B：

B = (K0-1)/(L+t)

在第二步以及后面的每一步，我需要这样计算B：

B = -(1/t)*ln((1+L*B)/K0)

当B的前一个值和当前值之间的相对差异不超过1%时，迭代就停止了。最终得到的B应该让第一个方程的右边等于1。

我该如何用python来实现这个呢？我听说过scipy中的零点查找方法，但我更希望用一些普通的编码方式（这样我能更好地理解）。我尝试过使用while循环。我可以写一个循环来迭代，并在第一个方程中的K接近1.0时停止迭代：

kinf = 1.123456e+00

tau = 2.832995e+01

L2 = 3.745903e+00

i = 1

b2 = (kinf-1)/(L2+tau)

def iterate():
    b = b2
    i = 1
    print "iteration no.{:.1f}; B^2 = {:.06e}; delta = {:.03f}".format(i, b, kinf*((exp(-tau*b))/(1+L2*b)) - 1.0)
    while abs(kinf*((exp(-tau*b))/(1+L2*b))-1.0)>0.0001:
        b = -(1/tau)*log((1+L2*b)/kinf)
        i+=1
        print "iteration no.{:.1f}; B^2 = {:.06e}; delta = {:.03f}".format(i, b, kinf*((exp(-tau*b))/(1+L2*b)) - 1.0)

但是我不明白，如何比较B的前一个值和当前值。我想这个问题是经典问题之一，但我非常感谢任何帮助。

更新： 谢谢你的帮助！我现在做得对吗？

def iterate():
b0 = (kinf-1)/(L2+tau)
bold = b0
b = -(1/tau)*log((1+L2*b0)/kinf)
bnew = b
diff = ((bnew-bold)/bnew)*100
while abs(diff)>=0.01:
    print 'previous B^2 = {:.06e}'.format(bold)
    bnew = -(1/tau)*log((1+L2*bold)/kinf)
    print 'B^2 = {:.06e}'.format(bnew)
    diff = ((bnew-bold)/bnew)*100
    print 'delta = {:.06e}'.format(diff)
    bold = bnew

Answer 1

为了在循环中实现这个功能，你可以创建一个叫 previous_b 的变量，并把它初始化为 None（表示没有值）。然后在你的循环里，如果 previous_b 是 None 或者差值大于你设定的阈值，就继续执行循环。

kinf = 1.123456e+00

tau = 2.832995e+01

L2 = 3.745903e+00

i = 1

b2 = (kinf-1)/(L2+tau)

def iterate():
    previous_b = None
    b = b2
    i = 1
    print "iteration no.{:.1f}; B^2 = {:.06e}; delta = {:.03f}".format(i, b, kinf*((exp(-tau*b))/(1+L2*b)) - 1.0)
    while previous_b is None or abs(kinf*((exp(-tau*b))/(1+L2*b))-1.0)>0.0001:
        previous_b = b
        b = -(1/tau)*log((1+L2*b)/kinf)
        i+=1
        print "iteration no.{:.1f}; B^2 = {:.06e}; delta = {:.03f}".format(i, b, kinf*((exp(-tau*b))/(1+L2*b)) - 1.0)

Answer 2

在这一行不要重新定义 b（这样会丢失 b 原来的值）：

b = -(1/tau)*log((1+L2*b)/kinf) # the old value of b gets lost here

相反，你可以在你的 while 循环中这样做：

b_new = -(1/tau)*log((1+L2*b)/kinf)
b_delta = b_new - b
# do whatever you want with b_delta and finally override the old b with the new b
b = b_new