I got a blank figure while trying to plot the phase space

1 Ansicht (letzte 30 Tage)
I got a blank figure while trying to plot the phase space of this code.
This is the function file
function dxdt = IP(t,x,B,C,D,r,z,A,K)
dxdt = zeros(4,1);
dxdt(1) = x(3);
dxdt(2) = x(4);
dxdt(3) = (z*D*exp(-z*x(1))*(2*A*exp(2*z*r)*exp(-z*x(1)))+(B*exp(z*r))*(1+exp(-(z*x(1))))+2*C)/((1-exp(-(z*x(1))))^3)+K*(x(2)-x(1));
dxdt(4) = (z*D*exp(-z*x(2))*(2*A*exp(2*z*r)*exp(-z*x(2)))+(B*exp(z*r))*(1+exp(-(z*x(2))))+2*C)/((1-exp(-(z*x(2))))^3)+K*(x(1)-x(2));
end
this is the script:
close all
D=0.02;
r=0.04;
A=2.5;
B=0.001;
z=3.0;
C=3.25;
q1=0;
q2=0;
p1=0.632;
E=0.2;
p2=sqrt(2*E- p1^2);
dt=1/5;
tspan=[10:dt:100];
K=[0 0.050 0.5 1.5];
m=length(K);
for i=1:m;
[t x]=ode45(@(t,x)IP(t,x,B,C,D,r,z,A,K(i)),[0 100],[q1 q2 p1 p2]);
figure(1)
subplot(2,2,(i))
plot(x(:,1),x(:,3),'-b',x(:,2),x(:,4),'-r ')
xlabel('q');
ylabel('p');
end

Akzeptierte Antwort

VBBV
VBBV am 18 Mai 2023
Choose better initial conditions e.g. for q1 and q2 they are both zero
close all
D=0.02;
r=0.04;
A=2.5;
B=0.01;
z=3.0;
C=3.25;
%>>>>>>>>> choose hetter intial coditions
q1=1;
q2=0.5;
%-->>>>>>>>>
p1=0.632;
E=0.2;
p2=sqrt(2*E- p1^2);
dt=1/5;
tspan=[10:dt:100];
K=[0 0.050 0.5 1.5];
m=length(K);
figure(1)
hold on
for i=1:m;
[t,x]=ode45(@(t,x)IP(t,x,B,C,D,r,z,A,K(i)),[0 100],[q1 q2 p1 p2]);
subplot(2,2,(i))
plot(x(:,1),x(:,3),'-b',x(:,2),x(:,4),'-r ')
end
xlabel('q');
ylabel('p');
function dxdt = IP(t,x,B,C,D,r,z,A,K)
dxdt = zeros(4,1);
dxdt(1) = x(3);
dxdt(2) = x(4);
dxdt(3) = (z*D*exp(-z*x(1))*(2*A*exp(2*z*r)*exp(-z*x(1)))+(B*exp(z*r))*(1+exp(-(z*x(1))))+2*C)/((1-exp(-(z*x(1))))^3)+K*(x(2)-x(1));
dxdt(4) = (z*D*exp(-z*x(2))*(2*A*exp(2*z*r)*exp(-z*x(2)))+(B*exp(z*r))*(1+exp(-(z*x(2))))+2*C)/((1-exp(-(z*x(2))))^3)+K*(x(1)-x(2));
end
  3 Kommentare
VBBV
VBBV am 18 Mai 2023
yes, xlabel and ylabel commands need to be present inisde the loop
OGUNGBEMI Ezekiel
OGUNGBEMI Ezekiel am 18 Mai 2023
Thanks for the assistance. the result is not producing a phase space plot.

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (0)

Kategorien

Mehr zu Mathematics finden Sie in Help Center und File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by