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

Question

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

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

例如：

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

我尝试将上面 link 中的 Lorenz 方程修改为 n 个粒子（这是一个非常非常粗略的尝试，因为我实际上不知道如何去做）通过尝试扩展 u 和 du 到二维数组。

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

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

替换你的洛伦兹定义

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) 就地版本。在您的函数中的各个振荡器上放置一个循环，您就快完成了。

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

Solve systems of ODEs in Julia with initial condition as 2D array

arrays

differential-equations

julia