ex12.m

function ex12()
% Loesung von Aufgabe 12

l = 1.8;  % laenge in m
a = 1.5;  % Federangriffspunkt in m
m = 10;   % in kg
c = 100;  % in N/m
c2 = 50;
g = 9.81; % in 

% Kreisfrequenz
omega = sqrt(a^2*c/(l^2*m) - g/l)
T = 2*pi/omega

% anderes c
omega2 = a^2*c2/(l^2*m) - g/l

% Teil c
% starke Feder

% Gleichgewichtslage
phig = acos(m*g*l/(a^2*c));
phi0 = [0.1, 0.4, phig-0.01, phig+0.01]; % Anfangswerte

func = @(t,y) dgl12(t, y, g/l, a^2*c/(m*l^2));  % Parameter in dgl12 einsetzen
[t1, y1] = ode45(func, [0 2*T], [phi0(1) 0]);
[t2, y2] = ode45(func, [0 2*T], [phi0(2) 0]);
[t3, y3] = ode45(func, [0 2*T], [phi0(3) 0]);
[t4, y4] = ode45(func, [0 2*T], [phi0(4) 0]);

plot(t1, y1(:,1), t2, y2(:,1), t3, y3(:,1), t4, y4(:,1));
hold on;
plot([0 2*T], [0 0], 'k-');         % t-Achse
plot([0 2*T], [phig phig], 'k:');   % Gleichgewichtslage 
plot([T T], [-phig 2*phig], 'k:');  % Schwingungsdauer, linearisiert 
hold off;

% Beschriftungen
title('Schwingungen bei verschiedenen Auslenkungen', 'FontSize', 16);
ylabel('phi', 'FontSize', 14);
xlabel('t [s]', 'FontSize', 14);
axis([0, 2*T, -phig, 2*phig]);
legend('0.1', '0.4', 'phig-0.01','phig+0.01', 'Location', 'Best')

F = getframe(gcf);
imwrite(F.cdata, 'bild71.png');

% Teil c
% schwache Feder

func = @(t,y) dgl12(t, y, g/l, a^2*c2/(m*l^2));  % Parameter in dgl12 einsetzen
[t5, y5] = ode45(func, [0 20], [0.1 0]);
plot(t5, y5(:,1), 'b-', [0 20], [0 0], 'k-');

% Beschriftungen
title('Bewegung bei schwacher Feder', 'FontSize', 16);
ylabel('phi', 'FontSize', 14);
xlabel('t [s]', 'FontSize', 14);

F = getframe(gcf);
imwrite(F.cdata, 'bild72.png');

%----------------------------------------------------------------------
function dydt = dgl12(t, y, a, b)
% rechte Seite der DGL bei Aufgabe 12
% a = g/l, b = a^2*c/(m*l^2)

dy1 = y(2);
dy2 = a*sin(y(1)) - b*sin(y(1)).*cos(y(1));
dydt = [dy1 dy2]';