在 Julia 中求解 ODE 系统，初始条件为二维数组

Question

http://diffeq.sciml.ai/latest/tutorials/ode_example.html

我正在尝试使用 Julia (DifferentialEquations.jl) 中的 ODE 求解器来求解 n 个相互作用粒子的系统。假设系统是二维的，每个粒子的运动方程由其位置相对于时间的二阶 ODE 来描述。然后每个粒子需要四个变量，两个用于位置，两个用于速度。然后需要声明 4n 个变量。有没有一种概括的方法，这样就不需要把所有 4n 个方程一一列出来？

例如：

http://diffeq.sciml.ai/latest/tutorials/ode_example.html#Example-2:-Solving-Systems-of-Equations-1

我尝试通过扩展u 和du 将上面链接中的洛伦兹方程稍微修改为 n 个粒子（这是一个非常非常粗略的尝试，因为我实际上不知道该怎么做）到二维数组。

using DifferentialEquations
using Plots
n = 4
function lorenz(du,u,p,t,i)
    du[i,1] = 10.0*(u[i,2]-u[i,1])*sum(u[1:n,1])
    du[i,2] = (u[i,1]*(28.0-u[i,3]) - u[i,2])*sum(u[1:n,1])
    du[i,3] = (u[i,1]*u[i,2] - (8/3)*u[i,3])*sum(u[1:n,1])
end

u0 = hcat([1.0;0.0;0.0], [0.0;1.0;0.0], [0.0;0.0;1.0])
tspan = (0.0,100.0)
prob = ODEProblem(lorenz,u0,tspan)
sol = solve(prob)

这毫无疑问是行不通的，但我希望你能明白我想要做什么。无论如何，ODE 求解器是否可以将 u 求解为二维数组（或其他可以实现类似目的的方法？）

Answer 1

问题不在于二维数组。例如用

替换你的 lorenz 定义

function lorenz(du,u,p,t)
    du[1,1] = 10.0*(u[1,2]-u[1,1])
    du[1,2] = (u[1,1]*(28.0-u[1,3]) - u[1,2])
    du[1,3] = (u[1,1]*u[1,2] - (8/3)*u[1,3])
end

会起作用的。

问题在于函数签名，不支持附加的i。如果您想解决洛伦兹振荡器网络，您必须使用具有相同签名的函数对其进行编码，例如lorenz_network!(du, u, p, t) 用于就地版本。在您的函数中的各个振荡器上放置一个循环，您就快到了。