DifferentialEquations.jl Matrix ODE error

General Usage
1

Hi. I’m trying to use DifferentialEquations.jl, but I’m running into an elementary problem of setting up the ODE (i’m probably incompetent, as i’m new to Julia). I have as my test function a simple system of ODEs. I’m running Julia 1.4.2 and DifferentialEquations v6.14.0. I get as an error MethodError: no method matching *(::Adjoint{Float64,Array{Float64,2}}, ::DiffEqBase.NullParameters). Here is my MWE:

using DifferentialEquations
using LinearAlgebra

# Parameters
K1	= [1.0 0.0; 0.0 1.0]
H1	= reshape([.16 0.0 0.0 .25 .09 0.0 0.0 .04],(2,2,2)) # CUBE
rho1 = [1.0;0.5]

# ODE inputs
tspan = (0.0,10.0)
u0    = [1.0,0.0]

# Linear function of cube multiplication
cube_multiply_dd(cube::Array{Float64,3},u) = dropdims(mapslices(i-> u'*i*u, cube; dims=[1,2]);dims=(1,2)) 
# mapslices maps the quadratic function i-> u'*i*u to each slice of cube (sliced along 3rd dimension, thus each slice is a square matrix), then dropdims kills the singleton dimensions (1,2)

# ODE 
function ODE_beta2!(du,u,t)
	du = (rho1 - K1'*u - 1/2*cube_multiply_dd(H1,u));
end

cube_multiply_dd(H1,u0) # checking that cube_multiply works
ODE_beta2!(u0,u0,0) # checking that ODE_beta2! works
prob2 = ODEProblem(ODE_beta2!,u0,tspan)
sol2 = solve(prob2)

In the “system of equations” documentation of DifferentialEquations.jl, the function refers to each element of du (i.e. du[1] du[2] …) directly (see Lorenz example), whereas I just write abstractly du, so maybe i’m departing from recommended syntax but i can’t seem to find where it would say not to write it in my way.

2

Your function needs to have an input for parameters even if you’re not using them:

ODE_beta2!(du,u,p,t)

3

Also, if the function is in place, you’ll need do use broadcasted assignment like:

du .= (rho1 - K1'*u - 1/2*cube_multiply_dd(H1,u))

4

Thanks so much – that solved it!

i had written it before not using an in-place function (and thus wouldn’t have needed the broadcasting) and with an included parameter p, but i must have messed up something else.

Two further questions: (1) does the parameter p have to be in a certain position in the definition of the ODE function? (2) is my understanding correct that an in-place function is most likely faster (but needs the broadcasting operation .=)?

5

Yes, it needs to be in the second-to-last position.

Mostly yeah. But things like StaticArrays are really fast out-of-place for small arrays.