ODE solving using RK4 method (Predator prey)
显示 更早的评论
function Q3
clc
clear all
close all
y = zeros(1,1);
x = zeros(1,1);
%Initial conditions
t = 0.0;
x(1) = 2;
y(1) = 1;
%Parameters
h = 0.001;
tfinal = 30.0;
nsteps=(tfinal-t)/h;
%store intial values
tout(1) = t;
xout(1) = x(1);
yout(1) = y(1);
for i=2:nsteps
k1 = getf(t,y,x);
k2 = getf(t+h/2,y+h/2*k1,x+h/2*k1);
k3 = getf(t+h/2,y+h/2*k2, x+h/2*k2);
k4 = getf(t+h,y+h*k3, x+h/2*k3);
y = y+h/6*(k1+2*k2 + 2*k3 + k4);
x = x+h/6*(k1+2*k2 + 2*k3 + k4);
t = t+h;
xout(i) = x(1);
yout(i) = y(1);
tout(i) = t;
end
plot(tout,yout,"x",tout,xout,"-r")
function f = getf(t,y,x)
alpha = 1.2;
beta = 0.6;
gamma = 0.8;
delta = 0.3;
f(1) = alpha*x(1)-beta*x(1)*y(1);
f(2) = -gamma*x(1)+delta*y(1)*x(1);
There is something wrong with this because the solution isnt right
5 个评论
John D'Errico
2020-4-3
Please do not remove your question, just saying the solution is not right. In fact, it is far more likely tjhat you just do not understand what was said. Of course, nobody can help you now, since we cannot see your question. This insults the person who tried to help you, and it makes your question useless for others who might try to understnd the problem later.
Diya Saha
2020-4-3
Stephen23
2020-4-3
Original question from Google Cache:
"ODE solving using RK4 method (Predator prey)"
function Q3
clc
clear all
close all
y = zeros(1,1);
x = zeros(1,1);
%Initial conditions
t = 0.0;
x(1) = 2;
y(1) = 1;
%Parameters
h = 0.001;
tfinal = 30.0;
nsteps=(tfinal-t)/h;
%store intial values
tout(1) = t;
xout(1) = x(1);
yout(1) = y(1);
for i=2:nsteps
k1 = getf(t,y,x);
k2 = getf(t+h/2,y+h/2*k1,x+h/2*k1);
k3 = getf(t+h/2,y+h/2*k2, x+h/2*k2);
k4 = getf(t+h,y+h*k3, x+h/2*k3);
y = y+h/6*(k1+2*k2 + 2*k3 + k4);
x = x+h/6*(k1+2*k2 + 2*k3 + k4);
t = t+h;
xout(i) = x(1);
yout(i) = y(1);
tout(i) = t;
end
plot(tout,yout,"x",tout,xout,"-r")
function f = getf(t,y,x)
alpha = 1.2;
beta = 0.6;
gamma = 0.8;
delta = 0.3;
f(1) = alpha*x(1)-beta*x(1)*y(1);
f(2) = -gamma*x(1)+delta*y(1)*x(1);
There is something wrong with this because the solution isnt right
Diya Saha
2020-4-3
Rena Berman
2020-5-14
(Answers Dev) Restored edit
采纳的回答
James Tursa
2020-4-3
Your basic problem is that you have two states, x and y, but your function arguments are inconsistent with this. Take this code:
k1 = getf(t,y,x);
getf( ) returns a 2-element vector. So k1 is a 2-element vector. Yet you turn around and treat it like it is a scalar in the very next line:
k2 = getf(t+h/2,y+h/2*k1,x+h/2*k1);
I.e., you end up passing vectors into getf in those 2nd and 3rd input aguments when they should be scalars.
It is unclear to me which element is supposed to be x and which one is y. Since y is first in your argument list, let's assume f(1) is y and f(2) is x. Then the calculation for k2 should be something like this instead:
k2 = getf(t+h/2,y+h/2*k1(1),x+h/2*k1(2));
So all of your function calls will have to be fixed for this.
Having said all of that, I would advise carrying one state vector variable that contains both y and x, instead of carrying two separate y and x variables. That would simplify your code and also put your functions in a form that could be used by calls to the MATLAB supplied integrators such as ode45( ).
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Programming 的更多信息
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
Translated by 
