I couldn't understand the problem in it? and I did't know where I should put the "o" ?
1 vue (au cours des 30 derniers jours)
Afficher commentaires plus anciens
clc
ti = 0;
tf = 1E-3;
tspan=[ti tf];
P = 1.27;
a = 0.1;
tc = 230E-6;
tp = 5.4E-9;
o = rand(1,5)*100E3;
k = 1E-3;
k1 = k;
k2 = k;
k3 = k;
k4 = k;
k5 = k;
f = @(t,y) [
(P - y(1).*((abs(y(6)))^2 +1))./tc;
(P - y(2).*((abs(y(7)))^2 +1))./tc;
(P - y(3).*((abs(y(8)))^2 +1))./tc;
(P - y(4).*((abs(y(9)))^2 +1))./tc;
(P - y(5).*((abs(y(10)))^2 +1))./tc;
(y(1)-a).*((y(6))./tp) + (k1./tp).*(y(7)).*cos(y(11)) + (k5./tp).*(y(10))*cos(y(15));
(y(2)-a).*((y(7))./tp) + (k2./tp).*(y(8)).*cos(y(12)) + (k1./tp).*(y(6))*cos(y(11));
(y(3)-a).*((y(8))./tp) + (k3./tp).*(y(9)).*cos(y(13)) + (k2./tp).*(y(7))*cos(y(12));
(y(4)-a).*((y(9))./tp) + (k4./tp).*(y(10)).*cos(y(14)) + (k3./tp).*(y(8))*cos(y(13));
(y(5)-a).*((y(10))./tp) + (k5./tp).*(y(9)).*cos(y(14)) +(k4./tp).* (y(6))*cos(y(15));
o(1,1) - (k1./tp).*(y(6)/y(7)).*(sin(y(11))) - (k1./tp).*(y(10)/y(6)).*(sin(y(11))) + (k2./tp).*(y(8)./y(7))*sin(y(12)) + (k5./tp).*(y(10)/y(6)).*sin(y(15));
o(1,2) - (k2./tp).*(y(7)/y(8)).*(sin(y(12))) - (k2./tp).*(y(8)/y(7)).*(sin(y(12))) + (k3./tp).*(y(9)./y(8))*sin(y(13)) + (k1./tp).*(y(6)/y(7)).*sin(y(11));
o(1,3) - (k3./tp).*(y(8)/y(9)).*(sin(y(13))) - (k3./tp).*(y(9)/y(8)).*(sin(y(13))) + (k4./tp).*(y(10)./y(9))*sin(y(14)) + (k2./tp).*(y(7)/y(8)).*sin(y(12));
o(1,4) - (k4./tp).*(y(9)/y(10)).*(sin(y(14))) - (k4./tp).*(y(10)/y(9)).*(sin(y(14))) + (k5./tp).*(y(6)./y(10))*sin(y(15)) + (k3./tp).*(y(8)/y(9)).*sin(y(13));
o(1,5) - (k5./tp).*(y(10)/y(6)).*(sin(y(15))) - (k5./tp).*(y(6)/y(10)).*(sin(y(15))) + (k1./tp).*(y(7)./y(6))*sin(y(11)) + (k4./tp).*(y(9)/y(10)).*sin(y(14));
];
[time,Y] = ode45(f,tspan./tp,[0; 0; 0; 0; 0; 1; 1; 1; 1; 1; 0.01; 0.01; 0.01; 0.01; 0.01],o);
subplot 211
plot(time,Y(:,11));
xlabel('time')
ylabel('phase')
subplot 212
plot(time,Y(:,12));
xlabel('time')
ylabel('Amplitude')
1 commentaire
Ankit
le 1 Sep 2022
clc
ti = 0;
tf = 1E-3;
tspan=[ti tf];
P = 1.27;
a = 0.1;
tc = 230E-6;
tp = 5.4E-9;
o = rand(1,5)*100E3;
k = 1E-3;
k1 = k;
k2 = k;
k3 = k;
k4 = k;
k5 = k;
f = @(t,y) [
(P - y(1).*((abs(y(6)))^2 +1))./tc;
(P - y(2).*((abs(y(7)))^2 +1))./tc;
(P - y(3).*((abs(y(8)))^2 +1))./tc;
(P - y(4).*((abs(y(9)))^2 +1))./tc;
(P - y(5).*((abs(y(10)))^2 +1))./tc;
(y(1)-a).*((y(6))./tp) + (k1./tp).*(y(7)).*cos(y(11)) + (k5./tp).*(y(10))*cos(y(15));
(y(2)-a).*((y(7))./tp) + (k2./tp).*(y(8)).*cos(y(12)) + (k1./tp).*(y(6))*cos(y(11));
(y(3)-a).*((y(8))./tp) + (k3./tp).*(y(9)).*cos(y(13)) + (k2./tp).*(y(7))*cos(y(12));
(y(4)-a).*((y(9))./tp) + (k4./tp).*(y(10)).*cos(y(14)) + (k3./tp).*(y(8))*cos(y(13));
(y(5)-a).*((y(10))./tp) + (k5./tp).*(y(9)).*cos(y(14)) + (k4./tp).* (y(6))*cos(y(15));
o(1,1) - (k1./tp).*(y(6)/y(7)).*(sin(y(11))) - (k1./tp).*(y(10)/y(6)).*(sin(y(11))) + (k2./tp).*(y(8)./y(7))*sin(y(12)) + (k5./tp).*(y(10)/y(6)).*sin(y(15));
o(1,2) - (k2./tp).*(y(7)/y(8)).*(sin(y(12))) - (k2./tp).*(y(8)/y(7)).*(sin(y(12))) + (k3./tp).*(y(9)./y(8))*sin(y(13)) + (k1./tp).*(y(6)/y(7)).*sin(y(11));
o(1,3) - (k3./tp).*(y(8)/y(9)).*(sin(y(13))) - (k3./tp).*(y(9)/y(8)).*(sin(y(13))) + (k4./tp).*(y(10)./y(9))*sin(y(14)) + (k2./tp).*(y(7)/y(8)).*sin(y(12));
o(1,4) - (k4./tp).*(y(9)/y(10)).*(sin(y(14))) - (k4./tp).*(y(10)/y(9)).*(sin(y(14))) + (k5./tp).*(y(6)./y(10))*sin(y(15)) + (k3./tp).*(y(8)/y(9)).*sin(y(13));
o(1,5) - (k5./tp).*(y(10)/y(6)).*(sin(y(15))) - (k5./tp).*(y(6)/y(10)).*(sin(y(15))) + (k1./tp).*(y(7)./y(6))*sin(y(11)) + (k4./tp).*(y(9)/y(10)).*sin(y(14));];
[time,Y] = ode45(f,tspan./tp,[0; 0; 0; 0; 0; 1; 1; 1; 1; 1; 0.01; 0.01; 0.01; 0.01; 0.01],o)
subplot 211
plot(time,Y(:,11));
xlabel('time')
ylabel('phase')
subplot 212
plot(time,Y(:,12));
xlabel('time')
ylabel('Amplitude')
Unfortunately your code is taking too much time to run.. I didnt check in details
But your Problem related to dimension is removed
Réponses (1)
Sam Chak
le 1 Sep 2022
Hi @SAHIL SAHOO
I guess that is an error in your ode argument, things that are inside the brackets of ode45().
I believe 'o' stands for options in your case. Try something like this:
o = odeset('RelTol', 1e-2, 'AbsTol', 1e-5);
[time,Y] = ode45(f, tspan./tp, [0; 0; 0; 0; 0; 1; 1; 1; 1; 1; 0.01; 0.01; 0.01; 0.01; 0.01], o);
Also check out this:
0 commentaires
Voir également
Catégories
En savoir plus sur Interactive Control and Callbacks dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!