what's the problem in defining the function in for loop?
1 vue (au cours des 30 derniers jours)
Afficher commentaires plus anciens
SAHIL SAHOO
le 31 Juil 2022
Commenté : Walter Roberson
le 31 Juil 2022
clc
ti = 0; %inital time
tf = 70E-9;% final time
tspan=[ti tf];
tp =1E-12;
T = 2E3;
p = 0.05;
k = (0.62).*10^(4);
c = 3E8;
N = 3.0;
n = k*(c/N)*tp
a = 5;
f = @(t,y) [
(y(4).*y(1) - n.*y(2).*sin(y(3)));
(y(5).*y(2) + n.*y(1).*sin(y(3)));
(-a.*(y(5)-y(4)) + n.*((y(1)./y(2)) - (y(2)./y(1))).*cos(y(3)));
(p - y(4)-(1+2.*y(4)).*(y(1)^2))./T;
(p - y(5)-(1+2.*y(5)).*(y(2)^2))./T;
];
[time,Y] = ode45(f,tspan./tp,[sqrt(p);sqrt(p);0;0;0.01]);
O = linspace(-10,10,100);
U = zeros(length(O),1) ;
for i = 1:length(O)
U(i) = -a.*(Y(:,5) - Y(:,4)).*O(i) + n.*((Y(:,1)./Y(:,2)) - (Y(:,2)./Y(:,1)))*sin(O(i));
end
subplot 311
plot(time,Y(:,3));
xlabel('time')
ylabel('phase')
subplot 312
plot(time,Y(:,2));
xlabel('time')
ylabel('Amplitude')
subplot 312
plot(O,U);
xlabel('phase')
ylabel('potential')
0 commentaires
Réponse acceptée
Walter Roberson
le 31 Juil 2022
tspan=[ti tf];
[time,Y] = ode45(f,tspan./tp,[sqrt(p);sqrt(p);0;0;0.01]);
When you use a vector of length 2 for the tspan, ode45 generates time outputs whenever it feels like it, so time will be a vector whose length is not easily predictable ahead of time. The Y output will have as many rows as there are entries in time and will have as many columns as there are entries in the initial conditions.
U(i) = -a.*(Y(:,5) - Y(:,4)).*O(i) + n.*((Y(:,1)./Y(:,2)) - (Y(:,2)./Y(:,1)))*sin(O(i));
The left hand side, U(i), is a scalar location.
On the right hand side, you use Y(:,1), Y(:,2), Y(:,4), and Y(:,5), each of which is a column vector with as many entries as the number of time values that were returned by ode45(). So the right hand side is a column vector of results, and you are trying to fit the column vector into a scalar location.
3 commentaires
Walter Roberson
le 31 Juil 2022
ti = 0; %inital time
tf = 70E-9;% final time
tspan = linspace(ti, tf, 1000);
%tspan = [ti, tf];
tp =1E-12;
T = 2E3;
p = 0.05;
k = (0.62).*10^(4);
c = 3E8;
N = 3.0;
n = k*(c/N)*tp
a = 5;
f = @(t,y) [
(y(4).*y(1) - n.*y(2).*sin(y(3)));
(y(5).*y(2) + n.*y(1).*sin(y(3)));
(-a.*(y(5)-y(4)) + n.*((y(1)./y(2)) - (y(2)./y(1))).*cos(y(3)));
(p - y(4)-(1+2.*y(4)).*(y(1)^2))./T;
(p - y(5)-(1+2.*y(5)).*(y(2)^2))./T;
];
[time,Y] = ode45(f,tspan./tp,[sqrt(p);sqrt(p);0;0;0.01]);
O = linspace(-0.08,0.08,10);
U = -a.*(Y(:,5) - Y(:,4)).*O + n.*((Y(:,1)./Y(:,2)) - (Y(:,2)./Y(:,1))).*sin(O);
subplot 311
plot(time,Y(:,3));
xlabel('time')
ylabel('phase')
subplot 312
plot(time,Y(:,2));
xlabel('time')
ylabel('Amplitude')
subplot 313
plot(O,U);
xlabel('phase')
ylabel('potential')
Plus de réponses (0)
Voir également
Catégories
En savoir plus sur Loops and Conditional Statements 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!