Keeping values of function when solving ODE
2 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
GCats
le 6 Mai 2021
Commenté : Bjorn Gustavsson
le 6 Mai 2021
Hello everyone!
I have a function describing a SDOF state space with a sinusoidal input y in one script as follows:
function dz = system(t,z,option,K,M,fe,y0)
dz=zeros(2,1);
y = y0 * sin(2*pi*fe*t); %base excitation
dz(1)=z(2);
dz(2)=(K*(y-z(1))/M
end
while in a separate script, I solve the ODE for different values of fe in a loop as follows:
clc
clear all
close all
% Linear Parameters________________________________________________________
K=5000;% linear stiffness
M=1; % system's mass
%Input Parameters__________________________________________________________
y0 = 1; % amplification of base input
T = 10;
fe = [1 2 3 4];
n_freq = length(fe);
for kk=1:n_freq
[tn,z] = ode15s('systeme1DDL_Dahl_bal',[t(1) t(end)],K,M,y0,fe(kk));
end
How can I save the value of y (the base excitation input) for each iteration?
Thank you in advance!
0 commentaires
Réponse acceptée
Bjorn Gustavsson
le 6 Mai 2021
As your ODE is written above we have to ask: Why bother? y is deterministically given provided your input parameters and the time. It should be trivial to calculate it for the times in you tn output. In case you have a slightly more complicated driving function you could perhaps do something like this:
function dz = system(t,z,option,K,M,fe,y_fcn)
dz=zeros(2,1);
y = y_fcn(t); %base excitation
dz(1)=z(2);
dz(2)=(K*(y-z(1))/M;
end
This would let you send in "any" function to your ODE:
% Linear Parameters________________________________________________________
K=5000;% linear stiffness
M=1; % system's mass
%Input Parameters__________________________________________________________
y0 = 1; % amplification of base input
T = 10;
fe = [1 2 3 4];
n_freq = length(fe);
z0 = 0; % you'll have to set this
for kk=1:n_freq
f_current = @(t) y0*sin(2*pi*fe(kk)*t);
[tn,z] = ode15s(@(t,z)system(t,z,K,M,fe(kk),f_current),[t(1) t(end)],z0);
y_exc{kk} = f_current(tn);
end
HTH
2 commentaires
Plus de réponses (0)
Voir également
Catégories
En savoir plus sur Ordinary Differential Equations 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!