Integration of a function that is numerical solution of differential equation

Question

I've obtained a numerical solution of a differential equation in a form of a vector (i.e., M(170,1)) by using ode45 (MATLAB) and tspan=0:500 (time range). How can I now integrate $M(t)$ with respect to time when I am not sure how time sequences are distributed? I know, that I could put tspan=0:0.01:500, but I would like to save some time now.

Answer 1

Consider augmented ODE of the form $$ M'(t) = F(t, M(t))\\ J'(t) = M(t)\\ M(0) = M_0\\ J(0) = 0 $$ Then $J(T)$ would be the $$ J(T) = \int_0^T M(t) dt. $$

Answer 2

Since you're using Matlab, you could use cumtrapz, which implements trapezoidal rule quadrature cumulatively:

Mint = cumtrapz(tout,M);

If you think you want something more, you could try combining cubic interpolation (or another interp1 option) with one of Matlab more advanced quadrature methods, e.g., integral:

Mint = zeros(length(tout),1);
for i = 2:length(tout)
    Mint(i) = integral(@(x)interp1(tspan,M,x,'cubic'),tout(1),tout(i))
    % Or maybe try:
    %Mint(i) = Mint(i-1)+integral(@(x)interp1(tout,M,x,'cubic'),tout(i-1),tout(i))
end

Both of these options will return the integral of your ode45 solution at each of the time points in the output time vector tout.