用python模拟矢量ODE
提问于 2019-09-17 09:49:50
我正在学习如何使用盖科 python包。作为第一步，我想模拟一个简单的向量ODE: dx/dt = A*x，其中A是矩阵，x是向量。我为GEKKO看到的所有ODE示例都是标量ODEs，在线文档中的数组示例没有说明在声明方程时如何合并.dt()方法。我知道，在声明方程时，可以使用列表，所以我认为这样的方法是可行的：
import numpy as np
from gekko import GEKKO
m=GEKKO()
m.time=np.linspace(0.,1.,10)
N=5
A=np.ones([N,N])
x=np.ones(N)
x=m.Var(value=x)
A=m.Param(value=A)
for i in range(N):
for j in range(N):
m.Equation(x[i].dt() += A[i][j] * x[j])
m.options.IMODE=4
m.solve()但是这段代码将失败有两个原因: 1) +=不是方程方法的有效比较，2) python抱怨xi.dt()不是x我的有效属性。那么，如果有可能的话，我将如何在GEKKO中模拟这首歌呢？
回答 1
Stack Overflow用户
回答已采纳
发布于 2019-09-17 20:42:04
模拟模型的一种方法是：
dx/dt =dx
将x声明为数组，并使用np.dot()对每一行A进行矩阵乘法。
import numpy as np
from gekko import GEKKO
m=GEKKO()
m.time=np.linspace(0.,1.,10)
N=5
A=np.ones((N,N))
ic = array([1., 1., 1., 1., 1.])
x=m.Array(m.Var,N,value=0.) #initialize to zero
for i in range(N):
x[i].value = ic[i] #set to some initial condition
m.Equations([x[i].dt()==np.dot(A[i,:],x) for i in range(N)])
m.options.IMODE=4
m.solve()
import matplotlib.pyplot as plt
for i in range(N):
plt.plot(m.time,x[i].value)
plt.show()另一种方法是使用Gekko中的状态空间物体
dx/dt =dx+B
y =C+D
包括B=0，C=0和D=0。这两种方法都得到了相同的结果。
import numpy as np
from gekko import GEKKO
m=GEKKO()
m.time=np.linspace(0.,1.,10)
N=5
A=np.ones((N,N))
B=np.zeros((N,1))
C=np.zeros((1,N))
x,y,u = m.state_space(A,B,C,D=None)
for i in range(N):
x[i].value=1
m.options.IMODE=4
m.solve()
import matplotlib.pyplot as plt
for i in range(N):
plt.plot(m.time,x[i].value)
plt.show()页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接：
https://stackoverflow.com/questions/57979580
复制
腾讯云开发者
