Doppelpendel (einfache Animation)

Aufgabe

Ein Doppelpendel wird definiert durch die beiden Pendelmassen m₁ und m₂, die auf die jeweiligen Schwerpunkte bezogenen Massenträgheitsmomente J_S₁ und J_S₂, die Schwerpunktabstände von den Drehpunkten s₁ und s₂ und den Abstand l₁ der beiden Drehpunkte voneinander.

Die Bewegung soll durch die Funktionen φ₁(t) und φ₂(t) beschrieben werden, die für das Zeitintervall t = 0 ... 10 s zu berechnen sind.

Die Aufgabe wird in den Kapiteln "Prinzipien der Mechanik" und "Verifizieren von Computerrechnungen" behandelt.

Die gegebenen Werte gelten für zwei schlanke Stäbe gleicher Masse und gleicher Länge. Sie sollen aus der nebenstehend skizzierten Anfangslage ohne Anfangsgeschwindigkeiten freigelassen werden, so dass folgende Anfangsbedingungen gelten:

Matlab-Script

function DoppelPendelAnimation (phi1Anf , phi2Anf)

if nargin == 2                          % Start mit vorgegebenen Anfangswerten?
   x0 = [phi1Anf ; phi2Anf ; 0 ; 0] ;   % Anfangswerte  [phi1 ; phi2 ; omega1 ; omega2]
else
   x0 = [pi/2 ; 0 ; 0 ; 0] ;            % Default-Anfangswerte
end

% Parameter (global definiert, damit Wertaenderung nur an einer Stelle erforderlich):
global mdm  J1  J2  s1  s2  gdl
mdm = 1 ; J1 = 1./12 ; J2 = 1./12 ; s1 = .5 ; s2  = .5 ; l1 = 1 ; l2 = 1 ; gdl = 9.81/l1 ;

tspan = [0  10] ;                       % Zeitintervall

% Integration des Anfangswertproblems:
clf
options = odeset('Mass' , @Massenmatrix , 'MaxStep' , 0.01 , ...
                 'OutputFcn' , @odeplot , 'OutputSel' , [2]) ;
% ... zeichnet die Funktion phi2(t) waehrend der ode45-Berechnung in das Graphik-Fenster
[t x]   = ode45 (@RechteSeite , tspan , x0 , options) ;

% Animation der Bewegung, zunaechst "Sortieren der Ergebnisse":
phi1i = x(:,1) ;    % x ist die Matrix der Ergebnisse,
phi2i = x(:,2) ; ;  % phi1i und phi2i sind Vektoren

clf
tt = tspan(1) ;
xx = [0  l1*sin(phi1i(1))  l1*sin(phi1i(1))+l2*sin(phi2i(1))] ;
yy = [0 -l1*cos(phi1i(1)) -l1*cos(phi1i(1))-l2*cos(phi2i(1))] ;
plot ([0 l1/8 -l1/8 0],[0 l1/6 l1/6 0] , 'r-') ;
hold on
p=plot (xx,yy,'o-','EraseMode','xor') ;
axis (1.1*[-(l1+l2) (l1+l2) -(l1+l2) (l1+l2)]) ;
zeit = title (sprintf('t = %8.2f s' , tt)) ;
set (gca , 'userdata' , zeit)

deltatt  = 1/25     ;    % 25 Bilder pro Sekunde
zeitlupe = 1        ;    % Verzögerungsfaktor
cstart   = cputime  ;
Zeitschritte = length (t)

for ii=1:Zeitschritte    % Schleife über alle berechneten Zeitschritte
   while (t(ii) >= tt)   % Ausgabe jeweils nur nach deltatt
       xx = [0  l1*sin(phi1i(ii))  l1*sin(phi1i(ii))+l2*sin(phi2i(ii))] ;
       yy = [0 -l1*cos(phi1i(ii)) -l1*cos(phi1i(ii))-l2*cos(phi2i(ii))] ;
       set(p , 'XData' , xx , 'YData' , yy)
       while (cputime-cstart < tt*zeitlupe)  % warteschleife
       end
       zeit = get (gca , 'userdata');
       set(zeit , 'string' , sprintf('t = %8.2f s' , tt))
       drawnow
       tt = tt + deltatt ;
   end
end

% ============================================================
% Funktion, die die "Massenmatrix" des Dgl.-Systems definiert:
function M = Massenmatrix (t , xv)

global mdm  J1  J2  s1  s2  gdl

phi1   = xv(1) ;
phi2   = xv(2) ;
omega1 = xv(3) ;
omega2 = xv(4) ;

c = cos(phi1-phi2);
M = [ 1  0      0             0      ;
      0  1      0             0      ;
      0  0  s1^2+J1+mdm    mdm*s2*c  ;
      0  0    mdm*s2*c   mdm*s2^2+J2 ] ;

% ============================================================
% Funktion, die die "rechte Seite" des Dgl.-Systems definiert:
function xvp = RechteSeite (t , xv)

global mdm  J1  J2  s1  s2  gdl

phi1   = xv(1) ;
phi2   = xv(2) ;
omega1 = xv(3) ;
omega2 = xv(4) ;

s = sin(phi1-phi2);

phi1p   = omega1 ;
phi2p   = omega2 ;
omega1p = -mdm*s2*omega2^2*s - (s1+mdm)*gdl*sin(phi1) ;
omega2p =  mdm*s2*omega1^2*s -   mdm*s2*gdl*sin(phi2) ;

xvp = [phi1p ; phi2p ; omega1p ; omega2p] ;

Zu den hellblauen Zeilen: Die Funktion "DoppelPendelAnimation" kann (aus dem Command Window) aufgerufen werden mit Angabe der Anfangsauslenkungen der beiden Stäbe oder ohne Parameter. Im letztgenannten Fall werden Default-Werte (das sind die Werte der Aufgabenstellung) verwendet.
Zu den weißen Zeilen: Mit der Option "Mass" wird angekündigt, dass eine Massenmatrix zu berücksichtigen ist, der folgende Parameter gibt den Namen der dafür bereitgestellten Funktion an. Die Option "MaxStep" ist zwingend, um über ein Zeitintervall von 10 Sekunden ausreichend genaue Ergebnisse zu erzielen (vgl. hierzu die Seite "Doppelpendel - Matlab (mit Kontrollen)".

Mit den beiden Keywords "OutputFcn" und "OutputSel" wird veranlasst, dass schon während der Berechnung die Funktion φ₂(t) in das Graphik-Fenster gezeichnet wird (damit etwas passiert, während die Ergebnisse berechnet und gespeichert werden).
Zu den grünen Zeilen: Dort wird die Animation mit den vorab berechneten Werten φ₁(t) und φ₂(t) realisiert. Mit dem Parameter "zeitlupe" kann die Bewegung verlangsamt dargestellt werden.

Download

Das oben gelistete Matlab-Script steht als DoppelpendelAnimation.m zum Download zur Verfügung.

Aufgabe

Bewegungs-Differenzialgleichungen, Vorbereitung der Matlab-Rechnung

Matlab-Script

Ergebnisse

Download