10. Differential-Algebraische Gleichungen (2)

In dieser Übung wird das Schwingen eines Pendels simuliert. Dabei ist es möglich, die Bewegung als ODE- oder auch als DAE-System zu beschreiben. Wie unterscheiden sich das jeweils resultierende System und dessen Ergebnisse?

Situation

Die Bewegung eines Pendels kann sowohl als gewöhnliche Differentialgleichung (ODE) sowie auch als Differential-Algebraisches System (DAE) beschrieben werden.

1. Programmimplementierung

Lösen Sie die Pendelbewegung mithilfe von ein_pendel_dae.m.

Pendel
%% file ein_pendel_dae.m
% if your Octave kernel does not plot, write `graphics_toolkit("gnuplot")`
global g m
g = 9.81; l = 0.15; m = 0.05;
y0 = [l 0 0 0 1]; tspan = [0, 3];
M = eye(5); M(5,5) = 0;
options = odeset('Mass',M,'RelTol',1.0e-5,'Stats','on');
tic
[t,y] = rodasp2(@dgl, tspan, y0,options);
toc
axis equal; hold on;
lges=1.2*l; r1 = 0.005; d1=2*r1;
axis ([-lges lges -lges 0.1*lges])
linie1 = line([0 y0(1)],[0 y0(3)]);
kugel1 = rectangle ('Position',[y0(1)-r1 y0(3)-r1 d1 d1],...
'Curvature', [1 1], 'Facecolor', 'r');
for i = 1:length(t)
    set(kugel1, 'Position', [-r1+y(i,1) -r1+y(i,3) d1 d1]);
    set(linie1, 'Xdata', [0 y(i,1)], 'Ydata', [0 y(i,3)]);
    drawnow
    if(i<length(t)) pause ((t(i+1)-t(i)));
    end
end

function dz = dgl(t,z)
    global g m
    dz(1,1) = z(2);
    dz(2,1) = 2*z(5)*z(1)/m;
    dz(3,1) = z(4);
    dz(4,1) = -g + 2*z(5)*z(3)/m;
    % dz(5,1) = z(1)^2 + z(3)^2-l^2; %Index = 3
    % dz(5,1) = z(1)*z(2) + z(3)*z(4); %Index = 2
    dz(5,1) = 2*(z(2)^2+z(1)*dz(2,1)+z(4)^2+z(3)*dz(4,1)); %Index
end

Die Funktion rodasp2 finden Sie hier.

2. Pendellänge

Leiten Sie aus den Simulationsergebnissen die Pendellänge her und überprüfen Sie, ob sie sich mit der Zeit verändert.

Plotten Sie die eventuelle Längenänderung über die Simulationszeit.

% if your Octave kernel does not plot, write `graphics_toolkit("gnuplot")`
% your code here

3. Systemdämpfung

Bringen Sie eine Dämpfung in das System ein.

4. Pendel als ODE

Das Pendel kann auch als ODE beschrieben werden (Drehbewegung). Der Drehwinkel definiert sich dann wie folgt:

\[\ddot{\varphi}=-\frac{g}{L}\sin(\varphi)\]

Dabei ist der Winkel so festgelegt, dass \(\varphi(t)=0\) für ein senkrecht nach unten zeigendes Pendel gilt. Implementieren Sie diese Gleichung und vergleichen Sie die Lösungen der beiden Ansätze.

% your code here

Lösung

- Hinweis: Dieses Video ist von YouTube aus eingebunden und nicht Teil des frei lizenzierten Materials! -

Zusatzaufgabe: Doppelpendel

Leiten Sie die Gleichungen eines Doppelpendels her und implementieren Sie diese.