我是一个年长的Fortran程序员,需要帮助一个年轻人使用Heun的方法数值求解一个常微分方程系统。他只知道Python,所以我必须学习最基本的python才能做到这一点。
这是我想出来的。测试代码是一个简单的两自由度系统,每个自由度都有指数增长。
我得到的错误是:
Traceback (most recent call last):
File "program.py", line 28, in <module>
heun (imax, y, dt, t)
File "program.py", line 7, in heun
rhs(y, ydot)
NameError: global name 'ydot' is not defined代码如下:
# Test RHS of ODE system:
def rhs (y, ydot):
ydot[0] = y[0]
ydot[1] = y[1]
return;
# Does one step of Heun's method:
def heun (imax, y, dt, t):
rhs(y, ydot)
for i in range(0, imax):
y_tilde[i] = y[i] + dt * ydot[i]
rhs(y_tilde, ydot_at_tilde)
for i in range(0, imax):
y[i] = y[i] + dt/2 * (ydot[i] + ydot_at_tilde[i])
t = t + dt
return;
# Initial condition
y = [0, 0]
t = 0
dt = 0.01
nsteps = 100
imax = 1
istep = 1
while istep <= nsteps:
heun (imax, y, dt, t)
print istep, y[0], y[1]
istep = istep + 1问:为什么python认为例程heun中的对象ydot是全局的?即使它是全局的,为什么我不能将其作为参数传递?
发布于 2017-07-24 15:36:14
好的,谢谢大家。这是一个代码,它与注释一起工作,以表明我学到了什么:
def rhs (y, ydot):
ydot[0] = y[0]
ydot[1] = y[1]
return ydot;
def heun (ndof, dt, y, t):
# These initializations of local arrays take the place of Fortran declarations:
ydot = [0] * (ndof)
y_tilde = [0] * (ndof)
ydot_at_tilde = [0] * (ndof)
ydot = rhs(y, ydot)
# Note: In python range means:
# range (first element, upto but not including last element)
for i in range(0, ndof):
y_tilde[i] = y[i] + dt * ydot[i]
ydot_at_tilde = rhs(y_tilde, ydot_at_tilde)
for i in range(0, ndof):
y[i] = y[i] + dt/2 * (ydot[i] + ydot_at_tilde[i])
t = t + dt
# Note: This lists the output arguments:
return y, t;
# Initial condition
y = [1, 1]
t = 0
dt = 0.01
nsteps = 100
ndof = 2
istep = 1
while istep <= nsteps:
# Note: This is how you get the output arguments:
y, t = heun (ndof, dt, y, t)
istep = istep + 1
print t, y[0], y[1]https://stackoverflow.com/questions/45272514
复制相似问题