首页
学习
活动
专区
圈层
工具
发布
首页
学习
活动
专区
圈层
工具
MCP广场
社区首页 >问答首页 >需要对python参数有基本的了解,以便能够编写ODE解决方案的代码系统

需要对python参数有基本的了解,以便能够编写ODE解决方案的代码系统
EN

Stack Overflow用户
提问于 2017-07-24 12:14:02
回答 4查看 151关注 0票数 1

我是一个年长的Fortran程序员,需要帮助一个年轻人使用Heun的方法数值求解一个常微分方程系统。他只知道Python,所以我必须学习最基本的python才能做到这一点。

这是我想出来的。测试代码是一个简单的两自由度系统,每个自由度都有指数增长。

我得到的错误是:

代码语言:javascript
运行
复制
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

代码如下:

代码语言:javascript
运行
复制
# 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是全局的?即使它是全局的,为什么我不能将其作为参数传递?

EN

Stack Overflow用户

发布于 2017-07-24 15:36:14

好的,谢谢大家。这是一个代码,它与注释一起工作,以表明我学到了什么:

代码语言:javascript
运行
复制
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]
票数 1
EN
查看全部 4 条回答
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接:

https://stackoverflow.com/questions/45272514

复制
相关文章

相似问题

领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档