Aufgabe
Ein Doppelpendel wird definiert durch die beiden Pendelmassen
m1
und m2, die auf die jeweiligen
Schwerpunkte bezogenen
Massenträgheitsmomente JS1
und JS2, die
Schwerpunktabstände
von den Drehpunkten s1
und s2 und den Abstand
l1 der beiden
Drehpunkte voneinander.
Die Bewegung soll durch die Funktionen
φ1(t) und
φ2(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:
Bewegungs-Differenzialgleichungen, Vorbereitung der Matlab-Rechnung
Unter Verwendung der in der Aufgabenstellung skizzierten Koordinaten
werden die Bewegungs-Differenzialgleichungen im Kapitel
"Prinzipien der Mechanik" hergeleitet. Man findet sie
einschließlich der Vorbereitung für die Matlab-Rechnung
auf der Seite
"Doppelpendel - Berechnung mit Matlab".
Hier soll nur eine Modifikation des dort gelisteten Matlab-Scripts
gezeigt werden, die mit einfachsten Mitteln eine Animation
erzeugt.
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
φ2(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
φ1(t) und
φ2(t)
realisiert. Mit dem Parameter "zeitlupe" kann die Bewegung
verlangsamt dargestellt werden.
Ergebnisse
Während der numerischen Lösung des nichtlinearen
Anfangswertproblems
wird die Funktion
φ2(t)
in das Graphik-Fenster gezeichnet. Diese wird nach dem Ende der Rechnung
gelöscht,
und die beiden Funktionen
φ1(t) und
φ2(t)
werden benutzt, um die Bewegung als Animation darzustellen.
Rechts sieht man einen Schnappschuss der Bewegung.
Die oben im Titelfeld des Graphik-Fensters angegebene Zeit
ist immer (auch bei Zeitlupen-Darstellung) der zur
gerade sichtbaren Stellung der beiden Pendel gehörende Wert.