My code won't run, related to ode function

2 vues (au cours des 30 derniers jours)
Dardenella Finegan
Dardenella Finegan le 24 Mai 2020
Commenté : Stephen23 le 14 Mai 2021
%E_as_a_function_of_time
m = 1; % mass [kg]
k = 4; % spring constant [N/m]
c = 1; % friction coefficient [Ns/m]
omega0 = sqrt(k/m); p = c/(2*m);
y0 =0.1; v0 = 0; % initial conditions
[t,Y] = ode45(@f,[0,10],[y0,v0],omega0,[],p); % solve for 0<t<15
y = Y(:,1); v = Y(:,2); % retrieve y, v from Y
figure(1); plot(t,y,'ro-',t,v,'b+-');% time series for y and v
grid on; axis tight;
%---------------------------------------------------
function dYdt = f(t,Y,omega0,p); % function defining the DE
y = Y(1); v = Y(2);
dYdt=[ v ; -omega0^2*y-2*p*v]; % fill-in dv/dt
end

Réponses (2)

David Goodmanson
David Goodmanson le 24 Mai 2020
Modifié(e) : David Goodmanson le 24 Mai 2020
Hi Dardenella,
try the same thing only with
[t,Y] = ode45(@(t,Y) f(t,Y,omega0,p), [0,10],[y0,v0]); % solve for 0<t<10 [not 15]
Sometimes it's less than obvious where the @'s should go. But the function call has to match the argument list it the function definition.

Poojana mathlab
Poojana mathlab le 14 Mai 2021
%mass [kg] %spring constant
k=4; %spring constant[N/m]
omega0=sqrt(k/m)
x0=0.1; v0=0 %initial conditions
[t,X]=ode45(@f,[0,10],[y0,v0],[],omega0); %solve for 0<t<10
x= X(:,1);v=X(:,2); %retrieve y,v from X
figure(1); plot(t,x,'b+-',t,v,'ro-'); %time series for x and v
%--------------------------------------
function dYdt= f(t,Y,omega0);x=X(1);v=X(2);
dYdt=[v;-omega062*x];end
doesnt work
  1 commentaire
Stephen23
Stephen23 le 14 Mai 2021
Works perfectly, following the answer and comment above and fixing all of the mismatching variable names:
m=1;%mass [kg]
k=4;%spring constant[N/m]
omega0=sqrt(k/m);
fun = @(t,Y) f(t,Y,omega0);
[t,Y]=ode45(fun,[0,10],[0.1,0]); %solve for 0<t<10
x=Y(:,1);
v=Y(:,2);
plot(t,x,'b+-',t,v,'ro-')
%--------------------------------------
function dYdt = f(t,Y,omega0)
x=Y(1);
v=Y(2);
dYdt=[v;-omega0*x];
end

Connectez-vous pour commenter.

Catégories

En savoir plus sur Programming 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!

Translated by