%   Script Frequg.m
%
%   Verhalten eines Feder-Masser-Schwingers bei
%   unterschiedlichen Frequenzen
%
% J. Lunze
% für 3. Auflage: 2.4.2001
% für 4. Auflage: 18.9.2003
% für 5. Auflage: 18.4.2005
% für 6. Auflage: 4.12.2006
% für 7. Auflage: 28.4.2008
% für 9. Auflage: 13.2.2012
% für 10. Auflage: 27.1.2014
% für 12. Auflage: 03.07.2019
%
%   
clear
close all
%
%          Feder-Masse-Schwinger aus Beispiel 4.2
%
%  Maßeinheiten: Zeit in Sekunden
%                Kraft in Newton
%                Entfernung in Metern
%
% u(t) = Erregung des Feder-Masse-Schwingers an der Feder
% y(t) = Position der Masse
%
m=0.2;              % Masse in kg
cF=1;                   % Federkonstante in N/m
dF=0.8;             % Dämpfungskonstante in Ns/m
a2=1;                   % Koeffizienten der Differenzialgleichung
a1=cF/dF;               % a2 yddot + a1 ydot + a0 y = b0 u
a0=cF/m;
b0=cF/m;
A=[0     1;         % Zustandsraummodell
   -a0 -a1];
b=[0; 1];
c=[b0 0];
d=0;
FMS=ss(A, b, c, d);
%% Erregung durch sinusförmige Eingangsgröße
%  Anfangsbedingung x0 so gewählt, dass das Übergangsverhalten
%  verschwindet
% 
figure(1)
subplot(3,1,1);
Tmax=50;
omega=0.4;
T=0:0.02/omega:Tmax/omega;
U=sin(omega*T);
Y=lsim(FMS, U, T);
T4=4*2*pi/omega;
plot(T,U,'b--', T,Y, 'b', [T4 T4],[-10 10],'k:', [T4 Tmax/omega],[0 0],'k:')
axis([0 Tmax/omega -2 2]);
xticks([0 20 40 60 80 100 120]);
yticks([-2 0 2]);
xlabel('');
latexylabel('$$y, u$$  in cm');
rotateY;
latextitle('Feder-Masse-Schwinger bei sinusf\"ormiger Erregung');
%
subplot(3, 1, 2);
omega=2.1;
T=0:0.1/omega:Tmax/omega;
U=sin(omega*T);
Y=lsim(FMS, U, T);
T4=4*2*pi/omega;
plot(T,U,'b--', T,Y, 'b', [T4 T4],[-10 10],'k:', [T4 Tmax/omega],[0 0],'k:')
axis([0 Tmax/omega -2 2]);
xticks([0 5 10 15 20 25]);
yticks([-2 0 2]);
xlabel('');
latexylabel('$$y, u$$  in cm');
rotateY;
%
subplot(3, 1, 3);
omega=5;
T=0:0.1/omega:Tmax/omega;
U=sin(omega*T);
Y=lsim(FMS, U, T);
T4=4*2*pi/omega;
plot(T,U,'b--', T,Y, 'b', [T4 T4],[-10 10],'k:', [T4 Tmax/omega],[0 0],'k:')
axis([0 Tmax/omega -2 2]);
xticks([0 5 10]);
yticks([-2 0 2]);
latexxlabel('$$t$$ in s');
latexylabel('$$y, u$$  in cm');
rotateY;


%%   Bodediagramm des Feder-Masse-Schwingers
%
figure(2)
z=[1];
n=[a2/b0, a1/b0, a0/b0];
FMS=tf(z, n)
[Betrag, Phase, Omega]=bode(FMS);
subplot(2, 1, 1)
semilogx(Omega, 20*log10(Betrag(:)), 'b');
hold on
AusgewFrequenzen=[0.4, 2.1, 5];
[BetragPunkte, PhasePunkte]=bode(FMS, AusgewFrequenzen);
plot(AusgewFrequenzen(1), 20*log10(BetragPunkte(1)), 'bo');
plot(AusgewFrequenzen(2), 20*log10(BetragPunkte(2)), 'bo');
plot(AusgewFrequenzen(3), 20*log10(BetragPunkte(3)), 'bo');
latextitle('Feder-Masse-Schwinger');
latexylabel('$$|G|$$ in dB');
rotateY;
xleer;
grid('on');
axis([0.1, 10, -40 20]);
yticks([-40 -20 0 20 40]);
hold off
subplot(2, 1, 2);
semilogx(Omega, Phase(:), 'b');
hold on;
plot(AusgewFrequenzen(1), PhasePunkte(1), 'bo');
plot(AusgewFrequenzen(2), PhasePunkte(2), 'bo');
plot(AusgewFrequenzen(3), PhasePunkte(3), 'bo');
latexylabel('$$\phi$$ in Grad');
rotateY;
latexxlabel('$$\omega$$ in rad/s');
xticks([0.1 1 10]);
yticks([-360 -270 -180 -90 0]);
axis([0.1, 10, -180, 0]);
grid('on');
hold off

%% Herstellung von .ps und .eps-Bildern
%
Epsbildc15('frequg');
