MATLAB is faster for my large system of ODEs?

function sys!(dx,x,p,t)
    mul!(dx, A, x)
    mul!(dx, B, u, 1, 1)
    nothing
end

this is the non-allocating version (depending a bit on where u is coming from).

If u is constant or piecewise constant, you can solve this differential equation by computing a single matrix exponential of exp([A B; 0 0]) and iterate the corresponding system

\begin{align} x(k+1) = A_d x(k) + B_d u(k) \end{align}
11 Likes