Usage of OrdinaryDiffEqLinear.jl for non-autonomous problems

Are you using MatrixOperator correctly? There is update_func and update_func! and the latter is the one you should be using I think.

using OrdinaryDiffEq
using Plots

alpha = 0.1; t0 = 0.0; t1 = 10.0; dt = 0.01
tspan = (t0, t1) 
u0 = ones(2) # initial condition
A0 = [0.0 1.0; 0.0 0.0] # time-independent part

function update_func!(A, u, p, t)
    A[1, 2] = -alpha * t 
    A[2, 1] = alpha * t
end

# Magnus solver
A = MatrixOperator(A0; update_func!)
prob = ODEProblem(A, u0, tspan)
sol = solve(prob, MagnusGL8(); dt)
plot(sol, label=["u[1] (magnus)" "u[2] (magnus)"], seriestype=:line, lw=7, opacity=0.3)

# exact solutions
F(t0,t1) = alpha * (t1^2 - t0^2) / 2
u_exact1(t) = cos(F(t0, t)) * u0[1] - sin(F(t0, t)) * u0[2]
u_exact2(t) = sin(F(t0, t)) * u0[1] + cos(F(t0, t)) * u0[2]
plot!(sol.t, u_exact1.(sol.t), label="u[1] (exact)")
plot!(sol.t, u_exact2.(sol.t), label="u[2] (exact)")

EDIT: Ah I just saw that you also wrote this in the linked issue. I’ll just leave this here for posterity.

1 Like