ode45 does not solve on the specified time interval. How do I fix this?

Hi, I have a system of differential equations that I want to solve. However, when trying to solve this DiffEq with ode45(@(t,y)odefun(t,y),tspan,y0), where I have put my odefun in a different file, it does produce results, but on a different time interval than I want. Namely, it solves for t = [0, 0.46]*10^-6, while I specified it should solve for t = [0,1].

My code: tspan = [0,1];

y0 = [0 110 -250 15];

[t,Xsolved] = ode45(@(t,y)odefun(t,y),tspan,y0);

odefun (very lengthy equations): function dXdt = odefun(t,y)

t_f = 30;

eq1 = -(82809797410786635*tan(log((- 74530005132491619276355244274867242431640625*y(4)^2*pi^2 - 8612356148643476394993094470262377676800000000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5) - 248801399849700440447422519369183536317307289600000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5) - 12261830142026331516929221358481899097130334070767616*y(2)^2*y(3)^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5))^(1/2)/(pi*y(4)*8633076226496069765625i + y(4)*pi*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)*498799959753106275696640000i - 110733148343331824175415296*y(2)*y(3)*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)))*1i))/(281474976710656*y(2)); eq2 = - (10019119168824577875*y(2)^2*(1361735765217477681/(28823037615171174400000*(9/10 - (17592186044416*y(2))/5918609519667053)^(27/10)) + (8450000000000*(8500259669165361/(9444732965739290427392*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)) + 13/250))/(250119*y(2)^4*cos(log((- 74530005132491619276355244274867242431640625*y(4)^2*pi^2 - 8612356148643476394993094470262377676800000000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5) - 248801399849700440447422519369183536317307289600000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5) - 12261830142026331516929221358481899097130334070767616*y(2)^2*y(3)^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5))^(1/2)/(y(4)*pi*8633076226496069765625i + y(4)*pi*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)*498799959753106275696640000i - 110733148343331824175415296*y(2)*y(3)*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)))*1i)^2) + 361283/50000000))/147573952589676412928 - (133588255584327705*((4332790137498830962146934784000*y(2)^2)/3184539876935769439309088518619 - (3982870920455782400*y(2))/5918609519667053 + 288000)*((1266637395197952*y(2))/29593047598335265 - (19652298123655411864023597056*y(2)^2)/175149693231467319161999868524045 + (397449804563656125325221541480305270980608*y(2)^3)/1036642641726526473848753494572130715682924789385 + (5520653160719109*y(4)*((4332790137498830962146934784000*y(2)^2)/3184539876935769439309088518619 - (3982870920455782400*y(2))/5918609519667053 + 288000))/731834939447705600000 + 33/2))/(590295810358705651712*((574685827824708321884380135181365943857027170818850816*y(2)^4)/557771182528674706092576514838123659160733159500324835031693855 - (1050791949051857975174900787749300236976128*y(2)^3)/1036642641726526473848753494572130715682924789385 + (119770698800860541596490268672*y(2)^2)/175149693231467319161999868524045 - (6350779162034176*y(2))/29593047598335265 + 229/5)); eq3 = 0; eq4 = (30*((1192349413690968375975664624440915812941824*y(2)^2)/1036642641726526473848753494572130715682924789385 - (39304596247310823728047194112*y(2))/175149693231467319161999868524045 + 1266637395197952/29593047598335265)*((1266637395197952*y(2))/29593047598335265 - (19652298123655411864023597056*y(2)^2)/175149693231467319161999868524045 + (397449804563656125325221541480305270980608*y(2)^3)/1036642641726526473848753494572130715682924789385 + (5520653160719109*y(4)*((4332790137498830962146934784000*y(2)^2)/3184539876935769439309088518619 - (3982870920455782400*y(2))/5918609519667053 + 288000))/731834939447705600000 + 33/2))/((574685827824708321884380135181365943857027170818850816*y(2)^4)/557771182528674706092576514838123659160733159500324835031693855 - (1050791949051857975174900787749300236976128*y(2)^3)/1036642641726526473848753494572130715682924789385 + (119770698800860541596490268672*y(2)^2)/175149693231467319161999868524045 - (6350779162034176*y(2))/29593047598335265 + 229/5) + (16561959482157327*y(4)*(300*y(2)^2*(36766865660871897387/(96970498370224996352000000*(9/10 - (17592186044416*y(2))/5918609519667053)^(37/10)) - (33800000000000*(8500259669165361/(9444732965739290427392*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)) + 13/250))/(250119*y(2)^5*cos(log((- 74530005132491619276355244274867242431640625*y(4)^2*pi^2 - 8612356148643476394993094470262377676800000000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5) - 248801399849700440447422519369183536317307289600000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5) - 12261830142026331516929221358481899097130334070767616*y(2)^2*y(3)^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5))^(1/2)/(y(4)*pi*8633076226496069765625i + y(4)*pi*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)*498799959753106275696640000i - 110733148343331824175415296*y(2)*y(3)*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)))*1i)^2) + 35851247107254511943359375/(86237027846621531955789824*y(2)^4*cos(log((- 74530005132491619276355244274867242431640625*y(4)^2*pi^2 - 8612356148643476394993094470262377676800000000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5) - 248801399849700440447422519369183536317307289600000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5) - 12261830142026331516929221358481899097130334070767616*y(2)^2*y(3)^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5))^(1/2)/(y(4)*pi*8633076226496069765625i + y(4)*pi*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)*498799959753106275696640000i - 110733148343331824175415296*y(2)*y(3)*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)))*1i)^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(28/5))) + 600*y(2)*(1361735765217477681/(28823037615171174400000*(9/10 - (17592186044416*y(2))/5918609519667053)^(27/10)) + (8450000000000*(8500259669165361/(9444732965739290427392*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)) + 13/250))/(250119*y(2)^4*cos(log((- 74530005132491619276355244274867242431640625*y(4)^2*pi^2 - 8612356148643476394993094470262377676800000000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5) - 248801399849700440447422519369183536317307289600000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5) - 12261830142026331516929221358481899097130334070767616*y(2)^2*y(3)^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5))^(1/2)/(y(4)*pi*8633076226496069765625i + y(4)*pi*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)*498799959753106275696640000i - 110733148343331824175415296*y(2)*y(3)*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)))*1i)^2) + 361283/50000000) + (((8665580274997661924293869568000*y(2))/3184539876935769439309088518619 - 3982870920455782400/5918609519667053)*((1266637395197952*y(2))/29593047598335265 - (19652298123655411864023597056*y(2)^2)/175149693231467319161999868524045 + (397449804563656125325221541480305270980608*y(2)^3)/1036642641726526473848753494572130715682924789385 + (5520653160719109*y(4)*((4332790137498830962146934784000*y(2)^2)/3184539876935769439309088518619 - (3982870920455782400*y(2))/5918609519667053 + 288000))/731834939447705600000 + 33/2))/((574685827824708321884380135181365943857027170818850816*y(2)^4)/557771182528674706092576514838123659160733159500324835031693855 - (1050791949051857975174900787749300236976128*y(2)^3)/1036642641726526473848753494572130715682924789385 + (119770698800860541596490268672*y(2)^2)/175149693231467319161999868524045 - (6350779162034176*y(2))/29593047598335265 + 229/5)))/73183493944770560000 - (30*((1149371655649416643768760270362731887714054341637701632*y(2)^3)/557771182528674706092576514838123659160733159500324835031693855 - (1576187923577786962762351181623950355464192*y(2)^2)/1036642641726526473848753494572130715682924789385 + (119770698800860541596490268672*y(2))/175149693231467319161999868524045 - 3175389581017088/29593047598335265)*((1266637395197952*y(2))/29593047598335265 - (19652298123655411864023597056*y(2)^2)/175149693231467319161999868524045 + (397449804563656125325221541480305270980608*y(2)^3)/1036642641726526473848753494572130715682924789385 + (5520653160719109*y(4)*((4332790137498830962146934784000*y(2)^2)/3184539876935769439309088518619 - (3982870920455782400*y(2))/5918609519667053 + 288000))/731834939447705600000 + 33/2)^2)/((574685827824708321884380135181365943857027170818850816*y(2)^4)/557771182528674706092576514838123659160733159500324835031693855 - (1050791949051857975174900787749300236976128*y(2)^3)/1036642641726526473848753494572130715682924789385 + (119770698800860541596490268672*y(2)^2)/175149693231467319161999868524045 - (6350779162034176*y(2))/29593047598335265 + 229/5)^2 - (82809797410786635*y(3)*tan(log((- 74530005132491619276355244274867242431640625*y(4)^2*pi^2 - 8612356148643476394993094470262377676800000000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5) - 248801399849700440447422519369183536317307289600000000*y(4)^2*pi^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5) - 12261830142026331516929221358481899097130334070767616*y(2)^2*y(3)^2*(9/10 - (17592186044416*y(2))/5918609519667053)^(46/5))^(1/2)/(pi*y(4)*8633076226496069765625i + y(4)*pi*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)*498799959753106275696640000i - 110733148343331824175415296*y(2)*y(3)*(9/10 - (17592186044416*y(2))/5918609519667053)^(23/5)))*1i))/(281474976710656*y(2)^2) - (96559323064259661442131689472*y(2))/35029938646293463832399973704809 - 1636073302130688/5918609519667053; dXdt = zeros(4,1); dXdt(1) = eq1*t_f; dXdt(2) = eq2*t_f; dXdt(3) = eq3*t_f; dXdt(4) = eq4*t_f; end

Thanks in advance!

[Torsten]

Don't you get a warning that the solver cannot continue integration ?

[Sam]

No I don't, strange enough

[Torsten]

I get this message. So you should again check your equations and take a look at the results up to the time when the solver gives up.

[Sam]

Will do, thanks for your effort.

[Jan]

The posted equation contains meaningless information. Remember that

12261830142026331516929221358481899097130334070767616

is rounded to a double internally, so it is the same as:

1.22618301420263e+52

If you are working numerically with the standard double precision, the equation can be simplified remarkably. Nevertheless, this is not the source of the problem.

[Sam]

Yes I know. The reason why the equations currently look like this is because they are the result of a series of (quite extensive) symbolic MATLAB operations.

[Jan]

I suggest to solve the actual problem of the integration interval at first. But then a massive simplification of the equation should be the next step, if runtime matters.

[Are Mjaavatten]

Your system is numerically highly unstable, and the solution depends heavily on the integration method and accuracy parameters. Using ode15s, which is probably more robust than ode45 in this case, the solution seems to have a near singularity at around t = 4.4e-6, where y(4) gets very close to zero.

Try

[t,u] = ode15s(@deWringer,[0,1e-5],y0);
plot(t,u,'.')

(I saved your odefun as deWringer.m.)

Most likely, there is something wrong in your derivation of odefun.