DDBanner02

Zur Startseite

TM-aktuell

Berechnung einer konischen Welle, Verifizieren der Ergebnisse mit dem RITZschen Verfahren mit MATLAB

Ein MATLAB-Beispiel
wird hier ausführlich
beschrieben

Die Ergebnisse für die Berechnung der Durchbiegung der stückweise konischen Welle (Elastizitätsmodul: E = 210000 N/mm2) nach der Finite-Elemente-Methode sollen mit dem RITZschen Verfahren verifiziert werden.

Die nachfolgend gelistete m-Datei lässt die Berechnung mit maximal 5 Polynom-Ansätzen zu, die jeweils für die gesamte Welle gelten:

v = a1 [- (z/l)2 + z/l] + a2 [-(z/l)3 + (z/l)2] + a3 [(z/l)4 - 2(z/l)3 + (z/l)2]
+ a4 [3(z/l)5 - 7(z/l)4 + 4(z/l)3] + a5 [4(z/l)6 - 9(z/l)5 + 5(z/l)4]  .    

Jede Funktion erfüllt für sich die geometrischen Randbedingungen (l ist die Gesamtlänge der Welle), die beiden letzten Funktionen auch die Bedingung der Momentenfreiheit an den Rändern.

% Ritzsches Verfahren für Biegeträger mit Polynom-Ansätzfunktionen,
% hier: Konische Welle

clear all

% Parameter:
tl1 = 100 ;
tl2 = 250 ;
tl3 = 250 ;
tl4 = 400 ;
tl5 = 100 ;
tl = tl1 + tl2 + tl3 + tl4 + tl5 ;
E  = 210000 ;
d1 = 20 ;
d2 = 50 ;
F1 = 2000 ;
F2 = 4000 ;

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

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  = 110 ;                 % 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

for i=1:nS                 % Durchmesser
   if     (zS(i) < tl1) d(i) = d1 ;
   elseif (zS(i) < tl1+tl2) d(i) = d1 + (d2-d1)/tl2 * (zS(i)-tl1) ;
   elseif (zS(i) < tl1+tl2+tl3) d(i) = d2 ;
   elseif (zS(i) < tl1+tl2+tl3+tl4) d(i) = d2 - (d2-d1)/tl4 * (zS(i)-tl1-tl2-tl3) ;
   else   d(i) = d1 ;
   end
end

EI = E * pi*d.^4/64 ;      % Biegesteifigkeiten

% 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)' ;

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

% Punkte der beiden Einzellasten:
pos2 = n * (tl1+tl2)/tl + 1 ;
pos3 = n * (tl1+tl2+tl3)/tl + 1 ;

for ii = 1:m                                         % Schleife über alle Gleichungen
   b(ii) = F1*vS(pos2,ii) + F2*vS(pos3,ii) ;        % Rechte Seite (Einzellastanteile)
  for jj = 1:m                                       % ii-te Zeile (Koeffizienten A)
   Summe = EI(1)*vddS(1,ii)*vddS(1,jj) + EI(nS)*vddS(nS,ii)*vddS(nS,jj) ;
   faktor = 4 ;                                     % Numerische Integration ...
   for k = 2:n                                      % ... nach Simpsonscher Regel ...
     Summe = Summe + EI(k) * vddS(k,ii) * vddS (k,jj) * faktor ;
     if   (faktor == 4) faktor = 2 ;                % ... für Anteil aus der ...
     else               faktor = 4 ;                % ... Biegesteifigkeit ...
     end ;
   end
   A(ii,jj) = Summe*dz/3 ;                         
  end 
end

% Lösen des Gleichungssystems
ai = A\b ;

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

% Spezielle Werte:
DurchbiegungLinks    = vSchlange (1)
DurchbiegungRechts   = vSchlange (nS)
MaximaleDurchbiegung = max(abs(vSchlange))
pos1 = n * (tl1)/tl + 1 ;
pos4 = n * (tl1+tl2+tl3+tl4)/tl + 1 ;
vPos1 = vSchlange(pos1)
vPos2 = vSchlange(pos2)
vPos3 = vSchlange(pos3)
vPos4 = vSchlange(pos4
)

Die am Ende der Datei farblich hervorgehobenen Zeilen sorgen für die graphische Ausgabe der Biegelinie (Bild unten) und die Ausgabe der Durchbiegung an den markanten Punkten (in das “Command Window”, Bild rechts), dies sind jeweils die Stellen des Übergangs vom konstanten Querschnitt auf einen veränderlichen Querschnitt bzw. umgekehrt.

Die Ergebnisse zeigen, dass die Qualität sich trotz der 5 Ansatzfunktionen in Grenzen hält, als Verifikation der mit der FEM berechneten Ergebnisse sind sie allerdings durchaus brauchbar.

Die oben angegebene m-Datei kann durch Modifikation der am Anfang farblich gekennzeichneten Zeile (m = 5) die Berechnung mit weniger Ansatzfunktionen ausführen. Die nachfolgende Tabelle zeigt, dass sich (wie zu erwarten) die Ergebnisse, die von “unbrauchbar” (eine oder zwei Ansatzfunktionen) bis zu “Abweichungen um etwa 10%” (fünf Ansatzfunktionen) reichen, mit zunehmender Anzahl der Ansatzfunktionen verbessert:

 

RITZsches Verfahren

FEM

m = 1

m = 2

m = 3

m = 4

m = 5

35 Elemente

vPos1 [mm]

1,4250

1,2455

2,2200

2,3731

2,3397

2,8192

vPos2 [mm]

3,7406

3,6351

4,9494

5,0325

5,0470

5,5267

vPos3 [mm]

4,2750

4,5723

5,8808

5,9336

5,9577

6,4245

vPos4 [mm]

1,4250

1,7469

2,9314

2,8067

2,7982

3,1986

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

Homepage TM-aktuell

www.JürgenDankert.de

D

nkert.de