How to store a variable in a vector?

Hola, he tenido una serie de problemas con le siguiente código:
%Función para el cálculo del periodo de un sistema L-V
function calculo_periodo(~)
%***********************************************************%
%Termina cuando y retorna al punto incial (donde su ángulo
%con mu es el mismo que el ángulo entre eta y mu)
function [g,isterm,dir] = pits(t,y)
sig = sign(eta(1)-mu(1)+realmin);
theta1 = atan2(sig*(y(2)-mu(2)),sig*(y(1)-mu(1)));
theta0 = atan2(sig*(eta(2)-mu(2)),sig*(eta(1)-mu(1)));
g = theta1 - theta0;
isterm = t > 1;
dir = 1;
end
%************************************************************%
%Ciclo anidado para el cálculo del periodo, donde j=C e i=eta1
%for j=205:220
for i=0.001:0.001:0.999999
mu = [1 1]; % Punto de Equilibrio.
eta = [i 1]; % Condiciones Iniciales.
F = @(t,y) [(1-(y(2))/mu(2))*(y(1));
-220*(1-(y(1))/mu(1))*(y(2))];
opts = odeset('reltol',1.e-8,'events',@pits);
[t,~,~] = ode45(F,[0 Inf],eta,opts);
A = [i,t(end)]
%************************************************************%
%Generación de gráficos
hold on
plot(i,t(end),'.b')
hold off;
%end
end
A
end
the idea of the code is to find the period function of a lotka volterra system. For this, it is necessary to store the value of t (end) inside a vecto. How can I save t (end) in a vector?
If someone can help me, I appreciate it.
Greetings.
Pedro.

 采纳的回答

iVec = 0.001:0.001:0.999999;
result = zeros(1, numel(iVec));
for k = 1:numel(iVec)
i = iVec(k);
... your code
result(k) = t(end);
end

5 个评论

Hello Jan, thanks for your response. When you put "my code" is all my code, or is it a part of my code?
Thanks.
%Función para el cálculo del periodo de un sistema L-V
function calculo_periodo(~)
%************************************************************%
%Ciclo anidado para el cálculo del periodo, donde j=C e i=eta1
iVec = 0.001:0.001:0.999999;
result = zeros(1, numel(iVec));
for k = 1:numel(iVec)
i = iVec(k);
mu = [1 1]; % Punto de Equilibrio.
eta = [i 1]; % Condiciones Iniciales.
F = @(t,y) [(1-(y(2))/mu(2))*(y(1));
-220*(1-(y(1))/mu(1))*(y(2))];
opts = odeset('reltol',1.e-8,'events',@pits);
[t,~,~] = ode45(F,[0 Inf],eta,opts);
A = [i,t(end)]
%************************************************************%
%Generación de gráficos
hold on
plot(i,t(end),'.b')
hold off;
%end
result(k) = t(end);
end
A
end
%***********************************************************%
%Termina cuando y retorna al punto incial (donde su ángulo
%con mu es el mismo que el ángulo entre eta y mu)
function [g,isterm,dir] = pits(t,y)
sig = sign(eta(1)-mu(1)+realmin);
theta1 = atan2(sig*(y(2)-mu(2)),sig*(y(1)-mu(1)));
theta0 = atan2(sig*(eta(2)-mu(2)),sig*(eta(1)-mu(1)));
g = theta1 - theta0;
isterm = t > 1;
dir = 1;
end
Excuse me Jan, I have tried to parameterize (with exponentials) the differential equation with m and n=0.5, I show you:
function periodo
function [g,isterm,dir] = pits(t,y)
sig = sign(eta(1)-mu(1)+realmin);
theta1 = atan2(sig*(y(2)-mu(2)),sig*(y(1)-mu(1)));
theta0 = atan2(sig*(eta(2)-mu(2)),sig*(eta(1)-mu(1)));
g = theta1 - theta0;
isterm = t > 1;
dir = 1;
end
iVec = 0.001:0.001:0.999999;
result = zeros(1, numel(iVec));
for m =0.1:0.1:0.9
for k = 1:numel(iVec)
i = iVec(k);
mu = [1 1]; % Punto de Equilibrio.
eta = [i 1]; % Condiciones Iniciales.
F = @(t,y) [(1-(y(2))^(m)/mu(2))*(y(1))^(0.5);
-10.3*(1-(y(1))^(0.5)/mu(1))*(y(2))^(m)];
opts = odeset('reltol',1.e-8,'events',@pits);
[t,~,~] = ode45(F,[0 Inf],eta,opts);
result(k) = t(end);
%************************************************************%
end
hold on
plot(0.001:0.001:0.999999,result)
hold off
end
%A(:,1) = 0.001:0.001:0.999999;
%A(:,2) = result;
%dlmwrite('data_2.txt',A,',');
end
but, I receive the following error message:
Error using atan2
Inputs must be real.
Error in periodo/pits (line 7)
theta1 = atan2(sig*(y(2)-mu(2)),sig*(y(1)-mu(1)));
Error in odezero (line 28)
[vnew,isterminal,direction] = feval(eventfun,tnew,ynew,eventargs{:});
Error in ode45 (line 406)
odezero(@ntrp45,eventFcn,eventArgs,valt,t,y,tnew,ynew,t0,h,f,idxNonNegative);
Error in periodo (line 26)
[t,~,~] = ode45(F,[0 Inf],eta,opts);
the function spits is for determinate the time of the revolution (the period). How can I improve my code to do not have de error mensaje?
Again, thank you for yout help!
@Pedro: Use the debugger to find out, why the values are not real. Insert some code to check, if the trajectories are real in the function to be integrated. Then set a breakpoint to stop Matlab.
The debugging is impeded massively, if you use anonymous functions. Prefer functions instead:
function dy = F(t,y)
dy = [(1-(y(2))^(m)/mu(2))*(y(1))^(0.5);
-10.3*(1-(y(1))^(0.5)/mu(1))*(y(2))^(m)];
if ~isreal(dy)
disp('not real') % Set a break point here
end
end
Hello Jan, again, thank you for yout help!

请先登录,再进行评论。

更多回答(1 个)

类别

Community Treasure Hunt

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

Start Hunting!

Translated by