code for eulers method help

1 visualización (últimos 30 días)
Michelle
Michelle el 5 de Oct. de 2020
Comentada: Walter Roberson el 15 de Nov. de 2020
The derivative equation is y'=r(1-(y/L))*y-((p*y^2)/(q+y^2)) where L can be 5.4,8.1,16.3
Heres the code I have so far but it is not giving me a value.
p=1.2;
q=1;
r=0.65;
L1=5.4;
L2=8.1;
L3=16.3;
h=0.01; %step size
t=0:h:60; %0 days to 60 days
y=zeros(size(t));
y(1)=100; %inital amount of 100
n=numel(y);
for i=1:n-1
f=@(y)r(1-(y/L))*y-((p*y^2)/(q+y^2));
y(i+1)=y(i)+h*f(y)
end
  2 comentarios
Rik
Rik el 5 de Oct. de 2020
This code returns an error, because L is missing. After that it errors because r is indexed with 1-(y/L), which is unlikely to only return 1 (which is the only value it is allowed to have, because r is a scalar. Once you get past that, there are many products and divisions with y, and none of them are element-wise operations, which is probably not how it should be. After that I stopped looking for issues with your code.
Walter Roberson
Walter Roberson el 15 de Nov. de 2020
Michelle, if you feel that the question is unclear, then as you are the person who posted the question, you should be the one who clarifies it to make it more clear.

Iniciar sesión para comentar.

Respuestas (1)

James Tursa
James Tursa el 5 de Oct. de 2020
This needs to be outside of and prior to your loop, and r needs a multiply operator:
f=@(y)r*(1-(y/L))*y-((p*y^2)/(q+y^2));
This needs to have y indexed when passed into f:
y(i+1)=y(i)+h*f(y(i));
Finally, you should define the L you are going to use prior to creating the f function handle. E.g.,
L = L1;
f=@(y)r*(1-(y/L))*y-((p*y^2)/(q+y^2));
You can change L manually and run your code three times, or perhaps you can put all of the L selection stuff in a loop around your other code. E.g.,
for L=[L1,L2,L3]
f = @(y)r*(1-(y/L))*y-((p*y^2)/(q+y^2));
% your other code
end

Categorías

Más información sobre Programming en Help Center y File Exchange.

Community Treasure Hunt

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

Start Hunting!

Translated by