TMCu, Seite 261, MATLAB, Bandmatrix

Ein MATLAB-Beispiel
wird hier ausführlich
beschrieben

Beispiel 1 auf den Seiten 260 bis 262 (Einstiegsbeispiel zum Differenzenverfahren), Lösung mit MATLAB, Ausnutzung der Bandstruktur der Koeffizientenmatrix

Zunächst wird das Gleichungssystem mit 9 Unbekannten

in einem MATLAB-Script erzeugt, um die Speicherung der Koeffizientenmatrix als Bandmatrix mit den Bandweiten IBWL=4 und IBWR=3 zu demonstrieren:

% Gleichungssystem auf Seite 261, Matrix A
% als Bandmatrix mit IBWL=4 und IBWR=3
% (Hauptdiagonalelemente von A in Spalte 4)

clear all

A = [0 0 0 0 0 1 ; % Matrix A ...
   0 0 0 -1 0 1 ; % ist als ...
   0 1 -4 6 -4 1 ; % Bandmatrix ...
   0 1 -4 6 -4 1 ; % gespeichert, ...
   0 1 -4 6 -4 1 ; % b muss kein ...
   0 1 -4 6 -4 1 ; % Spaltenvektor ...
   0 1 -4 6 -4 1 ; % sein
   -1 2 0 -2 1 0 ;
   1 -2 1 0 0 0] ;

b = [0 0 1/256 1/256 1/256 1/256 1/256 0 0] ;

v = gabamp (A , b , 4)

An die Stelle der 9*9 -Matrix tritt die 9*6-Bandmatrix mit den Hauptdiagonalelementen in der 4. Spalte, was gabamp beim Aufruf mitgeteilt werden muss.

Nach dem Abarbeiten des links zu sehenden Scripts S261Band1.m findet man das Ergebnis im “Command Window” (Bild rechts).

Das folgende Script S261Band2.m zeigt die Definition des beliebig großen Gleichungssystems in den wenigen hellblauen Zeilen. Die Möglichkeit unterschiedlich feiner Unterteilung des Trägers wird nur über die Angabe des Wertes für n_A (in der weißen Zeile am Anfang) gesteuert.

% Differenzenverfahren, Beispiel auf den Seiten 260 bis 262,
% Ausnutzung der Bandstruktur der Koeffizientenmatrix

clear all

L = 1 ;
q0 = 1 ;
EI = 1 ;
nA = 1000 ; % Anzahl der Abschnitte (nur dieser Wert muss geaendert werden)
n = nA + 5 ; % Anzahl der Gleichungen
h = L / nA ;

hdsp = 4 ; % Hauptdiagonalelemente stehen in der 4. Spalte der Bandmatrix A

A = zeros (n,7) ; % Nullmatrix fuer Aufnahme der Bandmatrix
b = zeros (n,1) ; % Nullvektor

A(1:2,:) = [ 0 0 0 0 0 1 0 ;
   0 0 0 -1 0 1 0 ] ; % Randbedingungen links: Einspannung
A(n-1:n,:) = [-1 2 0 -2 1 0 0 ;
   1 -2 1 0 0 0 0 ] ; % Randbedingungen rechts: Freier Rand

for i = 3:n-2
   A(i,:) = [0 1 -4 6 -4 1 0] ; % Standardgleichungen
   b(i) = q0*h^4/EI ;
end

% Lösen des Gleichungssystems
tic ;
v = gabamp (A , b , hdsp) ; % Berechnung der Durchsenkung
toc ;

format long
vEnd = v(n-2)
phiEnd = (-v(nA+2) + v(nA+4)) / (2*h)
Mblinks = -EI/h^2 * (v(2) - 2*v(3) + v(4))

for i = 3:n-2
   Mb(i) = - EI/h^2 * (v(i-1) - 2*v(i) + v(i+1)) ;
end

z = 0 : h : L ;
subplot (2 , 1 , 1) ; plot (z , v (3:n-2)) , axis ij , title ('Verschiebung')
subplot (2 , 1 , 2) ; plot (z , Mb(3:n-2)) , title ('Biegemoment')

Nach der Lösung des Gleichungssystems mit dem Script gabamp werden mit den Zeilen, die nicht mit einem Semikolon abgeschlossen sind, vier Ergebnisse in das “Command Window” ausgegeben: Neben der für die Lösung des Gleichungssystems benötigten Zeit die Verschiebung am Trägerende und der Biegewinkel an dieser Stelle, der mit Hilfe der Differenzenformel für die erste Ableitung aus den Verschiebungen der beiden Nachbarpunkte berechnet wird, sowie das auf entsprechende Weise mit der Näherung der zweiten Ableitung berechnete Biegemoment an der Einspannung. Schließlich wird in einer Schleife der gesamte Biegemomentenverlauf berechnet, und Durchbiegung und Biegemomentenverlauf werden in ein Graphik-Fenster gezeichnet. Für die oben gelistete Datei (mit nA = 1000, Gleichungssystem mit 1005 Gleichungen) sehen nach der Berechnng das “Command Window” und das Graphik-Fenster so aus:

Um die Werte der Tabelle auf Seite 262 zu berechnen, muss in der oben gelisteten M-Datei jeweils nur der Wer für n_A geändert werden.

Die oben gelisteten M-Dateien sind als S261Band1.m bzw. S261Band2.m zum Download verfügbar. Für die Rechnungen benötigt man zusätzlich das Script gabamp.m.

Homepage

www.D@nkert.de

nkert.de