Aufgabe 33-9a, Matlab

Das nachfolgend gelistete Matlab-Script löst die Aufgabenstellung a der Aufgabe 33-9. Es ist so geschrieben und wird nachfolgend entsprechend kommentiert, dass es als Muster-Script für entsprechende andere Aufgaben dienen kann.

Alle gelb eingetragen Zeilen können und sollten für die Lösung anderer Aufgaben unverändert übernommen werden. Zwingend ist nur die Anpassung der vier hellblau geschriebenen Bereiche, die von den Kommentarzeilen

% 11111111111 ANPASSEN AN DAS AKTUELLE PROBLEM VON HIER ... 111111111111111111111111111
...
% 11111111111 ... BIS HIER 111111111111111111111111111111111111111111111111111111111111

umschlossen sind. Eine Änderung der grün geschriebene Zeilen kann gegebenenfalls vorgenommen werden, ist aber nicht zwingend.

Im ersten anzupassenden Bereich werden zunächst die Zahlenwerte der Problemparameter definiert (Trägerlänge, Elastizitätsmodul, ...). Dabei ist darauf zu achten, dass für alle Größen einheitliche Dimensionen gewählt werden (in diesem Script wurden alle Größen auf mm und N eingestellt, damit dürfen auch für die "gemischten Dimensionen" nur N/mm², mm⁴, N/mm, Nmm verwendet werden).

Danach werden die Ritzschen Ansatzfunktionen (Vergleichsfunktionen) definiert. Das Script wurde für fünf Funktionen eingerichtet, was für die meisten Aufgabenstellungen ausreichend sein dürfte. Wenn - wie für die Aufgabenstellung 33-9 - nur drei Funktionen verwendet werden sollen, kann dies durch die Variable m gesteuert werden. Hier wurde m = 3 gesetzt, so dass nur die ersten drei Ansatzfunktionen in der nachfolgenden Rechnung verwendet werden.

In dem nach der Definition der Ansatzfunktionen folgenden Abschnitt werden die später benötigten ersten und zweiten Ableitungen der Ansatzfunktionen mit der Matlab-Function polyder berechnet. Diese Zeilen sollten nicht geändert werden.
Die Auswertung der Integrale für das Erzeugen der Elemente des aufzubauenden linearen Gleichungssystems erfolgt numerisch nach der Simpsonschen Regel (siehe "Einfache Formeln für die numerische Integration"). Die Anzahl der dafür erforderlichen Abschnitte wird mit der Variablen n festgelegt. Im vorliegenden Script gilt n = 100 (damit ergeben sich nS = n+1 Stützstellen für die Integration). Dieser Wert kann gegebenenfalls verändert werden (Bedingung ist ein geradzahliger Wert), dürfte aber in der Regel ausreichend sein.

Die Integrationsfomeln benötigen zur Berechnung der Werte für die Integranden die Funktionswerte der Biegesteifigkeit EI und der Linienlast q für alle nS Stützstellen. Im zweiten anzupassenden Bereich werden dafür die beiden Vektoren EIS bzw. qS aufgebaut.

Im Anschluss an diesen Bereich werden die Funktionswerte der Ansatzfunktionen und deren Ableitungen für alle Stützstellen mit der Matlab-Function polyval berechnet. Diese Zeilen sollten nicht geändert werden. Schließlich werden alle Ansatzfunktionen in ein Graphikfenster ausgegeben, so dass eine Sichtkontrolle der Zulässigkeit als Vergleichsfunktionen möglich ist.
Mit den beiden Vektoren EIS und qS werden die Integralanteile automatisch aufgebaut (auch dieser Bereich sollte nicht geändert werden), ergänzt werden müssen nur die diskreten Lasten und die diskreten Federn. Die diskreten Lasten werden im dritten anzupassenden Bereich eingebaut. Hier ist es die Einzelkraft F am rechten Rand (Punkt nS): + F*vS(nS,ii). Bei weiteren Einzelkräften kommen entsprechende zusätzliche Summanden hinzu, wobei nur die weiß geschriebenen Zeichen modifiziert werden müssen. Zu beachten ist, dass Einzelmomente mit der ersten Ableitung an dem entsprechenden Punkt (gespeichert in der Matrix vdS) multipliziert werden müssen. Ein Einzelmoment M in der Mitte des Trägers (Punkt iM = round (n/2+1)) würde mit dem Summanden + M*vdS(iM,ii) hinzugefügt werden müssen.
Im vierten anzupassenden Bereich werden schließlich noch diskrete Federn hinzugefügt, hier die Drehfeder am linken Rand (Punkt 1) mit dem Summanden + cT*vdS(1,ii)*vdS(1,jj). Während die Drehfederkonstante mit den ersten Ableitungen der Ansatzfunktionen multipliziert werden muss, würde eine der Verschiebung entgegen wirkende Feder mit der Federkonstanten c mit den Verschiebungswerten multipliziert werden, eine solche Feder in der Mitte des Trägers würde durch den zusätzlichen Summanden + c*vS(iM,ii)*vS(iM,jj) erfasst werden.

Der Rest des Scripts beinhaltet die Lösung des linearen Gleichungssystems (Berechnung der Ritzschen Ansatzparameter), die graphische Ausgabe der Biegelinie und des Momentenverlaufs und die Ausgabe spezieller Werte in das Command Window. Dieser Bereich muss gegebenenfalls an die entsprechende Aufgabenstellung angepasst werden.

% Ritzsches Verfahren für Biegetraeger mit Polynom-Ansatzfunktionen,
% Musterprogramm: Es koennen bis zu m=5 Ansatzfunktionen definiert werden,
% integriert wird nach der Simpson-Regel mit n Abschnitten (n kann beliebig
% - aber geradzahlig - gewaehlt werden).
%
% Wenn man - wie unten empfohlen - ein Feld mit den Werten fuer die Biegesteifigkeiten
% aller Stuetzstellen und ein weiteres Feld mit den Linienlastintensitaeten an den
% Stuetzstellen aufbaut, dann muss man nur die vier mit den Kommentarzeilen
% "ANPASSEN AN DAS AKTUELLE PROBLEM VON HIER ... ... BIS HIER"
% eingerahmten Bereiche an das aktuelle Problem anpassen.
%
% In dieser Datei sind in diesen Bereichen die Werte fuer die Aufgabe 33-9 eingetragen.

clear all

% 11111111111 ANPASSEN AN DAS AKTUELLE PROBLEM VON HIER ... 111111111111111111111111111
% Parameter:
tl = 1000 ;
E = 210000 ;
I0 = 200000 ;
q1 = 1 ;
F = 500 ;
cT = 100000000 ;

% Ansatzfunktionen und deren 1. und 2. Ableitungen:
m = 3 ; % Anzahl der Ansatzfuntionen (wenn m < 5 gesetzt ...
P1 = [1/tl 0] ; % ... wird, werden nur die ersten ...
P2 = [1/tl^2 0 0] ; % ... m Ansatzfunktionen verwendet)
P3 = [1/tl^3 0 0 0] ;
P4 = [1/tl^4 0 0 0 0] ;
P5 = [1/tl^5 0 0 0 0 0] ;
% 11111111111 ... BIS HIER 111111111111111111111111111111111111111111111111111111111111

P1D = polyder (P1);
P2D = polyder (P2);
P3D = polyder (P3);
P4D = polyder (P4);
P5D = polyder (P5);

P1DD = polyder (P1D) ;
P2DD = polyder (P2D) ;
P3DD = polyder (P3D) ;
P4DD = polyder (P4D) ;
P5DD = polyder (P5D) ;

% ------------------------------------------------------------------------------------
% Vorbereitung der numerischen Integration und der Ergebnisauswertung:
n = 100 ; % Anzahl der Abschnitte für numerische Integration
dz = tl / n ; % Breite eines Integrationsintervalls
nS = n + 1 ; % Anzahl der Stützstellen
zS = 0 : dz : tl ; % Koordinaten der Stützpunkte

EIS = zeros (nS , 1) ;
qS = zeros (nS , 1) ;

% 22222222222222 ANPASSEN AN DAS AKTUELLE PROBLEM VON HIER ... 222222222222222222222222
EIS(1:nS) = E * I0 ; % Biegesteifigkeit und ...
qS = q1 * (1 - zS/tl) ; % ... Linienlast
% 22222222222222 ... BIS HIER 222222222222222222222222222222222222222222222222222222222

% Funktionswerte und 1. und 2. Ableitung der Ansatzfunktionen an den Stützdtellen
vS (:,1) = polyval (P1 , zS)' ;
vS (:,2) = polyval (P2 , zS)' ;
vS (:,3) = polyval (P3 , zS)' ;
vS (:,4) = polyval (P4 , zS)' ;
vS (:,5) = polyval (P5 , zS)' ;
vdS (:,1) = polyval (P1D , zS)' ;
vdS (:,2) = polyval (P2D , zS)' ;
vdS (:,3) = polyval (P3D , zS)' ;
vdS (:,4) = polyval (P4D , zS)' ;
vdS (:,5) = polyval (P5D , zS)' ;
vddS(:,1) = polyval (P1DD , zS)' ;
vddS(:,2) = polyval (P2DD , zS)' ;
vddS(:,3) = polyval (P3DD , zS)' ;
vddS(:,4) = polyval (P4DD , zS)' ;
vddS(:,5) = polyval (P5DD , zS)' ;

subplot(311) , plot(zS , vS(:,1:m)) , axis ij , title('Ansatzfunktionen') ;

% ------------------------------------------------------------------------------------
% Aufbau des Gleichungssystems zur Bestimmung der Ansatzparameter:
K = zeros (m , m) ;
b = zeros (m , 1) ;

for ii = 1:m % Schleife über alle Gleichungen
   Summe = (qS(1)*vS(1,ii)+qS(nS)*vS(nS,ii)) ; % Rechte Seite ...
   faktor = 4 ;
   for k = 2:n % Numerische Integration ...
   Summe = Summe + qS(k)*vS(k,ii)*faktor ; % ... für Linienlastanteil ...
   if (faktor == 4) faktor = 2 ;
   else faktor = 4 ;
   end ;
   end
% 33333333333333 ANPASSEN AN DAS AKTUELLE PROBLEM VON HIER ... 3333333333333333333333333
   b(ii) = Summe*dz/3 + F*vS(nS,ii) ; % ... + Einzelkraftanteil,
   % + weitere Einzelkraftanteile + Anteile von aeusseren Momenten
% 33333333333333 ... BIS HIER ... 333333333333333333333333333333333333333333333333333333
for jj = 1:m % ii-te Zeile (Koeffizienten A)
   Summe = EIS(1)*vddS(1,ii)*vddS(1,jj) + EIS(nS)*vddS(nS,ii)*vddS(nS,jj) ;
   faktor = 4 ; % Numerische Integration ...
   for k = 2:n % ... nach Simpsonscher Regel ...
   Summe = Summe + EIS(k) * vddS(k,ii) * vddS (k,jj) * faktor ;
   if (faktor == 4) faktor = 2 ; % ... für Anteil aus der ...
   else faktor = 4 ; % ... Biegesteifigkeit ...
   end ;
   end
% 44444444444444 ANPASSEN AN DAS AKTUELLE PROBLEM VON HIER ... 4444444444444444444444444
   K(ii,jj) = Summe*dz/3 + cT*vdS(1,ii)*vdS (1,jj) ; % ... + Drehfederanteil
   % + weitere Anteile von Federn
% 44444444444444 ... BIS HIER ... 444444444444444444444444444444444444444444444444444444
end
end

% Loesen des Gleichungssystems
ai = zeros (5,1) ; % Maximal moegliche Zahl der Ansatzfznktionen: 5
aim = K\b ;
ai(1:m) = aim(1:m) ;

% Biegelinie graphisch darstellen:
vSchlange = vS * ai ;
subplot(312) , plot(zS , vSchlange) , axis ij , title('Durchbiegung') , ylabel('mm') ;

% Biegemoment graphisch darstellen:
v2Strich = vddS * ai ;
Mb = - EIS .* v2Strich;
subplot(313) , plot(zS , Mb) , title('Biegemoment') , ylabel('Nmm') ;

% Spezielle Werte:
DurchbiegungLinks = vSchlange(1)
DurchbiegungRechts = vSchlange(nS)
vStrich = vdS * ai ;
BiegewinkelLinks = vStrich (1)
BiegewinkelRechts = vStrich (nS)
BiegemomentLinks = Mb (1)
BiegemomentRechts = Mb (nS)

Die oben gelistete Datei ist als Aufgabe33_9a.m zum Download verfügbar .

Nachfolgend sind die Ergebnisse zu sehen, deren Qualität nur schwer zu beurteilen ist. Allerdings sieht man in der Graphik, dass das Biegemoment am rechten Rand doch deutlich von Null abweicht (der Wert im Command Window bestätigt das), so dass man sich mit dieser Rechnung mit 3 Ansatzfunktionen noch nicht zufrieden geben sollte. Das Script ist darauf vorbereitet, bis zu fünf Ansatzfunktionen zu verarbeiten, so dass eine Verbesserung schon durch das Hochsetzen des entsprechenden Parameters auf m = 5 erreicht wird. Dies wird bei der Behandlung der Aufgabenstellung b realisiert.

Aufgabe33_9aGraphik Aufgabe33_9aCW


	www.D@nkert.de


	D


	nkert.de