Error using horzcat Requested 36x28689400 (7.7GB) array exceeds
17 views (last 30 days)
Show older comments
Manish Kumar
on 22 Aug 2024 at 14:32
Edited: Manish Kumar
on 23 Aug 2024 at 13:48
clear
clc
% Define the time span for the solution
tspan = [0 1]; % For example, from 0 to 10 seconds
% Define the initial conditions for all 36 variables
y0 = 0.1 * ones(36, 1);
% Solve the differential equations using ode45
%options = odeset('RelTol',1e-4,'AbsTol',1e-6);
[t, y] = ode45(@odefun, tspan, y0);
% Plot the results (for example, plot y1 and y2 over time)
figure;
plot(t, y(:,1), '-o', t, y(:,2), '-x');
title('Solutions of the Differential Equations');
xlabel('Time (s)');
ylabel('Solution values');
legend('y1', 'y2');
grid on;
function dydt = odefun(t,y)
dydt = zeros(36,1);
dydt(1) = 26300/693 - (10000*y(2))/693 - (100*y(1))/693;
dydt(2) = (5000*y(1))/1731 - (5000*y(4))/(1731*y(2));
dydt(3) = 3141/50;
dydt(4) = (3141*y(14)*y(12))/1000 - (3141*y(4))/1000 + (3141*y(15)*y(13))/1000;
dydt(5) = (3141*y(15)*y(12))/1000 - (3141*y(14)*y(13))/1000 - (3141*y(5))/1000;
dydt(6) = 282215255899405/281474976710656 - y(12) - (171*y(5))/5000;
dydt(7) = -y(13);
dydt(8) = (3*y(14))/4 - y(10) - (171*y(5))/4000 + (39*y(6))/4000 - (5*y(12))/4 - (8019065281956203*y(13))/73786976294838206464 + 1411076279497025/1125899906842624;
dydt(9) = (3*y(15))/4 - y(11) + (39*y(7))/4000 + (8019065281956203*y(12))/73786976294838206464 - (5*y(13))/4;
dydt(10) = (6189965296924617*y(2))/2199023255552 - (5278537774523551167*y(5)*y(2))/54975581388800000 + (7349676656256685*y(8)*y(2))/2147483648 + (92605925868834231*y(14)*y(2))/54975581388800 + (1203877036294845003*y(6)*y(2))/54975581388800000;
dydt(11) = (7349676656256685*y(9)*y(2))/2147483648 + (92605925868834231*y(15)*y(2))/54975581388800 + (1203877036294845003*y(7)*y(2))/54975581388800000 + (4950753100254522323750085935885*y(2)*y(12))/20282409603651670423947251286016 + y(10)*((9954226338292172243339743750363*y(2))/633825300114114700748351602688 - 3141/100) - y(11)*((30868641956278077*y(2))/13743895347200 + 13800/187) - y(13)*((30868641956278077*y(2))/10995116277760 + 1469935331251337/137438953472);
dydt(12) = (4965280110982659*y(10))/17179869184 - (4965280110982659*y(14))/17179869184 + (3141*y(13))/100;
dydt(13) = (4965280110982659*y(11))/17179869184 - (4965280110982659*y(15))/17179869184 - (3141*y(12))/100;
dydt(14) = (3141*y(15))/100 + (2571957024037427*y(12))/2199023255552 - y(14)*((638718513552193224704*cos(y(3))^2)/7885983091510833 + (638718513552193224704*sin(y(3))^2)/7885983091510833 + 140/57) + (638718513552193224704*y(31)*cos(y(3)))/7885983091510833 + (638718513552193224704*y(32)*sin(y(3)))/7885983091510833;
dydt(15) = (2571957024037427*y(13))/2199023255552 - (3141*y(14))/100 - y(15)*((638718513552193224704*cos(y(3))^2)/7885983091510833 + (638718513552193224704*sin(y(3))^2)/7885983091510833 + 140/57) + (638718513552193224704*y(32)*cos(y(3)))/7885983091510833 - (638718513552193224704*y(31)*sin(y(3)))/7885983091510833;
dydt(16) = 26300/693 - (10000*y(17))/693 - (100*y(16))/693;
dydt(17) = (5000*y(16))/1731 - (5000*y(19))/(1731*y(17));
dydt(18) = (3*y(4))/10 - (3*y(19))/10 + 3141/50;
dydt(19) = (3141*y(29)*y(27))/1000 - (3141*y(19))/1000 + (3141*y(30)*y(28))/1000;
dydt(20) = (3141*y(30)*y(27))/1000 - (3141*y(29)*y(28))/1000 - (3141*y(20))/1000;
dydt(21) = 282215255899405/281474976710656 - y(27) - (171*y(20))/5000;
dydt(22) = -y(28);
dydt(23) = (3*y(29))/4 - y(25) - (171*y(20))/4000 + (39*y(21))/4000 - (5*y(27))/4 - (8019065281956203*y(28))/73786976294838206464 + 1411076279497025/1125899906842624;
dydt(24) = (3*y(30))/4 - y(26) + (39*y(22))/4000 + (8019065281956203*y(27))/73786976294838206464 - (5*y(28))/4;
dydt(25) = (6189965296924617*y(17))/2199023255552 - (5278537774523551167*y(20)*y(17))/54975581388800000 + (7349676656256685*y(23)*y(17))/2147483648 + (92605925868834231*y(29)*y(17))/54975581388800 + (1203877036294845003*y(21)*y(17))/54975581388800000 - (4950753100254522323750085935885*y(17)*y(28))/20282409603651670423947251286016 - y(26)*((9954226338292172243339743750363*y(17))/633825300114114700748351602688 - 3141/100) - y(25)*((30868641956278077*y(17))/13743895347200 + 13800/187) - y(27)*((30868641956278077*y(17))/10995116277760 + 1469935331251337/137438953472);
dydt(26) = (7349676656256685*y(24)*y(17))/2147483648 + (92605925868834231*y(30)*y(17))/54975581388800 + (1203877036294845003*y(22)*y(17))/54975581388800000 + (4950753100254522323750085935885*y(17)*y(27))/20282409603651670423947251286016 + y(25)*((9954226338292172243339743750363*y(17))/633825300114114700748351602688 - 3141/100) - y(26)*((30868641956278077*y(17))/13743895347200 + 13800/187) - y(28)*((30868641956278077*y(17))/10995116277760 + 1469935331251337/137438953472);
dydt(27) = (4965280110982659*y(25))/17179869184 - (4965280110982659*y(29))/17179869184 + (3141*y(28))/100;
dydt(28) = (4965280110982659*y(26))/17179869184 - (4965280110982659*y(30))/17179869184 - (3141*y(27))/100;
dydt(29) = (3141*y(30))/100 - (140*y(29))/57 + (2571957024037427*y(27))/2199023255552 - y(14)*((638718513552193224704*sin(y(3))*sin(y(18)))/7885983091510833 + (638718513552193224704*cos(y(3))*cos(y(18)))/7885983091510833) - y(15)*((638718513552193224704*cos(y(3))*sin(y(18)))/7885983091510833 - (638718513552193224704*cos(y(18))*sin(y(3)))/7885983091510833) + (638718513552193224704*y(31)*cos(y(18)))/7885983091510833 + (638718513552193224704*y(32)*sin(y(18)))/7885983091510833;
dydt(30) = (2571957024037427*y(28))/2199023255552 - (140*y(30))/57 - (3141*y(29))/100 + y(14)*((638718513552193224704*cos(y(3))*sin(y(18)))/7885983091510833 - (638718513552193224704*cos(y(18))*sin(y(3)))/7885983091510833) - y(15)*((638718513552193224704*sin(y(3))*sin(y(18)))/7885983091510833 + (638718513552193224704*cos(y(3))*cos(y(18)))/7885983091510833) + (638718513552193224704*y(32)*cos(y(18)))/7885983091510833 - (638718513552193224704*y(31)*sin(y(18)))/7885983091510833;
dydt(31) = (2758738413894493*y(33))/274877906944 - (5518743656403597*y(31))/549755813888 + (3141*y(32))/100 + (2758738413894493*y(14)*cos(y(3)))/274877906944 - (2758738413894493*y(29)*cos(y(18)))/274877906944 - (2758738413894493*y(15)*sin(y(3)))/274877906944 + (2758738413894493*y(30)*sin(y(18)))/274877906944;
dydt(32) = (2758738413894493*y(34))/274877906944 - (5518743656403597*y(32))/549755813888 - (3141*y(31))/100 + (2758738413894493*y(15)*cos(y(3)))/274877906944 - (2758738413894493*y(30)*cos(y(18)))/274877906944 + (2758738413894493*y(14)*sin(y(3)))/274877906944 - (2758738413894493*y(29)*sin(y(18)))/274877906944;
dydt(33) = (3788118419079005*y(31))/2199023255552 - (3788118419079005*y(35))/2199023255552 + (3786794628959991*y(33))/2199023255552 + (3788118419079005*y(32))/2199023255552 + (3141*y(34))/100 - (3788118419079005*y(15)*cos(y(3)))/2199023255552 - (3788118419079005*y(14)*sin(y(3)))/2199023255552;
dydt(34) = (3788118419079005*y(32))/2199023255552 - (3788118419079005*y(36))/2199023255552 - (3141*y(33))/100 + (1893728262009749*y(34))/549755813888 - (3788118419079005*y(30)*cos(y(18)))/2199023255552 - (3788118419079005*y(29)*sin(y(18)))/2199023255552;
dydt(35) = (277*y(31))/6 - (6659819921464029*y(35))/140737488355328 + (277*y(33))/6 + (3141*y(36))/100;
dydt(36) = (277*y(32))/6 - (6659819921464029*y(36))/140737488355328 - (3141*y(35))/100 + (277*y(34))/6;
%dydt=[dydt(1);dydt(2);dydt(3);dydt(4);dydt(5);dydt(6);dydt(7);dydt(8);dydt(9);dydt(10);dydt(11);dydt(12);dydt(13);dydt(14);dydt(15);dydt(16);dydt(17);dydt(18);dydt(19);dydt(20);dydt(21);dydt(22);dydt(23);dydt(24);dydt(24);dydt(25);dydt(26);dydt(27);dydt(28);dydt(29);dydt(30);dydt(31);dydt(32);dydt(33);dydt(34);dydt(34);dydt(35);dydt(36)];
end
error
Error using horzcat
Requested 36x28689400 (7.7GB) array exceeds
maximum array size preference (7.7GB). This might
cause MATLAB to become unresponsive.
Error in ode45 (line 476)
yout = [yout, zeros(neq,chunk,dataType)]; %#ok<AGROW>
Error in sol_ode (line 11)
[t, y] = ode45(@odefun, tspan, y0);
Related documentation
0 Comments
Accepted Answer
Steven Lord
on 22 Aug 2024 at 17:39
That suggests to me that your system of ODEs is stiff and so ode45 requires a lot of time steps. Try a stiffer solver like ode23s.
tspan = [0 1]; % For example, from 0 to 10 seconds
% Define the initial conditions for all 36 variables
y0 = 0.1 * ones(36, 1);
% Solve the differential equations using ode45
%options = odeset('RelTol',1e-4,'AbsTol',1e-6);
[t, y] = ode23s(@odefun, tspan, y0);
whos t y
So even with the stiffer solver, it seems like there's a singularity or some sort of problem well before you reach the end of the time span. I've moved the location of the legend on the following plot to the west side of the plot because for small values of t the function looks pretty flat. The "interesting" behavior take place at later values of t.
plot(t, y)
legend(Location="west")
So one of the components (from the color it looks like probably the 5th, 12th, or 19th component) plummets at the time where the solver reports it can't proceed.
subplot(1, 3, 1)
plot(t, y(:, 5))
title('y(:, 5)')
subplot(1, 3, 2)
plot(t, y(:, 12))
title('y(:, 12)')
subplot(1, 3, 3)
plot(t, y(:, 19))
title('y(:, 19)')
Looks like y(:, 12) is the likely culprit, though y(:, 5) and y(:, 19) are also showing signs of wildly oscillatory behavior. I'd make sure those equations in your odefun function are doing what you expect them to do and what you want them to do.
function dydt = odefun(t,y)
dydt = zeros(36,1);
dydt(1) = 26300/693 - (10000*y(2))/693 - (100*y(1))/693;
dydt(2) = (5000*y(1))/1731 - (5000*y(4))/(1731*y(2));
dydt(3) = 3141/50;
dydt(4) = (3141*y(14)*y(12))/1000 - (3141*y(4))/1000 + (3141*y(15)*y(13))/1000;
dydt(5) = (3141*y(15)*y(12))/1000 - (3141*y(14)*y(13))/1000 - (3141*y(5))/1000;
dydt(6) = 282215255899405/281474976710656 - y(12) - (171*y(5))/5000;
dydt(7) = -y(13);
dydt(8) = (3*y(14))/4 - y(10) - (171*y(5))/4000 + (39*y(6))/4000 - (5*y(12))/4 - (8019065281956203*y(13))/73786976294838206464 + 1411076279497025/1125899906842624;
dydt(9) = (3*y(15))/4 - y(11) + (39*y(7))/4000 + (8019065281956203*y(12))/73786976294838206464 - (5*y(13))/4;
dydt(10) = (6189965296924617*y(2))/2199023255552 - (5278537774523551167*y(5)*y(2))/54975581388800000 + (7349676656256685*y(8)*y(2))/2147483648 + (92605925868834231*y(14)*y(2))/54975581388800 + (1203877036294845003*y(6)*y(2))/54975581388800000;
dydt(11) = (7349676656256685*y(9)*y(2))/2147483648 + (92605925868834231*y(15)*y(2))/54975581388800 + (1203877036294845003*y(7)*y(2))/54975581388800000 + (4950753100254522323750085935885*y(2)*y(12))/20282409603651670423947251286016 + y(10)*((9954226338292172243339743750363*y(2))/633825300114114700748351602688 - 3141/100) - y(11)*((30868641956278077*y(2))/13743895347200 + 13800/187) - y(13)*((30868641956278077*y(2))/10995116277760 + 1469935331251337/137438953472);
dydt(12) = (4965280110982659*y(10))/17179869184 - (4965280110982659*y(14))/17179869184 + (3141*y(13))/100;
dydt(13) = (4965280110982659*y(11))/17179869184 - (4965280110982659*y(15))/17179869184 - (3141*y(12))/100;
dydt(14) = (3141*y(15))/100 + (2571957024037427*y(12))/2199023255552 - y(14)*((638718513552193224704*cos(y(3))^2)/7885983091510833 + (638718513552193224704*sin(y(3))^2)/7885983091510833 + 140/57) + (638718513552193224704*y(31)*cos(y(3)))/7885983091510833 + (638718513552193224704*y(32)*sin(y(3)))/7885983091510833;
dydt(15) = (2571957024037427*y(13))/2199023255552 - (3141*y(14))/100 - y(15)*((638718513552193224704*cos(y(3))^2)/7885983091510833 + (638718513552193224704*sin(y(3))^2)/7885983091510833 + 140/57) + (638718513552193224704*y(32)*cos(y(3)))/7885983091510833 - (638718513552193224704*y(31)*sin(y(3)))/7885983091510833;
dydt(16) = 26300/693 - (10000*y(17))/693 - (100*y(16))/693;
dydt(17) = (5000*y(16))/1731 - (5000*y(19))/(1731*y(17));
dydt(18) = (3*y(4))/10 - (3*y(19))/10 + 3141/50;
dydt(19) = (3141*y(29)*y(27))/1000 - (3141*y(19))/1000 + (3141*y(30)*y(28))/1000;
dydt(20) = (3141*y(30)*y(27))/1000 - (3141*y(29)*y(28))/1000 - (3141*y(20))/1000;
dydt(21) = 282215255899405/281474976710656 - y(27) - (171*y(20))/5000;
dydt(22) = -y(28);
dydt(23) = (3*y(29))/4 - y(25) - (171*y(20))/4000 + (39*y(21))/4000 - (5*y(27))/4 - (8019065281956203*y(28))/73786976294838206464 + 1411076279497025/1125899906842624;
dydt(24) = (3*y(30))/4 - y(26) + (39*y(22))/4000 + (8019065281956203*y(27))/73786976294838206464 - (5*y(28))/4;
dydt(25) = (6189965296924617*y(17))/2199023255552 - (5278537774523551167*y(20)*y(17))/54975581388800000 + (7349676656256685*y(23)*y(17))/2147483648 + (92605925868834231*y(29)*y(17))/54975581388800 + (1203877036294845003*y(21)*y(17))/54975581388800000 - (4950753100254522323750085935885*y(17)*y(28))/20282409603651670423947251286016 - y(26)*((9954226338292172243339743750363*y(17))/633825300114114700748351602688 - 3141/100) - y(25)*((30868641956278077*y(17))/13743895347200 + 13800/187) - y(27)*((30868641956278077*y(17))/10995116277760 + 1469935331251337/137438953472);
dydt(26) = (7349676656256685*y(24)*y(17))/2147483648 + (92605925868834231*y(30)*y(17))/54975581388800 + (1203877036294845003*y(22)*y(17))/54975581388800000 + (4950753100254522323750085935885*y(17)*y(27))/20282409603651670423947251286016 + y(25)*((9954226338292172243339743750363*y(17))/633825300114114700748351602688 - 3141/100) - y(26)*((30868641956278077*y(17))/13743895347200 + 13800/187) - y(28)*((30868641956278077*y(17))/10995116277760 + 1469935331251337/137438953472);
dydt(27) = (4965280110982659*y(25))/17179869184 - (4965280110982659*y(29))/17179869184 + (3141*y(28))/100;
dydt(28) = (4965280110982659*y(26))/17179869184 - (4965280110982659*y(30))/17179869184 - (3141*y(27))/100;
dydt(29) = (3141*y(30))/100 - (140*y(29))/57 + (2571957024037427*y(27))/2199023255552 - y(14)*((638718513552193224704*sin(y(3))*sin(y(18)))/7885983091510833 + (638718513552193224704*cos(y(3))*cos(y(18)))/7885983091510833) - y(15)*((638718513552193224704*cos(y(3))*sin(y(18)))/7885983091510833 - (638718513552193224704*cos(y(18))*sin(y(3)))/7885983091510833) + (638718513552193224704*y(31)*cos(y(18)))/7885983091510833 + (638718513552193224704*y(32)*sin(y(18)))/7885983091510833;
dydt(30) = (2571957024037427*y(28))/2199023255552 - (140*y(30))/57 - (3141*y(29))/100 + y(14)*((638718513552193224704*cos(y(3))*sin(y(18)))/7885983091510833 - (638718513552193224704*cos(y(18))*sin(y(3)))/7885983091510833) - y(15)*((638718513552193224704*sin(y(3))*sin(y(18)))/7885983091510833 + (638718513552193224704*cos(y(3))*cos(y(18)))/7885983091510833) + (638718513552193224704*y(32)*cos(y(18)))/7885983091510833 - (638718513552193224704*y(31)*sin(y(18)))/7885983091510833;
dydt(31) = (2758738413894493*y(33))/274877906944 - (5518743656403597*y(31))/549755813888 + (3141*y(32))/100 + (2758738413894493*y(14)*cos(y(3)))/274877906944 - (2758738413894493*y(29)*cos(y(18)))/274877906944 - (2758738413894493*y(15)*sin(y(3)))/274877906944 + (2758738413894493*y(30)*sin(y(18)))/274877906944;
dydt(32) = (2758738413894493*y(34))/274877906944 - (5518743656403597*y(32))/549755813888 - (3141*y(31))/100 + (2758738413894493*y(15)*cos(y(3)))/274877906944 - (2758738413894493*y(30)*cos(y(18)))/274877906944 + (2758738413894493*y(14)*sin(y(3)))/274877906944 - (2758738413894493*y(29)*sin(y(18)))/274877906944;
dydt(33) = (3788118419079005*y(31))/2199023255552 - (3788118419079005*y(35))/2199023255552 + (3786794628959991*y(33))/2199023255552 + (3788118419079005*y(32))/2199023255552 + (3141*y(34))/100 - (3788118419079005*y(15)*cos(y(3)))/2199023255552 - (3788118419079005*y(14)*sin(y(3)))/2199023255552;
dydt(34) = (3788118419079005*y(32))/2199023255552 - (3788118419079005*y(36))/2199023255552 - (3141*y(33))/100 + (1893728262009749*y(34))/549755813888 - (3788118419079005*y(30)*cos(y(18)))/2199023255552 - (3788118419079005*y(29)*sin(y(18)))/2199023255552;
dydt(35) = (277*y(31))/6 - (6659819921464029*y(35))/140737488355328 + (277*y(33))/6 + (3141*y(36))/100;
dydt(36) = (277*y(32))/6 - (6659819921464029*y(36))/140737488355328 - (3141*y(35))/100 + (277*y(34))/6;
%dydt=[dydt(1);dydt(2);dydt(3);dydt(4);dydt(5);dydt(6);dydt(7);dydt(8);dydt(9);dydt(10);dydt(11);dydt(12);dydt(13);dydt(14);dydt(15);dydt(16);dydt(17);dydt(18);dydt(19);dydt(20);dydt(21);dydt(22);dydt(23);dydt(24);dydt(24);dydt(25);dydt(26);dydt(27);dydt(28);dydt(29);dydt(30);dydt(31);dydt(32);dydt(33);dydt(34);dydt(34);dydt(35);dydt(36)];
end
1 Comment
More Answers (0)
See Also
Categories
Find more on Ordinary Differential Equations in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!