what's the problem in defining the function in for loop?

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

n = 0.6200

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

Unable to perform assignment because the left and right sides have a different number of elements.

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')

Answer 1

[Walter Roberson]

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.

[SAHIL SAHOO]

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);

now i change this part it takes lot of time to run.

then how should I modify my U ?

[Torsten]

O = linspace(-0.08,0.08,10);
for  k = 1:numel(O)
    U(:,k) = -a*(Y(:,5) - Y(:,4))*O(k) + n.*((Y(:,1)./Y(:,2)) - (Y(:,2)./Y(:,1)))*sin(O(k));
end

[Walter Roberson]

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

n = 0.6200

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')