Ode45 error: Unable to meet integration tolerances without reducing step size below the smallest value allowed

I have the following function:

function V = fun_name(t,y)

b=1;
c=1;
u=3;

V = [2*y(2)*((y(1)).^2 + c^2) - b*sin(t); -2*y(1)*(y(2)).^2 + u*cos(t)];

end

And I have the following script:

tspan = 0:0.01:40;
initial_cond=[1 ; 0.5 ];
[tv,Yv]=ode45('fun_name',tspan,initial_cond);
plot(tv,Yv(:,2),'x');

However the graph plotted does not go beyond t=0.7, and I receive the error message:

Warning: Failure at t=6.867282e-01.  Unable to meet integration tolerances
without reducing the step size below the smallest value allowed (1.776357e-15)
at time t.

From the plot, it does not look like the function is encountering a singularity. I have also tried using ode15s and ode23s instead of ode45 but I get the same error.

I can obtain graphs by splitting the time intervals up by using tspan = 0:0.01:0.7, then tspan = 0.7:0.01:2.4, tspan = 2.4:0.01:3.8 etc, where the function plotted stops again at t=2.4, t=3.8, etc. However I was wondering if there is a way I can obtain the plot for the full range i.e. tspan = 0:0.01:40 without encountering this error.

Answer 1

[Torsten]

Integrating up to t=0.68, my plot of the solution clearly indicates a singularity in the solution for y(2).

Best wishes

Torsten.

[Guestsaint]

Hello,

I integrated up to t=0.68 and got the following plot:

I don't see y(2) approaching an infinite value around t=0.68 from the plot, is this the same form of graph you obtained?

Many thanks.

[Torsten]

Sorry, the above comment concerns y(1).

I get values of about 30000 for y(1) at t=0.68.

Best wishes

Torsten.

[Cristian Vargas]

Hello

I'm having the same problem with ode15s and I haven't able to fix the problem. Could you please help me? I have tried to fix the problem with odeset as option but I don't really understand.

The message error is the following:

Warning: RelTol has been increased to 2.22045e-14.

> In odearguments (line 130)

In ode15s (line 150)

Warning: Failure at t=1.830877e-05. Unable to meet integration tolerances without reducing the

step size below the smallest value allowed (5.421011e-20) at time t.

> In ode15s (line 730)

21113 successful steps

175 failed attempts

42374 function evaluations

1 partial derivatives

658 LU decompositions

42357 solutions of linear systems

%This is the code, ode15s calls a function name conoxozono

tspan=linspace(0,3600);

PMNO2=46;PMNO=30;PMO3=48;PMO=16;PMCO=28;

Vmolecular=22.41*(288/273)*(1013/755);%Litros

Cono2=100*(PMNO2/Vmolecular)*(6.022*10^23)/((100^3)*1000*1000*PMNO2);%molec/cm3

CoCo=604*(PMCO/Vmolecular)*(6.022*10^23)/((100^3)*1000*1000*PMCO);%molec/cm3

options = odeset('RelTol',5.421011e-21,'AbsTol',5.421011e-21,'OutputFcn',@odeplot,'Stats','on');

[t,x]=ode15s(@conoxozono,tspan,[Cono2 0 0 0 0 0 CoCo 0 0 0],options);

function dxdt=conoxozono(t,x)

O2=(0.745127/(0.082057*288))*6.022*10^23*0.21*(1000/(100^3));%Moléculas/cm3

M=(0.745127/(0.082057*288))*(6.022*10^23)*(1000/(100^3));%Moléculas/cm3

HR=70;%%

denvaps=12.83; %densidad de vapor de saturación a 15°C (g/m3)

H2O=(HR*denvaps/100)*6.022*10^23/(18*100^3);

k1=0.0001*2.7+0.0088;%s-1

k2=(6*10^-34)*((288/300)^-2.3);%cm^6/(molec^2*s)

k3=(2.2*10^-12)*exp(-1430/288);%cm^3/(molec*s)

k4= 0.0028*k1;%2*10^-6*2.7+5*10^-5; %s-1

k5=2.9*10^-11;%cm3/molec*sec

k6=2.2*10^-10;%cm3/molec*sec

k7=2.2*10^-13;%cm3/molec*sec

k8=3.7*10^-12*exp(240/288);%cm3/molec*sec

ko=2.6*10^-30*(300/288)^2.9;

kinf=6.70*10^-11*(300/288)^0.6;

fc=0.43;

k9=(kinf*ko*M/(kinf+ko*M))*fc^((1+(log(ko*M/kinf)^2))^-1);

dxdt=zeros(10,1);%zeros(10,1);

dxdt(1)=-k1*x(1)+k3*x(2)*x(3)+k8*x(4)*x(3)-k9*x(5)*x(1);%dno2

dxdt(2)=k1*x(1)-k3*x(2)*x(3)-k8*x(4)*x(3);%dno

dxdt(3)=k1*x(1)-k2*x(6)*O2*M +k5*x(7)*M;%do

dxdt(4)=k2*x(6)*O2*M-k3*x(2)*x(3)-k4*x(2);%do3

dxdt(5)=k4*x(2)-k5*x(7)*M-k6*x(7)*H2O;%d(o1d)

dxdt(6)=2*k6*x(7)*H2O-k7*x(8)*x(5)+k8*x(4)*x(3)-k9*x(5)*x(1);%doh

dxdt(7)=-k7*x(8)*x(5);%dco

dxdt(8)=k7*x(8)*x(5);%dco2

dxdt(9)=k7*x(8)*x(5)-k8*x(4)*x(3);%dho2

dxdt(10)=k9*x(5)*x(1);%dhno3

end

[darova]

I reduced tspan

tspan = [0 1e-5];

No need such precision. I reduced steps of integration

options = odeset('RelTol',5.421011e-6,'AbsTol',5.421011e-6);%'OutputFcn',@odeplot,'Stats','on');
[t,x]=ode15s(@conoxozono,tspan,[Cono2 0 0 0 0 0 CoCo 0 0 0],options);
plot(t,x)

The result