scipy中的复杂ODE系统
Posted
技术标签:
【中文标题】scipy中的复杂ODE系统【英文标题】:complex ODE systems in scipy 【发布时间】:2012-03-15 23:50:50 【问题描述】:我在求解光学布洛赫方程时遇到了问题,这是一个具有复数值的一阶 ODE 系统。我发现 scipy 可以解决这样的系统,但是他们的网页提供的信息太少,我几乎看不懂。
我有 8 个耦合的一阶 ODE,我应该生成如下函数:
def derv(y):
compute the time dervative of elements in y
return answers as an array
然后做complex_ode(derv)
我的问题是:
-
我的 y 不是列表而是矩阵,我怎样才能给出正确的输出
适合 complex_ode() 吗?
complex_ode()
需要一个雅可比,我不知道如何开始构建一个
应该是什么类型?
我应该把初始条件放在哪里,就像在正常颂歌和
时间线空间?
这是 scipy 的 complex_ode 链接: http://docs.scipy.org/doc/scipy/reference/generated/scipy.integrate.complex_ode.html
谁能给我更多的信息,以便我能学到更多。
【问题讨论】:
【参考方案1】:我认为我们至少可以为您指明正确的方向。光学的 布洛赫方程是一个科学界很好理解的问题 社区,虽然不是我的:-),所以互联网上已经有解决方案 针对这个特殊问题。
http://massey.dur.ac.uk/jdp/code.html
但是,为了满足您的需求,您谈到了使用 complex_ode,我想这是 很好,但我认为简单的 scipy.integrate.ode 也可以正常工作 根据他们的文档:
from scipy import eye
from scipy.integrate import ode
y0, t0 = [1.0j, 2.0], 0
def f(t, y, arg1):
return [1j*arg1*y[0] + y[1], -arg1*y[1]**2]
def jac(t, y, arg1):
return [[1j*arg1, 1], [0, -arg1*2*y[1]]]
r = ode(f, jac).set_integrator('zvode', method='bdf', with_jacobian=True)
r.set_initial_value(y0, t0).set_f_params(2.0).set_jac_params(2.0)
t1 = 10
dt = 1
while r.successful() and r.t < t1:
r.integrate(r.t+dt)
print r.t, r.y
你还有一个额外的好处,那就是更老的更成熟和更好的
记录的功能。我很惊讶你有 8 个而不是 9 个耦合 ODE,但我
确保您比我更了解这一点。是的,您是正确的,您的功能
应该是 ydot = f(t,y)
的形式,你称之为 def derv()
但你是
需要确保你的函数至少有两个参数
喜欢derv(t,y)
。如果您的y
在矩阵中,没问题!只需“重塑”它
derv(t,y)
函数如下:
Y = numpy.reshape(y,(num_rows,num_cols));
只要num_rows*num_cols = 8
,你的ODE数量应该没问题。然后
在计算中使用矩阵。完成后,请务必返回
一个向量,而不是像这样的矩阵:
out = numpy.reshape(Y,(8,1));
Jacobian 不是必需的,但它可能会允许计算继续进行 快得多。如果您不知道如何计算,您可能需要咨询 ***或微积分教科书。这很简单,但可能很耗时。
就初始条件而言,您可能应该已经知道那些应该是什么 是,无论它是复杂的还是真正有价值的。只要您选择的值 在合理的范围内,这应该无关紧要。
【讨论】:
非常感谢您的出色回答,一定花费了您很多时间。您提供的链接绝对完美,正是我正在寻找的。是的,你说得对,应该是 9 点 :) 好的,今天有时间,如果我在某个地方遇到了麻烦,我可以回到你身边吗?再次感谢。以上是关于scipy中的复杂ODE系统的主要内容,如果未能解决你的问题,请参考以下文章
使用 SymPy 表达式和 SciPy 求解器求解一阶 ODE 系统
通过Scipy和Numpy使用Python将数据拟合到ODE系统