Hello! I am new to Julia. I want to solve the variational equation of an ODE without handwriting the Jacobi matrix. Here is an MME. Suppose that the ODE equations are the followings:
\dot{x}=y\\
\dot{y}=x-x^3-0.1y+0.1\cos(t)
Set f(x,y,t)=(y,x-x^3-0.1y+0.1\cos(t))^T. Given a initial condition (x_0,y_0) at t=0, suppose the solution is (x(t),y(t)) for t\in[0,2\pi]. Now I want to solve the following matrix ODE:
\dot{X}=d_{(x,y)}f(x(t),y(t),t)X
with initial condition X(0)=id. I can solve the above matrix ODE by using the following codes:
using OrdinaryDiffEq, StaticArrays
function vectorfield(u, p, t)
SA[u[2], u[1]-u[1]^3-0.1*u[2]+0.1*cos(t)] # vector field of the system
end
function timemap(u0)
tspan = (0.0, 2 * pi)
prob = ODEProblem(vectorfield, u0, tspan)
sol = solve(prob, Vern7())
function dv(X, p, t)
SA[0.0 0.1; 1-3*(sol(t)[1])^2 -0.1] * X # variational equation in matrix form
end
X0 = SA[1.0 0.0; 0.0 1.0] # initial condition of variational equation
probJ = ODEProblem(dv, X0, tspan)
solJ = solve(probJ, Vern7(), save_everystep=false)
(sol[end], solJ[end])
end
However, when the ODE has a large dimension, hand calculation of the Jacobi matrix is not a good choice. Can anyone show me how to solve the matrix ODE like the above example without hand calculation of the Jacobi matrix?