Info

This question is closed. Reopen it to edit or answer.

why plots are not matching with ode

1 view (last 30 days)
shiv gaur
shiv gaur on 12 Mar 2022
Closed: Jan on 12 Mar 2022
This question was flagged by Jan
% Set parameters
Tstart = 0.0;
Tend = 200.0;
Nt = 200000;
dT = (Tend-Tstart)/Nt;
X0 = 1;
Y0 = 1;
Z0 = 0;
N = 20;
SIGMA = 0.2;
R = 0.2;
B = 5.7;
%
% Initialize coefficient arrays
T = zeros(Nt+1,1);
X = zeros(Nt+1,1);
Y = zeros(Nt+1,1);
Z = zeros(Nt+1,1);
a = zeros(N+1,1);
b = zeros(N+1,1);
c = zeros(N+1,1);
T(1) = 0.0;
X(1) = X0;
Y(1) = Y0;
Z(1) = Z0;
for j = 2:Nt+1
a(1) = X(j-1);
b(1) = Y(j-1);
c(1) = Z(j-1);
for k = 1:N
SB = 0.0;
for i= 0:k-1
SB = SB + a(i+1)*c(k-i);
end
a(k+1) = ((-b(k) - a(k)))/k+1;
b(k+1) = ((a(k) +R*b(k)))/k+1 ;
c(k+1) =( SIGMA-B*c(k) + SB)/k+1 ;
end
x = a(1);
y = b(1);
z = c(1);
for k = 2:N+1
x = x + a(k)*dT^(k-1);
y = y + b(k)*dT^(k-1);
z = z + c(k)*dT^(k-1);
end
T(j) = T(j-1) + dT;
X(j) = x;
Y(j) = y;
Z(j) = z;
end
figure(1)
plot(T,X,'Color','red')
hold on
figure(2)
plot(T,Y,'Color','red')
hold on
figure(3)
plot(T,Z,'Color','red')
hold on
figure(4)
plot3(X,Y,Z,'Color','red')
hold on
%matching with ode
function kk2
close all; clear all;
% accuracy
%options = odeset('RelTol',1e-9,'AbsTol', 1e-16);
a=0.2;
b=0.2;
k=5.7;
[t,y]=ode45(@f,[0 200],[1;1;0.05] );
%plot(t,y(:,2),t,y(:,4),t,y(:,6))
%
plot(t,y(:,1),t,y(:,2),t,y(:,3))
xlabel('t');
ylabel('x,y,z')
%plot3(y(:,1),y(2,:),y(:,3))
function dy=f(t,y)
a=0.2;
b=0.2;
k=5.7;
dy=[(-y(2)-y(3));y(1)+a*y(2);b+y(3)*(y(1)-k)];
end
end
what is the problem above program
  3 Comments
Jan
Jan on 12 Mar 2022
Your question does not contain enough information to be answered. You state, that the "plots do not match the ODE", but of course the diagrams are the results of the equations. As said already, there is no chance to guess, why you expect different results, because you do not mention any reason.
Re-opening your question without adding any further useful details is a misusage of the forum. Stop this.

Answers (0)

Tags

Products


Release

R2021b

Community Treasure Hunt

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

Start Hunting!

Translated by