Info

此问题已关闭。 请重新打开它进行编辑或回答。

Numerical solution is off for system of ODEs

1 次查看(过去 30 天)
Thanks in advance for your help. I'm not looking for an explicit solution to my issue, but rather to be pointed in the right direction.
I have been plugging away at solving a system of non-linear, first order ODEs in MATLAB. The system was solved numerically in this study: http://web.math.ku.dk/~moller/e04/bio/ludwig78.pdf
I have done all of the work to understand and recreate the model from scratch. I presented the qualitative part for a class project. What I am doing now is taking that project a step farther and trying to first solve the system in MATLAB with runge-kutta (or any method that works). Finally, I want to dive into the theory behind the numerical analysis to find out why the chosed method converges.
ludwig.jpg
I have found that I can create a plot with roughly the same shape, but there are several problems:
1.) The time-scale over which the change occurs is three times that of the above plot.
2.) The range of function values is is vastly wrong.
3.) The desired shapes only occur if I tweak the initial conditions to be significantly different than what is shown near t=0 above.
So, finally, my question is simply to ask what might be causing the disparity I am seeing in my plots?
Code thus far:
% System Parameters:
r_b = 1.52;
k_b = 355;
alph = 1.11;
bet = 43200;
r_e = 0.92;
k_e = 1;
p = 0.00195;
r_s = 0.095;
k_s = 25440;
tspan = [0 200];
init = [1 1 1];
[t, Y] = ode45(@(t,y) odefcn(t, y, r_b, k_b, alph, bet, r_e, k_e, p, r_s, k_s), tspan, init);
subplot(3,1,1);
plot(t,Y(:,1),'b');
title('Budworm Density');
subplot(3,1,2)
plot(t,Y(:,2),'g');
title('Branch Density');
subplot(3,1,3);
plot(t,Y(:,3),'r');
title('Foliage Condition');
function dydt = odefcn(t, y, r_b, k_b, alph, bet, r_e, k_e, p, r_s, k_s)
dydt = [ r_b*y(1)*(1 - y(1)/(k_b*y(2))) - bet*(y(1)^2/((alph*y(2))^2 + y(1)^2));
r_s*y(2)*(1 - (y(2)*k_e)/(k_s*y(3)));
r_e*y(3)*(1 - (y(3)/k_e)) - p*y(1)/y(2)
];
end

回答(1 个)

SandeepKumar R
SandeepKumar R 2019-3-6
To be consistent use years also in your ode function (better to stick to year). Secondly play with tolerances of ode( Reltol,AbsTol). If that doesn't work change the ode solver (23,15s..etc..)

此问题已关闭。

标签

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by