Klein B. FEM: Grundlagen und Anwendungen der Finite-Element-Methode im Maschinen - und Fahrzeugbau
Подождите немного. Документ загружается.
7 Elementkatalog für elastostatische Probleme
128
»
»
»
»
¼
º
«
«
«
«
¬
ª
»
¼
º
«
¬
ª
34
33
22
11
yx
yx
yx
yx
ȟ1ȟ1ȟ1ȟ1
Ș1Ș1Ș1Ș1
4
1
J (7.81)
und für die entsprechende Determinante
>@
Șyxyxȟyxyxyxyx
8
1
det
412323411234342113242413
J
.
(7.82)
Hierbei wurden folgende Abkürzungen benutzt:
xxx yyy
ij i j ij i j
,
.
Somit können nach längerer Rechnung die Ableitungen der Formfunktionen gebildet wer-
den. Für
g
1
ergibt sich dann beispielsweise
Jdet8
Șyȟyy
x
g
324324
1
w
w
bzw. (7.83)
Jdet8
Șxȟxx
y
g
233442
1
w
w
,
wodurch letztlich die
B-Matrix gegeben ist. Dies ermöglicht es nun, auch die Steifigkeits-
matrix
dydxt
t
xy
BEBk
³³
aufzustellen, worin noch die Integration in dx dy ersetzt werden muss durch
dȘdȟdetdydx
J . (7.84)
Nach dem Ersetzen erhält man
K
³³
ddȟdett
1
1
1
1
t
JBEBk . (7.85)
Aus der Ausmultiplikation ergeben sich recht komplizierte Integrale der Form
K[
K[
K[K[
³³
dd
cba
cbacba
1
1
1
1
333
222111
,
7.2 Scheiben-Elemente
129
die effektiv nur noch numerisch gelöst werden können. Auf die numerische Integration fini-
ter Gebiete soll aber zu einem späteren Zeitpunkt eingegangen werden.
7.2.9 Isoparametrische Elemente
Wie zuvor schon erwähnt, muss es das Ziel jeglicher Modellierung sein, die äußere Geo-
metrie eines Bauteils so exakt wie möglich zu erfassen. In der Praxis wird dies mit
ausschließlich gerade berandeten Elemente nicht möglich sein, da technische Konturen in
der Regel nicht immer gerade, sondern viel öfter gekrümmt sind. Um hier bessere
Möglichkeiten der Abbildung zu haben, müssen Elemente mit flexibel anpassbaren Rändern,
so genannte
isoparametrische Elemente, definiert werden. Diese Aufgabenstellung ist nicht
trivial, da bei der Elementbeschreibung über das umschließende Gebiet integriert werden
muss. Man kann sich das Problem aber erleichtern, in dem man zur Randbeschreibung die
gleiche Funktion nimmt, mit der der Ansatz des Elements gebildet wird. Dies ist somit das
entscheidende Merkmal der isoparametrischen Elemente.
Seitens der Nomenklatur gibt es jetzt folgende Zuordnung:
Elementtyp
Zahl der
Zwischenknoten
Polynomgrad der
Seitengeometrie
linear
quadratisch
kubisch
parabolisch
parabolisch
0
1
2
3
4
#
n
lineare Funktion
quadratische Parabel
kubische Parabel
Parabel 4. Ordnung
Parabel 5. Ordnung
#
Parabel n + 1 Ordnung
Bild 7.18: Typisierung der isoparametrischen Elemente
Die Funktionsweise des isoparametrischen Konzeptes wollen wir jetzt an einfach nachvoll-
ziehbaren Elementen studieren.
Im Kapitel 7.2.3 haben wir zuvor schon die Darstellung des
Dreieck-Scheiben-Elements in
Flächenkoordinaten (siehe insbesondere Gl. (7.29)) diskutiert. Bei der Koordinatentransfor-
mation haben wir dabei von der Beziehung
332211
332211
yyyy
,xxxx
]]]
]
]
]
Gebrauch gemacht. Diese Beziehung ist von der Bauform völlig identisch zum Verschie-
bungsansatz (siehe Gl. (7.38))
7 Elementkatalog für elastostatische Probleme
130
,vvvv
,uuuu
332211
332211
]]] ]
]
]
]
]
sodass wir es zuvor schon mit einem linearen, isoparametrischen Element zu tun hatten.
Gleiches gilt für das
Viereck-Scheiben-Element. Wir hatten hier für die Koordinatentransfor-
mation
>@
>@
4321
4321
y11y11y11y11
4
1
y
,x11x11x11x11
4
1
x
K[K[K[K[
K[K[K[K[
gefunden. Nach Umformung erhielten wir für den Verschiebungsansatz
>@
>@
.v11v11v11v11
4
1
,v
,u11u11u11u11
4
1
,u
4321
4321
K[K[K[K[ K[
K[K[K[K[ K[
Verallgemeinert man die vorstehende Erkenntnis, so kann definiert werden:
Isoparametrische Elemente nutzen eine spezielle Transformation mit gleichartigen Funk-
tionen, um ein krummlinig berandetes Gebiete in gerade berandete „Einheitsgebiete“ zu
überführen.
Durch weiteren Vergleich, z. B. zwischen dem
Rechteck-Element in Kapitel 7.2.5 und dem
Viereck-Element in Kapitel 7.2.8, ist als weiterer Vorteil festzustellen, dass bei Benutzung
von natürlichen Koordinaten immer sofort der Zusammenhang zwischen den Knotenver-
schiebungen und dem Verschiebungsfeld im Elementinneren gegeben ist.
Auch führen die Formulierungen stets noch auf analytisch auswertbare Integrale. Insofern ist
es nahe liegend, diese Vorgehensweise auf alle krummlinig berandeten Elemente auszu-
weiten, da dies bei der Netzbildung sehr vorteilhaft anwendbar sind. Wir wollen nun dieses
Konzept an den beiden im umseitigen Bild 7.19
gezeigten quadratischen Scheiben-Ele-
menten kurz darlegen. Der hier eingeführte Seitenmittenknoten ist deshalb notwendig, um
mit einer quadratischen Funktion gekrümmte Ränder näherungsweise nachbilden zu können.
Entsprechend ist das Knotenbild zu erweitern, wenn noch höhergradige Polynome verwendet
werden sollen. Dies kann erforderlich sein, wenn auf einem Elementrand auch Krümmungs-
wechsel erforderlich sind.
Dies darf aber nicht zu dem Missverständnis führen, dass jetzt beliebige Konturen der exak-
ten Körpergeometrie nachgebildet werden können. Soll beispielsweise ein Bohrungsrand
modelliert werden, so wird mit isoparametrischen Elementen der Rand stückchenweise
durch ein Polynom erfasst; der Geometriefehler ist zwar klein, aber dennoch vorhanden. Mit
höhergradigem Polynom reduziert sich dieser Fehler jedoch immer mehr, sodass er praktisch
für die Steifigkeit und die Spannungsverteilung in einem Körper keine Rolle spielt. Diese
7.2 Scheiben-Elemente
131
Geometrieabweichung wird aber allgemein zu den Fehlerquellen der FEM gezählt, welches
den Näherungscharakter der Methode ausmacht.
2
4
6
5
3
1
]
1
= 0,5
]
1
= 0
]
3
= 0,5
]
3
= 0
]
2
= 0,5
2
4
6
5
3
1
]
2
= 0
]
1
= 0
]
3
= 0
(-1, -1) (1, -1)
(1, 1)(-1, 1)
K
[
21
8
4 3
6
5
7
y
x
x
(x, y)
K
[
1
8
4
6
5
7
2
3
]
2
= 0
x
1
x
2
x
1
x
2
x
5
y
1
x
1
x
4
x
2
y
1
y
3
y
6
Bild 7.19: Abbildung krummliniger isoparametrischer Elemente auf ein „Einheitselement“
Der maßgebende quadratische Ansatz für das
Dreieck-Element mit sechs Knoten lautet für
eine Richtung in Flächenkoordinaten:
>
@
ddg ]]]]]]]]]]]] ]
313221332211
t
444121212u.
(7.86)
Beim
Viereck-Element ist es hingegen günstiger, im normierten Schwerpunkt-Koordinaten-
system zu arbeiten. Die verallgemeinerte Ansatzfunktion für die acht Knoten lautet:
u gugugugugugugugu[K,
11 22 33 44 55 66 77 88
(7.87)
mit den speziellen Formfunktionen
7 Elementkatalog für elastostatische Probleme
132
gfüri
gfüri
gfüri
iiiii
ii
ii
1
4
11 1 1234
1
2
11 57
1
2
11 68
2
2
[[ KK[[KK
KK [
[[ K
,,,
,
,.
(7.88)
Sowohl der Ansatz von Gl. (7.86) sowie der mit Gl. (7.87) zu bildende Ansatz erfüllt die im
Kapitel 6 formulierten Steifigkeits-, Starrkörperverschiebungs- und Konstantdehnungsbedin-
gungen, insofern können auch die isoparametrischen Elemente verträglich formuliert wer-
den.
Da die Vorgehensweise beim
Viereck-Element für eine ganze Elementklasse repräsentativ
ist, sollen die weiteren Ausführungen nun auf das
Viereck-Element angepasst werden. Ge-
mäß den vorausgegangenen Darlegungen ist es zweckmäßig, den Verschiebungsansatz fol-
gendermaßen zu spezialisieren:
d
g0
0g
d,G,u
»
¼
º
«
¬
ª
K[
»
¼
º
«
¬
ª
K[
Șȟ,
Șȟ,
v
u
t
t
. (7.89)
Hierbei ist der Knotenverschiebungsvektor eingeführt als
> @
8765432187654321
t
vvvvvvvvuuuuuuuu d
.
Demgemäß ist die Formfunktionsmatrix zu partitionieren in die Zeilenvektoren
t
g . Für die
Geometrie des isoparametrischen Elements wird der gleiche Ansatz gemacht, und zwar
»
¼
º
«
¬
ª
»
¼
º
«
¬
ª
»
¼
º
«
¬
ª
y
x
g0
0g
Șȟ,
Șȟ,
y
x
t
t
mit beispielsweise
>@
321
t
xxx x . (7.90)
Damit kann jetzt die Elementsteifigkeitsmatrix berechnet werden aus
dȘdȟdett=dAt
1
1
1
1
t
t
A
JBEBBEBk
³³³
. (7.91)
Zur Auflösung dieser Gleichung muss aber noch die B-Matrix bestimmt werden. Diese er-
gibt zunächst formal zu
»
¼
º
«
¬
ª
»
»
»
»
»
»
¼
º
«
«
«
«
«
«
¬
ª
w
w
w
w
w
w
w
w
t
t
xy
y
0
0
x
g0
0g
GDB .
7.2 Scheiben-Elemente
133
Da aber die Differenziationen nach [ und K auszuführen sind, müssen gemäß des vorstehen-
den Zusammenhangs die Differenzialoperatoren ersetzt bzw. zuvor geschickt mit einer Zu-
weisungsmatrix zerlegt werden. Dies erfolgt mit der Operation
»
»
»
¼
º
«
«
«
¬
ª
»
»
»
»
»
»
»
»
»
»
¼
º
«
«
«
«
«
«
«
«
«
«
¬
ª
w
w
w
w
w
w
w
w
»
»
»
»
»
»
¼
º
«
«
«
«
«
«
¬
ª
t
t
y
0
x
0
0
y
0
x
0110
1000
0001
g0
0g
B . (7.92)
Wird dies nun berücksichtigt, so folgt weiter
»
»
»
»
»
»
»
»
»
»
¼
º
«
«
«
«
«
«
«
«
«
«
¬
ª
w
w
w
w
w
w
w
w
»
»
»
¼
º
«
«
«
¬
ª
»
»
»
»
»
»
¼
º
«
«
«
«
«
«
¬
ª
»
»
»
¼
º
«
«
«
¬
ª
»
»
»
»
»
»
»
»
»
»
¼
º
«
«
«
«
«
«
«
«
«
«
¬
ª
w
w
w
w
w
w
w
w
»
»
»
¼
º
«
«
«
¬
ª
»
»
»
»
»
»
¼
º
«
«
«
«
«
«
¬
ª
t
t
t
t
1
1
t
t
1
1
Ș
0
ȟ
0
0
Ș
0
ȟ
0110
1000
0001
Ș
0
ȟ
0
0
Ș
0
ȟ
0110
1000
0001
g
g
g
g
J0
0J
g0
0g
J0
0J
B
. (7.93)
Hierin muss noch die Jacobi-Matrix bzw. die invertierte Matrix eingesetzt werden zu
>@
yxJ
»
»
»
»
¼
º
«
«
«
«
¬
ª
wK
w
w[
w
, (7.94)
7 Elementkatalog für elastostatische Probleme
134
»
»
»
»
¼
º
«
«
«
«
¬
ª
w[
w
wK
w
w[
w
wK
w
xx
yy
Jdet
1
1
J mit (7.95)
w[
w
wK
w
wK
w
w[
w
yxyx
det J .
Der Algorithmus ist bis jetzt so allgemein gültig, dass mit g Formfunktionen für eine ganze
Elementklasse (1. oder 2. oder 3. Ordnung etc.) eingesetzt werden können. Insofern ist es
angebracht, die Gl. (7.91) numerisch auszuwerten. Im Vorgriff auf das nachfolgende Kapitel
gilt es also, mit einer Quadraturformel den folgenden Ausdruck auszuwerten:
¦¦
³³
{|
n
1i
t
iii
n
1i
i
1
1
1
1
dettwȘ,ȟFwdȘdȟȘȟ,F JBEB . (7.96)
Mit
ii
,
K
[
sind hier die Stützstellen und mit
w
i
die Wichtungsfaktoren der Quadraturformel
eingeführt worden.
7.2.10
Numerische Integration
Als eine wesentliche Schwierigkeit bei den vorhergehend eingeführten verallgemeinerten
Elementen ist die Integration der Steifigkeitsmatrix zu sehen, die wegen der teils komplizier-
ten Koeffizienten und des unregelmäßig berandeten Gebiets nur numerisch durchgeführt
werden kann. Wir wollen jetzt diese Problematik beispielhaft aufgreifen und für den ein- und
zweidimensionalen Fall kurz diskutieren.
Der eindimensionale Fall dient hier mehr der Anschauung. Unterstellen wir, es sei der im
umseitigen Bild 7.20
skizzierte Funktionsverlauf gegeben, und es sei die Fläche unter dieser
Funktion zu bestimmen.
Der einfachste Weg besteht dann darin, das Integrationsintervall in n äquidistante Abschnitte
i
n
ab
ax
i
(i = 0, 1, ..., n)
zu unterteilen und durch die n + 1 Stützstellen
i
xF
ein Interpolationspolynom n-ten
Grades zu legen und darunter zu integrieren. Für die Ordnung n = 1 führt dies zur Trapez-
regel und für die Ordnung n = 2 zur Simpson‘schen Regel. Die sich somit für unterschied-
liche Ordnungen (n = 1 bis n = 4) ergebenden Formeln werden als
Newton-Cotes-
Quadraturen bezeichnet. In FEM-Programmen werden aber bevorzugt die Gauß‘schen
Quadraturformeln verwendet, da sie in der Regel ein genaueres Ergebnis liefern.
7.2 Scheiben-Elemente
135
Bild 7.20:
Integration einer Funktion mit einem Interpolationspolynom 2. Grades
Der grundsätzliche Unterschied zu den einfachen Integrationsformeln besteht darin, dass bei
den Gauß‘schen Quadraturformeln ein gewichteter Ansatz gemacht wird und die Stützstellen
optimiert sind. Allgemein lautet der Ansatz:
i
n
1i
i
b
a
xFw
2
ab
dxxFI
¦
³
| . (7.97)
Es liegen somit 2 (n + 1) Freiwerte, nämlich n + 1 Stützstellen und n + 1 Gewichtskoeffi-
zienten
w
i
, vor. Demnach kann ein Polynom 2(n + 1)-ten Grades exakt integriert werden.
1
F
0
F
7 Elementkatalog für elastostatische Probleme
136
Die Hauptanwendung der numerischen Integration ist die Bestimmung der eingeschlossenen
Fläche bei isoparametrischen Elementen. Deshalb empfiehlt es sich, von einem normierten
Intervall (-1 d[d 1) auszugehen und das Integral
i
n
1i
i
1
1
FwdFI [ [[
¦
³
(7.98)
zu Grunde zu legen. Dies ist auch insofern vorteilhaft, da für das normierte Intervall bereits
die günstigsten Stützstellen und die Wichtungsfaktoren tabelliert vorliegen. Einen kleinen
Ausschnitt aus einer derartigen Tabelle zeigt Bild 7.21
.
0,0000 0000 0000
0,0000 0000 0000
± 0,5773 5026 9190
± 0,7745 9666 9241
2,0000 0000 0000
1,0000 0000 0000
0,8888 8888 8889
0,5555 5555 5556
n
1.
2.
3.
1-1
1-1
1-1
i
[
i
w
Bild 7.21:
Stützstellen und Wichtungskoeffizienten der eindimensionalen Gauß‘schen Qua-
draturformel (nach /DAN 77/) im Intervall [-1, +1]
Die Festlegung auf ein Einheitsintervall stellt aber keine Beschränkung der Allgemeingültig-
keit dar, da mittels der Transformation
2
ab
2
ab
x
ii
[
auf jedes andere Intervall umgerechnet werden kann.
Von weit größerer Wichtigkeit ist hier aber die Behandlung von Doppelintegralen für die
Integration von Dreieck- und Viereckbereichen. Diese Integrale liegen in der Form
³³
d
cy
b
ax
dydxy,xFI (7.99)
vor. Man erkennt, dass die Näherung durch zweimalige Anwendung der Quadraturregel ent-
steht und deshalb auch auf Flächen ausgedehnt werden kann. Die Anzahl der Integrations-
punkte kann dabei ohne Weiteres in beiden Richtungen unterschiedlich sein.
7.2 Scheiben-Elemente
137
Bei der praktischen Anwendung der Quadraturformel auf Rechteckbereiche empfiehlt es
sich, analog zum Einheitsintervall von einem Einheitsquadrat (s. Bild 7.22
) auszugehen.
x
h
y
d
c
a
b
(-1, -1) (1, -1)
(1, 1)(-1, 1)
[
Bild 7.22: Transformation eines Rechteckbereichs auf ein Einheitsquadrat
Durch die Transformation
2
cd
2
cd
y
,
2
ab
2
ab
x
ii
ii
K
[
(7.100)
kann nunmehr das Integral Gl. (7.99) für ein beliebiges Gebietsintervall angenähert werden
durch
ii
m
1j
n
1i
ij
d
c
b
a
,Fww
4
cdab
dydxy,xFI K[
|
¦¦
³³
(7.101)
oder für gleiche Stützstellen (m = n) im Einheitsquadrat
*)
angenähert werden durch
¦
K[
|
m,n
1i
iii
,FW
4
cdab
I . (7.102)
Für den häufig benutzten Fall mit zwei richtungsabhängigen Stützstellen n = m = 3 weist
Bild 7.23
die so genannte Neun-Punkte-Formel aus, die Polynome bis zur 5. Ordnung exakt
integrieren kann.
*)
Anmerkung: Bei finiten Flächenelementen werden an den Gauß-Punkten in der Regel auch die Element-
spannungen angegeben.