Issue integrating inside ode45 loop

Question

When I try to use this Matlab code it goes in an infinite loop. I am trying to perform integration inside ode45:

clear
clc
options = odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-4 1e-5]);
[T,Y] = ode45(@rigid,[0 12],[0 1 1],options);

plot(T,Y(:,1),'+',T,Y(:,2),'*',T,Y(:,3),'.')


function dy = rigid(t,y)
dy = zeros(3,1);    % a column vector

dy(1) = y(2) ;

dy(2) = -y(1) * y(3);
fun = @(t) exp(-t.^2).*log(t).^2+y(1);
q = integral(fun,0,Inf);
dy(3) = y(2) * y(3) + q;

Answer 1

There is no "infinite loop." Your function just takes a very long time to integrate. Try setting tspan to [0 1e-7]. It appears to be a high frequency oscillation, but I don't know if your equations are correct (that's a math question rather than a programming one). Such systems are hard to integrate accurately (ode15 might be a better choice), let alone quickly.

You also didn't bother to mention the important fact that the call to integral is generating a warning message:

Warning: Minimum step size reached near x = 1.75484e+22. There may be a
singularity, or the tolerances may be too tight for this problem. 
> In funfun/private/integralCalc>checkSpacing at 457
  In funfun/private/integralCalc>iterateScalarValued at 320
  In funfun/private/integralCalc>vadapt at 133
  In funfun/private/integralCalc at 84
  In integral at 88
  In rtest1>rigid at 17
  In ode15s at 580
  In rtest1 at 5 

Printing out warning messages on each iteration greatly slows down integration. There's a good reason for this warning. You do realize that the function that you're evaluating with integral from 0 to Inf is equivalent to the following, right?

sqrt(pi)*((eulergamma + log(4))^2/8 + pi^2/16) + Inf*y(1)

where eulergamma is psi(1) or double(sym('eulergamma')). Your integrand doesn't converge.

If you like, can try to avoid the warning message in one of two ways.

1. Turn off the the warning (being sure to re-enable it afterwards). You can do that with the following code:

...
warning('OFF','MATLAB:integral:MinStepSize');
[T,Y] = ode45(@rigid,[0 12],[0 1 1],options);
warning('ON','MATLAB:integral:MinStepSize');
...

You can obtain the the ID for a waring via the lastwarn function.

2. The other option might be to change your integration bounds and avoid the warning altogether, e.g.:

...
q = integral(fun,0,1e20);
...

This may or may not be acceptable, but integral is not be returning the correct solution either because the result doesn't converge.