% Skript wurf_2.m, in dem der Wurf einer kleinen Masse
% von einer großen schwingeneden Masse untersucht wird
clear;

% ------- Parameter des Systems
m1 = 5;     c1 = 5;     k1 = 5;
m2 = 0.1;   % Masse die ''geworfen'' wird
g = 9.89;
% ------- Anregungskraft
ampl = 150;    % Amplitude der Anregungskraft
fanr = 4;      % Frequenz der Anregungskraft
% ------- Simulation mit Euler-Verfahren
Tfinal = 100;
dt = 0.001;
t = 0:dt:Tfinal;    nt = length(t);
% Initialisierungen
x = zeros(1,nt);    vx = x;        ax = x;
y = x;              vy = x;
Fe = ampl*cos(2*pi*fanr*t);

for k = 1:nt-1
    ax(k+1) = (Fe(k)-c1*vx(k)-k1*x(k))/m1-g;
    vx(k+1) = vx(k) + dt*ax(k+1);
    x(k+1) = x(k) + dt*vx(k+1);
    if ax(k) > - g & y(k) <= x(k+1)
        vy(k+1) = vx(k+1);
        y(k+1) = x(k+1);
    elseif y(k) < x(k);     
        vy(k+1) = vx(k+1);
        y(k+1)  = x(k+1); 
    elseif y(k) >= x(k);
        vy(k+1) = vy(k) + dt*(-g);
        y(k+1)  = y(k) + dt*vy(k+1);
    else
        vy(k+1) = vy(k) + dt*(-g);
        y(k+1)  = y(k) + dt*vy(k+1);
    end;
end;

figure(1);    clf;
nperioden = 5;
Tnperioden = 5/fanr;     nT = fix(Tnperioden/dt);
nd = nt - nT:nt;
subplot(211), plot(t(nd), ax(nd),'r');
title(['Beschleunigung; Frequenz der Anregung    fanr = ',...
    num2str(fanr),'  Hz']);    xlabel('Zeit in s');
axis tight; grid on;   hold on;
La = axis;    plot([La(1), La(2)], [-g, -g],'r');
hold off;
subplot(212), plot(t(nd), x(nd), t(nd), y(nd));
title('Lage der Masse 1 und der Masse 2')
xlabel('Zeit in s');        grid on;      
legend('x', 'y');           axis tight;   

