I got two different graphs from my code like ODE89 and Eulers method. I need to compare the graphs in it.
2 次查看(过去 30 天)
显示 更早的评论
ODE89
clc
clear all
A0=1;
B0=3;
P0=0;
K=5*10^-5;
Yb=1;
Yp=0.15;
tspan = [0 43200];
[t,Y] = ode89(@(t,Y) odefun(t,Y, K, Yb, Yp),tspan, [A0;B0;P0]);
figure (1)
plot(t,Y(:,1))
figure (2)
plot(t,Y(:,2))
figure (3)
plot(t,Y(:,3))
function dYdt = odefun(t,Y,K,Yb,Yp)
dYdt = [(-K*Y(1)*Y(2));
(-Yb*(K*Y(1)*Y(2)));
(Yp*(K*Y(1)*Y(2)))];
end
Eulers
nsteps = 12;
t = zeros (nsteps,1);
A = zeros (nsteps,1);
B = zeros(nsteps, 1);
P = zeros(nsteps,1);
A(1) = 1;
B(1) = 3;
C(1) = 0;
K = 5*10^-5
for k = 2:13
t(k) = t(k-1)+3600
A(k) = A(k-1)+(-K*A(k-1)*B(k-1))*3600;
B(k) = B(k-1)+(-Yb*(K*A(k-1)*B(k-1)))*3600;
P(k)= P(k-1)+ Yp*(K*A(k-1)*B(k-1))*3600;
end
plot (t,A)
figure (1)
plot(t,A(:,1))
plot (t,B)
figure (2)
plot(t,B(:,1))
plot(t,P)
figure (3)
plot(t,P(:,1))
I need to compare the ODE89 graph of A,B and P with Euler's A,B and P
2 个评论
AndresVar
2022-3-14
They look pretty similar, so to make it easier to draw comparisons you can plot the corresponding solutions on the same axis.
Maybe you can plot the percent differce on the right-hand side axis also.
采纳的回答
Davide Masiello
2022-3-15
编辑:Davide Masiello
2022-3-15
clc
clear all
%% ODE89
A0=1;
B0=3;
P0=0;
K=5*10^-5;
Yb=1;
Yp=0.15;
tspan = [0 43200];
[x,Y] = ode89(@(t,Y) odefun(t,Y, K, Yb, Yp),tspan, [A0;B0;P0]);
%% Euler
nsteps = 12;
t = zeros (nsteps,1);
A = zeros (nsteps,1);
B = zeros(nsteps, 1);
P = zeros(nsteps,1);
A(1) = 1;
B(1) = 3;
C(1) = 0;
K = 5*10^-5;
for k = 2:13
t(k) = t(k-1)+3600;
A(k) = A(k-1)+(-K*A(k-1)*B(k-1))*3600;
B(k) = B(k-1)+(-Yb*(K*A(k-1)*B(k-1)))*3600;
P(k)= P(k-1)+ Yp*(K*A(k-1)*B(k-1))*3600;
end
figure(1)
subplot(1,3,1)
plot (t,A,x,Y(:,1))
xlabel('t')
ylabel('A')
subplot(1,3,2)
plot (t,B,x,Y(:,2))
xlabel('t')
ylabel('B')
subplot(1,3,3)
plot(t,C,x,Y(:,3))
xlabel('t')
ylabel('P')
legend('ODE89','Euler','Location','Best')
function dYdt = odefun(t,Y,K,Yb,Yp)
dYdt = [(-K*Y(1)*Y(2));
(-Yb*(K*Y(1)*Y(2)));
(Yp*(K*Y(1)*Y(2)))];
end
更多回答(0 个)
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!