function BolzenkraefteEv

global JSS          m1S         m2S       gS 
       JSS = 1/12 ; m1S = 0.5 ; m2S = 2 ; gS = 29.43 ;

maxstep = 0.01 ;                                                 % ... fuer MaxStep-Option
phi0 = [asin(1/3) ; 0] ; tspan = [0 2] ;                         % Anfangswerte, Zeitintervall
options = odeset ('MaxStep' , maxstep , 'Events' , @Ereignis) ;  % Event-Behandlung wird angekuendigt 

[t phiphip te phiphipe ereig] = ode45 (@Dgl , tspan , phi0 , options) ;
phi = phiphip(:,1) ; phip = phiphip(:,2) ;                       % Sortieren der Ergebnisse

subplot(1,2,1) ; plot (t , phi)  , grid on , title  ('\phi(t)')
subplot(1,2,2) ; plot (t , phip) , grid on , title  ('\omega(t)')
 
disp([' ']) ;
disp(['Rechnung mit MaxStep = ' , num2str(maxstep) , ':']) ;
for i=1:length(te)                                                         % Auswertung der "Events"
     disp(['omega(t = ' , num2str(te(i)) , ') = ' , num2str(phiphipe(i,2))]) ;
end
disp(['phi  (t_end = ' , num2str(t(end)) , ') = ' , num2str(phi(end))])  ; % Zusaetzliche Kontrolle:
disp(['omega(t_end = ' , num2str(t(end)) , ') = ' , num2str(phip(end))]) ; % Werte am Ende des Intervalls

% ===================================================================
function [value,isterminal,direction] = Ereignis (t,phiphip) 

phi  = phiphip(1) ;
value(1)      = phi ; % Ereignis 1: "Nulldurchgang"
isterminal(1) = 0 ;   % Kein Abbruch, aber registrieren ...
direction(1)  = 0 ;   % ... in jedem Fall (Hin- und Rueckweg) 

% ============================================================
% Function, die die "rechte Seite" des Dgl.-Systems definiert: 
function f = Dgl (t , phiphip)

global JSS          m1S         m2S       gS 

phi   = phiphip(1) ;
omega = phiphip(2) ;

phiP  = omega ;
zaehl = -((m2S-m1S)*omega^2*sin(2*phi)+(2*m1S+1)*gS*cos(phi)) ;
nenn  =  (0.5+2*m2S*(sin(phi))^2+2*m1S*(cos(phi))^2+2*JSS) ;

f = [phiP ; zaehl/nenn] ;