使用 Julia 中的 DifferentialEquations 包求解矩阵 ODE

Question

我想解决：

[\mathbf{M} \ddot{ \mathbf{U} }+ \mathbf{C} \dot{ \mathbf{U} }+ \mathbf{K} \mathbf{U} = \mathbf{P}(t) ]

或者，以状态空间形式：

[\dot{\mathbf{Y}}=f(\mathbf{Y},t)]

在哪里：

[\mathbf{Y} = \left[ \begin{array}{ c} \mathbf{U} \ \dot{ \mathbf{U} \end{array} \right] ]

和：

[f( \mathbf{Y} ,t)= \left[ \begin{array}{ c} \dot{ \mathbf{U} }\ \mathbf{M}^{-1} \mathbf{P} (t )- \mathbf{M} ^{-1} \mathbf{C} \dot{ \mathbf{U} }- \mathbf{M} ^{-1} \mathbf{K} \mathbf{U} \end{数组} \right] ]

我在 Julia 中尝试了以下代码，使用

\mathbf{M} = \left[ \begin{array}{ cc} 2&0\ 0&1 \end{array} \right];

\mathbf{C} = \left[ \begin{array}{ cc} 0&0\ 0& 0 \end{array} \right];

\mathbf{K} = \left[ \begin{array}{ cc} 96&-32\ -32& 32 \end{array} \right];

\mathbf{P} (t)= \left[ \begin{array}{ c} 10\ 10 \end{array} \right]

。

using DifferentialEquations
function eq(t,u,du)
    v=reshape([u...],Int(length(u)/2),2)
    du=reshape([v[:,2];-[10;10]-M\C*v[:,2]-M\K*v[:,1]],length(u))
end
u0=[0;0;0;0];
tspan=(0.0,10.0);
prob=ODEProblem(eq,u0,tspan)
sol=solve(prob)

但运行这些代码行会导致此错误：

ERROR: InexactError()

我正在使用朱莉娅版本。0.5.2。

请帮我。谢谢。

Answer 1

您的问题是 DifferentialEquations.jl 尊重您的输入类型。

u0=[0;0;0;0];

这是一个整数数组，因此这意味着您的问题将演变为一个整数数组。在第一步中，它发现计算返回浮点数，因此它不知道要做什么，u因为你说它必须是整数数组。解决这个问题的方法是说你的问题出在浮点数上：

u0=[0.0;0.0;0.0;0.0];

现在已经正确发展了。

但让我们再迈出一步。DifferentialEquations.jl 尊重您的输入类型，因此只需将初始条件设置为矩阵，DifferentialEquations.jl 就会使其成为矩阵问题。u0如果您制作矩阵，则无需重新整形：

u0=[0.0 0.0
    0.0 0.0]

一旦你这样做了，你只需编写一个直接适用于矩阵的方程。例子：

using DifferentialEquations
M = 1
K = 1
C = 1
function eq(t,u,du)
  du .= u[:,2] .- 10 - M./C.*u[:,2] - M.\K.*u[:,1]
end
u0=[0.0 0.0
    0.0 0.0]
tspan=(0.0,10.0);
prob=ODEProblem(eq,u0,tspan)
sol=solve(prob)

我不确定我是否得到了您试图解出的方程是否正确，因为它非常难以阅读，但这应该会让您非常接近您想要的结果。