Keeping values of function when solving ODE
2 views (last 30 days)
Show older comments
GCats
on 6 May 2021
Commented: Bjorn Gustavsson
on 6 May 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 Comments
Accepted Answer
Bjorn Gustavsson
on 6 May 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 Comments
More Answers (0)
See Also
Categories
Find more on Ordinary Differential Equations in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!