7 Der Arbeits- und Energiebegriff in der Elastostatikuser · 2 q q 6 4 d f 3 f A 1 0 2 2 a l l l l...

F. U. Mathiak 7-1

7 Der Arbeits- und Energiebegriff in der Elastostatik 7.1 Die Arbeit einer Kraft längs eines Verschiebungswe-ges

Abb. 7-1 Arbeit einer Kraft längs eines Verschiebungsweges Für die an einem starren Körper angreifende Kraft F, deren Angriffspunkt sich auf einer Bahnkurve C bewegt (Abb. 7-1), definieren wir die differentielle Arbeit längs des Verschie-bungsweges dr als das Skalarprodukt

α=α=⋅= cosFdr)(cosdddAa rrF(r)rF(r) Gl. 7-1

Die skalare Größe dAa ist das Produkt aus der Kraftkomponente in Wegrichtung, also

αcosF , und dem Verschiebungszuwachs dr, wenn Kraft- und Wegrichtung den Winkel α miteinander einschließen. Der Verschiebungszuwachs dr tangiert dabei an jeder Stelle r die Bahnkurve C. Auf dem endlichen Verschiebungsweg von r1 nach r2 verrichtet die Kraft die Arbeit

7-2 Der Arbeits- und Energiebegriff

rrFr

r

d)(A2

1

a ∫ ⋅= Gl. 7-2

Sonderfall: Hängt die Kraft zFeF = nicht vom Verschiebungswege r ab, und sind Kraft und

Verschiebungsdifferential zdzd er = für den gesamten Verschiebungsweg parallel, dann geht

Gl. 7-2 über in

)zz(FdzFd)(A 12

z

za

2

1

2

1

−==⋅= ∫∫ rrFr

r

Gl. 7-3

Die Arbeit kann sowohl positiv, negativ oder auch Null sein. Die Definition wurde gerade so gewählt, daß bei positiver Arbeit Aa die Kraft F Arbeit verrichtet, während bei negativer Ar-beit Aa Arbeit gegen die Kraft aufgewendet werden muß. Für rF d⊥ ist der differentielle Ar-

beitsanteil dAa gleich Null. Die Arbeit hat die Dimension:

[ ] 2

2

a )Zeit()Länge(MasseA ⋅

= Einheit: JNmskgm 22 ==−

Beispiel: 7-1

Gesucht wird die Arbeit einer linear veränderlichen Streckenlast (Abb. 7-2) an einem vorge-gebenen Verschiebungsweg

( )22 463f)x(w ξ+ξ−ξ= . ( )lx=ξ

Abb. 7-2 Kragträger mit linear veränderlicher Streckenlast

Wir führen die Lösung des Problems auf die Arbeit einer infinitesimalen Kraft dF längs des Verschiebungsweges w(x) zurück. Die Arbeit der differentiellen Kraft dF = q(x) dx am Ver-schiebungswege w(x) ist: )x(wdx)x(q)x(w)x(dFdAa ==

F. U. Mathiak 7-3

Für die Arbeit der gesamten Linienlast q(x) gilt dann: dx)x(w)x(qdAA0x

aa ∫∫=

==l

.

Beachten wir ξ∆+=−

+= qqxqqq)x(q rl

ll

l, dann liefert die Integration

( ) ( ) ⎟⎠⎞

⎜⎝⎛ ∆+=ξξ+ξ−ξξ∆+= ∫ q

4513q

52fd46qq

3fA

1

0

22a ll l

l

und für den Sonderfall Gleichstreckenlast erhalten wir mit 0r qqq == l und 0q =∆

fq4,0f5q2

A 00

a ll

== .

Achtung: Prinzipiell falsch wäre es, im Fall der Gleichstreckenlast mittels der Resultierenden

l0qR = und dem zugehörigen Angriffspunkt 2* l=ξ sowie der dortigen Verschiebung

f4817)(w * =ξ die Arbeit fq354,0f

48q17

)(RwA 00*

a ll

==ξ= zu ermitteln.

7.2 Die Arbeit eines Kräftepaares mit dem Moment M

Abb. 7-3 Arbeit eines Kräftepaares

Abb. 7-4 Momentane Drehachse

Die Arbeit eines Kräftepaares mit dem Moment FrM ×= (Abb. 7-3) leiten wir unmittelbar aus der Definition Gl. 7-1 her. Nach Euler kann nämlich die infinitesimale Lageänderung (Abb. 7-4) eines Punktes P des starren Körpers darstellt werden als die Hintereinanderschal-tung einer für alle Körperpunkte identischen Translation drA und einer Rotation mit dem dif-ferentiellen Drehwinkel dϕ um den Punkt A, also

APA ddd rrr ×+= ϕ Gl. 7-4


Dabei ist A ein beliebiger Punkt des Körpers und rAP der Verbindungsvektor von A nach P. Nach Gl. 7-1 ist die differentielle Arbeit des Kräftepaares:

( ) ( ) ( )[ ]( )[ ] ( ) ( ) ( ) ϕϕϕϕϕ

ϕϕdMdrFdrFrdFaadF

addradrFrrFdrFrdF

21

2A21

⋅=⋅×−=×⋅−=×⋅=−×⋅=×+−×+⋅=−⋅=⋅−+⋅= 1A21a ddddA

Der translatorische Anteil hebt sich offensichtlich heraus, und es verbleibt

ϕϕ d)M( ⋅=adA Gl. 7-5

Dreht sich der Körper mit dem Kräftepaar von 1ϕ bis 2ϕ , so wird die Arbeit

∫ ⋅=2

1

aAϕ

ϕ

ϕϕ d)M( Gl. 7-6

verrichtet.

Beispiel: 7-2

Gesucht wird die Arbeit eines Momentes M an der Verdrehung des Stabendes nach Gl. 7-7,

wenn die Verschiebung ( )22 463f)x(w ξ+ξ−ξ= vorgegeben ist.

Abb. 7-5 Balken mit Endmoment

Wird der Balken durch ein Moment M mit Drehrichtung um die negative y-Achse belastet, dann leistet dieses Moment Arbeit an der Tangentenneigung w', denn es gilt für kleine Ver-formungen w'(x)(x)(x)tan =ϕ≈ϕ und damit

'MwA2

1

a =ϕ⋅= ∫ϕ

ϕ

dM

Differentiation der Biegelinie nach x liefert: ( )l

ll 3

f4)x('w333

f4'w 2 ==→ξ+ξ−ξ=

Für die Arbeit des Momentes erhalten wir dann: l

lf

3M4)('MwAa ==

F. U. Mathiak 7-5

7.3 Das Potential einer Kraft Zur Auswertung des Integrals Gl. 7-2 ist in aller Regel die explizite Angabe der Bahnkurve C erforderlich, da sich mit der Lageänderung des Körpers im allgemeinen auch die Kraft F än-dert. Wir sprechen in diesem Fall von einem Kraftfeld F(r). In einem stationären1 Kraftfeld ist F(r) nur vom Ort r abhängig, in einem instationären Kraftfeld hängt F(r,t) zusätzlich noch von der Zeit t ab.

1

2(a) (b)

r1

r2

dr

F(r)

P

Abb. 7-6 Bewegung einer Kraft auf einer geschlossenen Bahnkurve

Betrachten wir Abb. 7-6, dann ist i.a. )b(21

)a(21 AA −− ≠ . Ist jedoch die Arbeit vom Weg unabhän-

gig, dann hängt sie nur vom Anfangs- und Endpunkt der Bahnkurve ab. Wir sprechen dann von einem konservativen2 Kraftfeld. Wegunabhängigkeit

)b(21

)a(21 AA −− =

oder auch

01

)b(2

2

)a(1

=⋅+⋅ ∫∫ drFdrF

ist gegeben, wenn gilt

0A)C(

a =⋅= ∫ drF Gl. 7-7

Die Arbeit verschwindet demnach längs eines beliebigen geschlossenen Weges. Allgemein kann gezeigt werden, daß für ein Kraftfeld, das Gl. 7-7 genügt, ein Potential U(r) existieren muß, aus dem durch Gradientenbildung

1 von lat. stationarius ›stillstehend‹, ›zum Standort gehörig‹ 2 zu lat. conservare ›bewahren‹, ›erhalten‹


⎟⎟⎠

⎞⎜⎜⎝

⎛∂∂

+∂∂

+∂∂

−=−∇=−= zyx zU

yU

xU)(U)(Ugrad eeerrF Gl. 7-8

das Kraftfeld F selbst gewonnen werden kann. Setzen wir nämlich Gl. 7-8 in Gl. 7-2, dann folgt unter Beachtung von

( )

dUdzzUdy

yUdx

xU

dzdydxzU

yU

xUU zyxzyx

=∂∂

+∂∂

+∂∂

=

++⋅⎟⎟⎠

⎞⎜⎜⎝

⎛∂∂

+∂∂

+∂∂

=⋅∇ eeeeeedr

[ ])(U)(U)(dU)(UA 12

r

r

r

r

r

r21

2

1

2

1

2

1

rrrrdrrdF(r) −−=−=⋅∇−=⋅= ∫∫∫− Gl. 7-9

Die Wegunabhängigkeit eines konservativen Kraftfeldes begründet sich aus dem Sachverhalt, daß die Arbeit allein aus der Potentialdifferenz der Orte r2 und r1 gewonnen werden kann. Die Komponentendarstellung von Gl. 7-8 hinsichtlich einer kartesischen Basis liefert

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

∂∂

−

∂∂

−

∂∂

−

=⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

zUyUxU

FFF

z

y

x

F Gl. 7-10

Die Rotation eines Vektorfeldes (hier des Kraftfeldes F) wird in der Vektoranalysis als Wir-bel des Vektorfeldes bezeichnet und symbolisch, unter Verwendung des Nablaoperators, in der Form FF ×∇=rot geschrieben. Kraftfelder, die ein Potential besitzen, sind demnach wirbelfrei, denn es gilt:

0eee =⎟⎟⎠

⎞⎜⎜⎝

⎛∂∂

∂−

∂∂∂

+⎟⎠⎞

⎜⎝⎛

∂∂∂

−∂∂

∂+⎟⎟

⎠

⎞⎜⎜⎝

⎛∂∂

∂−

∂∂∂

=∇×∇=xy

Uyx

Uzx

Uxz

Uyz

Uzy

UUgradUrot22

z

22

y

22

x

Das Potential U(x,y,z) des Kraftfeldes muß also den folgenden Integrabilitätsbedingungen genügen:

0xy

Uyx

U;0zx

Uxz

U;0yz

Uzy

U 222222

=∂∂

∂−

∂∂∂

=∂∂

∂−

∂∂∂

=∂∂

∂−

∂∂∂ Gl. 7-11

F. U. Mathiak 7-7

7.3.1 Das Potential einer Gewichtskraft Als Beispiel einer Kraft, der ein Potential zugeordnet werden kann, betrachten wir die Ge-wichtskraft G eines schweren Körpers in der Nähe der Erdoberfläche (Abb. 7-7).

Abb. 7-7 Arbeit der Gewichtskraft G Die Gewichtskraft zGeG −= hat in dem gewählten Koordinatensystem nur eine Komponen-

te. Mit dem Ortsvektordifferential zyx dzdydx eeedr ++= erhalten wir zunächst

dzGrdAa −=⋅= dG Gl. 7-12

Integrieren wir diesen Ausdruck längs des Verschiebungsweges von r1 nach r2, also

( ) ( )21

z

z1221 zzGzzGGdzA

2

1

2

1

−=−−=−=⋅= ∫∫−

r

r

drG

dann erhalten wir die Arbeit der Gewichtskraft G, die offensichtlich nur von den z- Koordina-ten der beiden Endpunkte abhängt. Das Potential der Gewichtskraft folgt aus Gl. 7-10 zu

adAGdzdUdzdUG −==→−=−

nach Integration zwischen den beiden Lagen r1 und r2

)zz(GUU 1212 −=− Gl. 7-13

Nehmen wir das Nullniveau bei z2 = 0 an (U2 = 0), dann hat der Körper mit der Gewichtskraft

G bezüglich der Ebene NN die Energie der Lage oder die potentielle Energie


hGU = Gl. 7-14

Die Dimension dieser skalaren Größe ist

[ ] JNmsmkgEinheit)Zeit(

)Länge(MasseU 222

2

==⋅

= −

Die potentielle Energie kann anschaulich gedeutet werden als diejenige Energie, die benötigt

wird, um G vom Nullniveau (z = 0) um z anzuheben. Sie ist > 0, wenn sich der Schwerpunkt

oberhalb des Nullniveaus befindet, Null, wenn der Schwerpunkt im Nullniveau liegt und < 0,

wenn er sich unterhalb desselben befindet.

7.3.2 Das Potential einer Federkraft

Abb. 7-8 Lineare Feder, Kraft-Verformungsdiagramm Als weiteres Beispiel für eine konservative Kraft betrachten wir die äußere Kraft xeF cx= ,

die eine lineare Feder mit der Federkonstanten c aus der ungespannten Lage (x = 0) in die Lage x auslenkt (Abb. 7-8). Nach Gl. 7-2 leistet die Kraft dabei die Arbeit

Fx21cx

21xdxcxd)x(FA 2

x

0x

x

0xa ==== ∫∫

==

Gl. 7-15

F. U. Mathiak 7-9

In Gl. 7-15 wurde xdxr ed = berücksichtigt. Die Federkraft fF leistet als Reaktionskraft dann

die innere Arbeit

Fx21cx

21xd)x(FA 2

x

0xf −=−=−= ∫

=

Gl. 7-16

Zur Berechnung des Potentials der Federkraft beachten wir Gl. 7-10 und erhalten

cxdx

dUcxdx

dUFF fffx =→−=−=≡

Integration liefert:

2f cx

21U = Gl. 7-17

Geometrisch entspricht das Potential der Federkraft Uf der in Abb. 7-8 schraffierten Dreiecks-fläche unterhalb der linearen Kraft-Verschiebungskurve. Auch das Potential Uf ist nur bis auf eine (additive) Konstante festgelegt. Entsprechende Beziehungen lassen sich für eine Drehfeder mit der Federkonstanten cd herlei-ten. Ist ydc eM ϕ= das äußere Moment, das die lineare Feder aus der ungespannten Lage

0=ϕ in die Lage ϕ auslenkt, dann folgt unter Beachtung von ϕd yd eϕ= die dabei vom

äußeren Moment geleistete Arbeit

ϕ=ϕ=ϕϕ=ϕϕ= ∫∫ϕ

=ϕ

ϕ

=ϕ

M21c

21dcd)(MA 2

d0

d0

a Gl. 7-18

Für das Federmoment MM −=f ist dann analog zu Gl. 7-17

2df c

21U ϕ= Gl. 7-19

Hinweis: Zu den Kräften, die sich nicht aus einem Potential ableiten lassen, gehören z.B. die geschwindigkeitsabhängigen Reibungskräfte, die dem Materialgesetz

v)v(f vR −= mit 0)v(f >

genügen. Unter Beachtung von dtdtdt

r vdrd == liefert Gl. 7-7


0dtv)v(fdtv

)v(frAa <−=⋅−=⋅= ∫∫∫ vvdR

eine Arbeit, die immer negativ ist. Wegen 0Aa ≠ lässt sich ein Potential nicht nachweisen.

Da diese Kräfte Arbeit zerstreuen, werden sie auch dissipative Kräfte1 genannt. Zur Berech-nung der Arbeit einer dissipativen Kraft muß deshalb der vollständige Verschiebungszustand des Kraftangriffspunktes bekannt sein.

7.4 Formänderungs- und Ergänzungsenergie für elasti-sche Körper

Wir betrachten in einem ersten Schritt einen elastischen Körper, dessen Ausgangszustand spannungs- und verzerrungsfrei ist. Dieser Körper sei einem einachsigen Spannungs- und Deformationszustand unterworfen. Die einzigen von Null verschiedenen Spannungs- und Verzerrungskomponenten sind dann z.B. σ xx und ε xx . Ist die Spannung σ xx eine allgemeine

Funktion von ε xx , also

)(f xxxx ε=σ Gl. 7-20

dann definieren wir als Differential der spezifischen2 Formänderungsenergie (Abb. 7-9)

xxxx)s( d)(fdW εε= Gl. 7-21

Nach Integration über den gesamten Verzerrungszustand erhalten wir die spezifische Form-änderungsenergie

∫ε

=ε

εε=εxx

xx 0xxxxxx

)s( d)(f)(W Gl. 7-22

Fassen wir umgekehrt die Dehnungen als Funktion der Spannungen auf, also

)(f xx*

xx σ=ε Gl. 7-23

dann definieren wir als Differential der spezifischen Ergänzungsenergie (Abb. 7-9)

1 zu lat. dissipare ›zerstreuen‹, ›verschwenden‹ 2 auf die Volumeneinheit bezogen

F. U. Mathiak 7-11

xxxx**)s( d)(fdW σσ= Gl. 7-24

Abb. 7-9 Spezifische Formänderungs- und Ergänzungsenergie

Nach Integration über den gesamten Spannungszustand erhalten wir die spezifische Ergän-zungsenergie

xx

xx

(s)* *xx xx xx

0

W ( ) f ( )dσ

σ =

σ = σ σ∫ Gl. 7-25

Aus Abb. 7-9 wird auch deutlich, warum *)s(W spezifische Ergänzungsenergie genannt wird. Sie ergänzt offensichtlich die spezifische Formänderungsenergie, die sich im eindimensiona-len Fall geometrisch als die Fläche unterhalb der Kurve )(f xxxx ε=σ interpretieren läßt, zu

einem Rechteck der Größe

(s) (s)*xx xxW W+ = σ ε Gl. 7-26

Die spezifische Formänderungsenergie und auch die Ergänzungsenergie haben Potentialcha-rakter, denn aus Gl. 7-21 und Gl. 7-24 folgen unmittelbar


(s)

xx xxxx(s)*

*xx xx

xx

dW f ( )d

dW f ( )d

= ε = σε

= σ = εσ

Gl. 7-27

oder in Worten:

1. Die spezifische Formänderungsenergie ist das Potential der Spannung. 2. Die spezifische Ergänzungsenergie ist das Potential der Verzerrung.

Abb. 7-10 Linear elastisches Material, Formänderungs- und Ergänzungsenergie

Bei linear elastischem Material sind im isothermen Fall (T = 0) die Maßzahlen für die spezifi-sche Formänderungsenergie und die spezifische Ergänzungsenergie gleich. Hier gilt nämlich nach Hooke: xxxx Eε=σ und damit

(s) (s)* 2 2xx xx xx xx

E 1 1W W2 2E 2

= = ε = σ = σ ε

In Verallgemeinerung auf den dreidimensionalen Fall erhalten wir die spezifische Formände-rungsenergie für Hookesches Material in Matrizenschreibweise

( ) εDεεDεεσ TTT(s)

21

21

21W === Gl. 7-28

Ausrechnen von Gl. 7-28 liefert in kartesischen Koordinaten

( ) ( )⎥⎦⎤

⎢⎣⎡ +++++

−+++= 2

zx2

yz2

xy2

zzyyxx2

zz2

yy2

xx(s) γγγ

21εεε

ν21νεεεGW Gl. 7-29

F. U. Mathiak 7-13

und unter Beachtung von Luε = gilt weiter

u)LDL(uDεε TTT ()21

21W(s) == Gl. 7-30

Die Auswertung der rechten Seite von Gl. 7-30 liefert

⎪⎪⎪

⎭

⎪⎪⎪

⎬

⎫

⎪⎪⎪

⎩

⎪⎪⎪

⎨

⎧

⎥⎥⎦

⎤

⎢⎢⎣

⎡⎟⎠⎞

⎜⎝⎛

∂∂

+∂∂

+⎟⎠⎞

⎜⎝⎛

∂∂

+∂∂

+⎟⎟⎠

⎞⎜⎜⎝

⎛∂∂

+∂∂

+

+⎟⎟⎠

⎞⎜⎜⎝

⎛∂∂

+∂∂

+∂∂

−+⎟

⎠⎞

⎜⎝⎛∂∂

+⎟⎟⎠

⎞⎜⎜⎝

⎛∂∂

+⎟⎠⎞

⎜⎝⎛∂∂

=222

2222

(s)

zx

xz

vz

zv

xv

yu

21

zw

yv

xu

ν21ν

zw

yv

xu

GW Gl. 7-31

Die Formänderungsenergie ergibt sich in jedem Fall durch Integration der spezifischen Form-änderungsenergie über das Gesamtvolumen des betrachteten Körpers.

∫=(V)

(s)dVWW Gl. 7-32

7.5 Formänderungsenergie für den geraden Balken

Wir beschränken uns auf linear elastisches Material. Für die folgenden Untersuchungen wird der isotherme Fall mit T = 0 zugrunde gelegt.

7.5.1 Schiefe Biegung mit Normalkraft Von den Verzerrungen verbleibt nur die Dehnung vywzuxx ′′−′′−′=ε . Alle anderen Verzer-

rungen sind Null. Die spezifische Formänderungsenergie Gl. 7-29 reduziert sich unter Beach-tung von 0=ν und G = E/2 auf

( )

( )22222

22xx

)s(

vyvwyz2wzvyu2wzu2u2E

vywzu2E

2EW

′′+′′′′+′′+′′′−′′′−′=

′′−′′−′=ε=

und mit dxdAdV = erhalten wir bei Bezugnahme auf Hauptzentralachsen


dxdAy)x(vdAyz)x(v)x(w2dAz)x(w

dAy)x(v)x(u2dAz)x(w)x(u2dA)x(u

2EdVWW

0x

I

)A(

22

0

)A(

I

)A(

22

0

)A(

0

)A(

A

)A(

2

)V(

)s(

zzyy

∫ ∫∫∫

∫∫∫

∫=

===

===

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

′′+′′′′+′′+

′′′−′′′−′

==l

3214321321

321321321

Nach Zusammenfassung verbleibt

[ ] dx)x(vEI)x(wEI)x(uEA21W

0x

2zz

2yy

2∫=

′′+′′+′=l

Gl. 7-33

Die spezifische isotherme Formänderungsenergie läßt sich unter Beachtung von

⎥⎥⎦

⎤

⎢⎢⎣

⎡−+=

σ=ε y

I)x(Mz

I)x(M

A)x(N

E1

E zz

z

yy

yxxxx

auch als Funktion der Schnittlasten N(x), My(x) und Mz(x) darstellen:

2

zz

z

yy

y2xx

)s( yI

)x(MzI

)x(MA

)x(NE21

2EW

⎥⎥⎦

⎤

⎢⎢⎣

⎡−+=ε=

Integration über das Stabvolumen führt auf

dxEI

)x(MEI

)x(MEA

)x(N21W

0x zz

2z

yy

2y

2

∫= ⎥

⎥⎦

⎤

⎢⎢⎣

⎡++=

l

Gl. 7-34

Beispiel: 7-3

Für den durch eine Normalkraft F und eine Linienlast q0 belasteten Kragträger ist die Formänderungsenergie zu berechnen. Lösung: Aus Gleichgewicht ergeben sich die Schnittlasten:

N = F = konst.; ( )20y x

2q)x(M −−= l

Mit Gl. 7-34 erhalten wir dann

F. U. Mathiak 7-15

( )yy

520

2

0x

4

yy

20

2

0x yy

2y

2

EI40q

2EAFdxx

4EIq

EAF

21dx

EI(x)M

EA(x)N

21W ll

lll

+=⎥⎥⎦

⎤

⎢⎢⎣

⎡−+=

⎥⎥⎦

⎤

⎢⎢⎣

⎡+= ∫∫

==

Hinweis: Die Verschiebung u eines Dehnstabes infolge Einzellast F am Stabende ist bekannt-

lich EAFu l

= . Die Arbeit der Kraft F an dieser Verschiebung ist nach Gl. 7-15

EA2FFu

21A

2l==

und damit identisch mit der im Körper gespeicherten Formänderungsenergie, die z.B. dazu genutzt werden kann, den Körper bei der Entlastung wieder in seinen Ausgangszustand zu bringen.

7.5.2 Querkraftbeanspruchung

Für den schubelastischen Balken gilt: ( )w21

21

yyxz ′+ω=ϕ=ε . Mit Gl. 7-29 finden wir dann

2y

)S( G21W ϕ=

und nach Integration über das Gesamtvolumen des Stabes

dx)x(GA21W

0x

2y∫

=

ϕ=l

Gl. 7-35

Unter Beachtung des Werkstoffgesetzes für die Querkraft yz GAQ ϕ= können wir mit Gl.

7-35 auch

dxGA

)x(Q21W

0x

2z∫

=

=l

Gl. 7-36

schreiben. Durch formale Erweiterung auf den zweidimensionalen Fall erhalten wir in Erwei-terung von Gl. 7-35

dx)x(GA21dx)x(GA

21W

0x

2z

0x

2y ∫∫

==

ϕ+ϕ=ll

Gl. 7-37

Und entsprechend von Gl. 7-36


dxGA

)x(Q21dx

GA)x(Q

21W

0x

2y

0x

2z ∫∫

==

+=ll

Gl. 7-38

Im Fall des schubstarren Balkens liegt kein Stoffgesetz für die Querkräfte vor. Wir gehen deshalb von den aus Gleichgewichtsbetrachtungen ermittelten Schubspannungen

)z(b)x(I)z(S)x(Q

yy

yzxz =σ

aus. Mit Gl. 7-25 erhalten wir 2

yy

yz2xz

)S()S(

)z(b)x(I)z(S)x(Q

G21

G21W*W

⎥⎥⎦

⎤

⎢⎢⎣

⎡=σ==

und damit

dx)x(A)x(I)x(A)x(Qdz

)z(b)z(S

dzdx)z(b)z(b)x(I)z(S)x(Q

G21dAdx

)z(b)x(I)z(S)x(Q

G21dV

G1

21W

0x2

2z

z

zz

2y

)V(22

2y

2z

)V(22

2y

2z

)V(

2xz

yy

u

o

yyyy

∫ ∫

∫∫∫

= = ⎥⎥⎦

⎤

⎢⎢⎣

⎡

⎟⎟⎠

⎞⎜⎜⎝

⎛=

==σ=

l

Die obige Beziehung können wir noch etwas kompakter schreiben, wenn wir den dimensions-losen Querschnittswert

∫=

=κu

o

2

z

zz

2y

yyz dz

)z(b)z(S

)x(I)x(A)x( Gl. 7-39

einführen. Wir erhalten dann die Formänderungsenergie

∫=

κ=l

0x

2z

z dx)x(A)x(Q)x(

G21W Gl. 7-40

Beispiel 7-1: Gesucht wird die Formänderungsenergie für einen Balken mit Rechteckquer-schnitt der Breite b und der Höhe h.

Lösung: ⎥⎥⎦

⎤

⎢⎢⎣

⎡⎟⎠⎞

⎜⎝⎛−=

22

y hz21

8bh)z(S . Für den Querschnittswert erhalten wir

( ) 56

120hb

hb144dz)z(S

b12bhbhdz

b)z(S

IA 52

52

2/h

2/hz

2y23

2/h

2/hz

2y

yyz 2 ====κ ∫∫

−=−=

und damit

F. U. Mathiak 7-17

∫=

=l

0x

2z dx)x(Q

GA53W

Beispiel 7-2: Gesucht wird die Formänderungsenergie eines Balkens mit Kreisquerschnitt (Ra-dius a).

2323

y

4

yy2

az1a

32)z(S;

4aI;aA

⎥⎥⎦

⎤

⎢⎢⎣

⎡⎟⎠⎞

⎜⎝⎛−=

π=π= ;

2

az1a2)z(b ⎟⎠⎞

⎜⎝⎛−= ; dz)z(bdA =

65a

0z

2525

a

az

2y a

725a

325a

94dz

az1a

94dz

)z(b)z(S

π=π=⎥⎥⎦

⎤

⎢⎢⎣

⎡⎟⎠⎞

⎜⎝⎛−= ∫∫

=−=

910a

725

4a

adA)z(b)z(S

IA 6

24

2

)A(2

2y

yyz 2 =π

⎟⎟⎠

⎞⎜⎜⎝

⎛ π

π==κ ∫

und damit

∫=

=l

0x

2z dx)x(Q

GA95W

Für den Fall der schiefen Biegung mit Qz und Qy werden die oben hergeleiteten Beziehungen sinnvoll erweitert. Wir erhalten entsprechend Gl. 7-40

∫∫==

κ+κ=ll

0x

2y

y0x

2z

z dxGA

)x(Q)x(

21dx

GA)x(Q)x(

21W Gl. 7-41

Bei dünnwandigen Walzprofilen (I-Profil, U-Profil usw.) wird die Querkraft Qz vorwiegend durch den Steg abgetragen. Bezeichnet A die Querschnittsfläche des Gesamtquerschnittes und ASteg die Querschnittsfläche des Steges, dann kann näherungsweise

Stegz A

A=κ Gl. 7-42

gesetzt werden.


7.5.3 Torsion Wir beschränken uns auf den wölbfreien Kreis- bzw. Kreisringquerschnitt. Für den Verschie-bungsvektor gilt

)]x(y);x(z;0[ xx ϑϑ−=u

mit den daraus resultierenden Gleitungen

)x(y)x(z

xxz

xxy

ϑ′=γ

ϑ′−=γ

Abb. 7-11 Torsion eines kreisförmigen Balkens Alle anderen Verzerrungen sind Null. Mit Gl. 7-29 erhalten wir

( ) )x()r

zy(2G

2GW 2

x

2

222xz

2xy

)s( ϑ′

=+=γ+γ=321

und nach Integration über das Stabvolumen

∫∫ ∫==

=

ϑ′=ϑ′=ll

321 0x

2xp

0x

2x

I

)A(

2 dxGI21dx)dAr(G

21W

p

Gl. 7-43

Beachten wir noch das Werkstoffgesetz p

xx GI

M)x()x(D =ϑ′≡ , dann lässt sich die in einem tor-

dierten Stab gespeicherte Formänderungsenergie auch durch das Schnittlastmoment Mx aus-drücken:

F. U. Mathiak 7-19

∫=

=l

0x p

2x dxGI

)x(M21W Gl. 7-44

Liegt eine kombinierte Beanspruchung vor, dann dürfen die Einzelbeanspruchungen zur Ge-samtlösung superponiert werden. Die Addition sämtlicher Formänderungsanteile liefert:

dxGI

)x(MGA

)x(QGA

)x(QEI

)x(MEI

)x(MEA

)x(N21W

0x p

2x

2y

y

2z

zzz

2z

yy

2y

2

∫= ⎥

⎥⎦

⎤

⎢⎢⎣

⎡+κ+κ+++=

l

Normalkraft Biegung Querkraft Torsion

∫=

′=l

0x

2 dx)x(uEA21W

∫

∫

=

=

′′+

′′=

l

l

0x

2zz

0x

2yy

dx)x(vEI21

dx)x(wEI21W

dx)x(GA21

dx)x(GA21W

0x

2z

0x

2y

∫

∫

=

=

ϕ+

ϕ=

l

l

(schubelastischer Balken)

∫=

ϑ′=l

0x

2xp dxGI

21W

∫=

=l

0x

2

dxEA

)x(N21W

∫

∫

=

=

+

=

l

l

0x zz

2z

0x yy

2y

dxEI

)x(M21

dxEI

)x(M21W

∫

∫

=

=

κ+

κ=

l

l

0x

2y

y

0x

2z

z

dxGA

)x(Q21

dxGA

)x(Q

21W

(schubstarrer Balken)

∫=

=l

0x p

2x dx

GI)x(M

21W

Tabelle 7-1 Formänderungsenergien für den geraden Stab bei linear elastischem Materialverhalten

Hinweis: Der Anteil der Formänderungsenergie aus Querkraft ist im Vergleich zu den übrigen Beanspruchungen von untergeordneter Bedeutung und wird deshalb in praktischen Berech-nungen oftmals vernachlässigt.

7.6 Die isotherme Formänderungsenergie für die Scheibe

Die Scheibenebene liege parallel zur x-y-Ebene. Es gilt der ebene Spannungszustand mit dem Verschiebungsvektor

)v,u(T =u und den daraus abgeleiteten Verzerrungen


[ ] ⎥⎦

⎤⎢⎣

⎡∂∂

+∂∂

∂∂

∂∂

=γεε=xv

yu,

yv,

xu,, xyyyxx

Tε

sowie der Materialmatrix

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

−−

=

2ν100

01ν

0ν1

ν1E

2ESD

Für die spezifische Formänderungsenergie erhalten wir dann

⎥⎦⎤

⎢⎣⎡ γ

ν−+ν++

−== 2

xyyyxx2

yy2

xx2(s)

21εε2εε

)ν1(2E

21W εDε ES

T Gl. 7-45

und unter Beachtung von Luε =

⎥⎥⎦

⎤

⎢⎢⎣

⎡⎟⎟⎠

⎞⎜⎜⎝

⎛∂∂

+∂∂ν−

+∂∂

∂∂

ν+⎟⎟⎠

⎞⎜⎜⎝

⎛∂∂

+⎟⎠⎞

⎜⎝⎛∂∂

−=

=

222

2

(s)

xv

yu

21

yv

xu2

yv

xu

)ν1(2E

21W (Lu)D(Lu) ES

T

Gl. 7-46

Die Formänderungsenergie der gesamten Körpers ist dann (Scheibendicke h)

( )∫∫∫∫ ==A(A)

(s) dA2hdAWhW (Lu)D(Lu) ES

T Gl. 7-47

7.7 Formänderungsenergie für die schubstarre Platte

Die Plattenmittelfläche liege in der x-y-Ebene. Die Verschiebung

⎟⎟⎠

⎞⎜⎜⎝

⎛∂∂

−∂∂

−==ywz,

xwz)v,u(Tu

eines Punktes P mit dem Abstand z von der Plattenmittelfläche läßt sich durch die Änderung des Verschiebungsfeldes w(x,y) ausdrücken. Daraus ergeben sich die Verzerrungen

F. U. Mathiak 7-21

[ ] zzyx

w2,yw,

xw,, T

2

2

2

2

2

xyyyxxT κε =⎥

⎦

⎤⎢⎣

⎡∂∂

∂−

∂∂

−∂∂

−=γεε=

Mit der Materialmatrix für den ebenen Spannungszustand

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

ν−ν

ν

ν−=

2100

0101

1E

2D

folgt für die Spannungen Dεσ = und damit die spezifische Formänderungsenergie

DκκDεεεσ T2

TT)s(

2z

21

21W === Gl. 7-48

Ausrechnen liefert

( )⎪⎭

⎪⎬⎫

⎪⎩

⎪⎨⎧

⎥⎦

⎤⎢⎣

⎡⎟⎟⎠

⎞⎜⎜⎝

⎛∂∂

∂−

∂∂

∂∂

ν−−∆ν−

=yx

wyw

xw)1(2w

)1(2EzW

2

2

2

2

22

2

2)s( Gl. 7-49

In jedem Fall ist die Formänderungsenergie der gesamten Platte

dA21dA

12h

21dVWW

)A(

T

)A(

T3

)V(

)s( ∫∫∫∫∫ === κDκDκκ PL Gl. 7-50

wobei

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

ν−ν

ν

ν−=

2100

0101

)1(12Eh

2

3

PLD

die Materialmatrix der Plattenbiegung bezeichnet.

F.U. Mathiak 8-1

8 Das Prinzip der virtuellen Verrü-ckung In der Statik bilden die beiden Axiome

− Kraftgleichgewicht und − Momentengleichgewicht

die Grundlage zur Berechnung mechanischer Systeme. Neben diesen Grundgesetzen wer-den außer dem Befreiungsprinzip und dem Schnittprinzip bei Systemen starrer Körper die kinematischen Beziehungen und bei deformierbaren Körpern zusätzlich geeignete Materialgesetze benötigt. Die Vorgehensweise bei der Lösung mechanischer Probleme ist dabei immer dieselbe. Systeme starrer Körper:

1. Befreiung des Systems von den Auflagern, Anbringen der Auflagerreaktionslasten nach dem Befreiungsprinzip,

2. Zerschneiden des Systems, Anbringen der Schnittlasten nach dem Schnittprinzip, 3. Anwendung der Gleichgewichtsbedingungen auf jedes freigeschnittene Teilsystem, 4. Formulierung der kinematischen Beziehungen für das Starrkörpersystem, 5. Berechnung der Auflager- und Schnittkräfte.

Deformierbare Körper: 1. Herausschneiden eines Volumenelements, 2. Anwendung der Gleichgewichtsbedingungen auf das Element, 3. Aufstellen der kinematischen Beziehungen, also den Beziehungen zwischen den

Verschiebungen und den Verzerrungen, 4. Wahl eines Materialgesetzes, 5. Elimination der Spannungen führt auf ein System partieller Differentialgleichungen

für die Verschiebungen, den Lamé-Navierschen Verschiebungsdifferentialglei-chungen, Elimination der Verschiebungen führt auf die Beltrami-Michellschen Spannungsdifferentialgleichungen,

6. Lösung des Gleichungssystems unter Berücksichtigung der Randbedingungen lie-fert die vollständige Deformation des Körpers.

In vielen Fällen kann sich diese Vorgehensweise jedoch als recht umständlich erweisen. Bei Systemen starrer Körper interessieren häufig die Schnittkräfte gar nicht. In diesen

8-2 Das Prinzip der virtuellen Verrückung

Fällen müssen, insbesondere bei großen Systemen, viele überflüssige Unbekannte mit berechnet werden. Bei deformierbaren Systemen ist oft der Deformationszustand des gesamten Systems nicht gefragt. Nicht selten reicht es aus, Zustandsgrößen an einzelnen Punkten zu berechnen (Verschiebungen, Auflagerkräfte bei statisch unbestimmten Sys-temen). In diesen Fällen wird bei der Integration meist komplizierter Differentialglei-chungen viel überflüssige Arbeit geleistet. Für die praktisch interessanten Fälle ist eine analytische Integration der anfallenden partiellen DGLn meist ohnehin nicht möglich. Für solche Fälle benötigen wir Näherungsverfahren, mit denen wir bei möglichst geringem Rechenaufwand eine gute Annäherung an die theoretisch exakte Lösung erreichen. Dazu dienen die Energieprinzipien der Mechanik, die folgende Vorteile bieten: 1. Angewandt auf starre Körper oder Systeme von starren Körpern erlauben sie eine

schnelle Herleitung der Gleichgewichtsbedingungen ohne Kenntnis der Schnittgrö-ßen.

2. Angewandt auf deformierbare Körper ermöglichen sie die Ermittlung von Kräften oder Verschiebungen an einzelnen Körperpunkten ohne Kenntnis der Lösung der Grundgleichungen (Sätze von Castigliano).

3. Außerdem bilden sie die Grundlage zur Herleitung von Näherungsverfahren, wie dem Verfahren von Ritz und der Methode der Finiten Elemente (FEM), die sehr eng mit dem Ritzschen Verfahren verbunden ist.

Das Prinzip der virtuellen1 Verrückung ist eine den Gleichgewichtsbedingungen äquivalente Energieaussage, das sowohl für starre als auch für deformierbare Kör-per formuliert werden kann.

Abb. 8-1 Beispiel einer virtuellen Verrückung am beidseitig drehbar gelagerter Balken

1 virtualis zu lat. virtus ›Tüchtigkeit‹, ›Mannhaftigkeit‹, fachsprachlich für nicht wirklich; scheinbar; ge-dacht.

F.U. Mathiak 8-3

Definition: Eine virtuelle Verrückung uδ ist eine dem aktuellen Deformationszustand zusätzlich überlagerte Verschiebung (Variation des Verschiebungszustandes) mit den folgenden Eigenschaften: 1. uδ ist geometrisch möglich, d.h. bei der virtuellen Verrückung wird der Zusammen-

hang des Körpers gewahrt. uδ ist verträglich (kompatibel) mit den Lagerungsbedin-gungen.

2. Die Verrückung uδ ist differentiell klein. Damit können in allen Rechnungen Terme

von höherer Ordnung als klein gestrichen werden. 3. Sämtliche inneren und äußeren Kraftgrößen werden bei der Durchführung der Variati-

on uδ konstant gehalten, d.h. nicht variiert. 4. Die virtuelle Verrückung ist eine gedachte Verrückung bei festgehaltener Zeit. Dabei

ist uninteressant, in welcher Zeit wir uns diese virtuelle Verrückung entstanden den-ken.

Das Prinzip der virtuellen Verrückungen für elastische Körper besagt (ohne Herleitung)

ia δAδA = Gl. 8-1

oder in Worten: Befindet sich ein Körper im Gleichgewicht, dann ist bei einer virtuellen Verrückung des Körpers die Arbeit der äußeren Kräfte gleich der Arbeit der inneren Kräfte. Werden bei der virtuellen Verrückung die Lagerungsbedingungen des Systems berück-

sichtigt, dann kann statt aδA auch (e)aδA (Arbeit der eingeprägten Kräfte1) geschrieben

werden, da dann die Reaktionskräfte bei der virtuellen Verrückung keine Arbeit leisten. Es kann gezeigt werden (hier ohne Beweis), daß die Variation der inneren Arbeit iden-tisch ist mit der Variation der Formänderungsenergie, also δWδAi = , so dass wir für Gl.

8-1 auch δWδAa = oder umgeformt

( ) 0δΠAWδδAδW aa =≡−=− Gl. 8-2

schreiben können. Der Ausdruck

1 Die äußeren Kräfte setzen sich bekanntlich aus eingeprägten Kräften und Reaktionskräften zusammen


)dOdVρ(WAWΠO(V)(V)

n

1kk

m

1jjja ∫∫∑∑

↓↓

=

↓

=

↓

+++−=−= usufMuF Ok ϕ

wird elastisches Potential genannt, und Gl. 8-2 heißt Satz vom Extremum des elasti-schen Potentials. In der obigen Gleichung wird im Ausdruck für die äußere Arbeit durch die aufgesetzten Pfeile angedeutet, daß allein die Verrückungsgrößen zu variieren sind. Das ist auch dann der Fall, wenn die Belastungen, z.B. bei Federkraftbelastungen, von den Verschiebungen abhängen. In Worten besagt Gl. 8-2: Von allen denkbaren Verschiebungszuständen eines elastischen Körpers ist Derjenige der wirklich Eintretende, für den die Energiegröße Π einen stationären Wert annimmt.

Abb. 8-2 Biegeträger mit Querbelastung

Ein Mittels des Prinzips der virtuellen Verrückung bereitgestelltes funktional unbestimm-tes Problem kann in eine Verschiebungsdifferentialgleichung übergeführt werden. Dazu betrachten wir den Biegeträger nach Abb. 8-2. Im Sinne des Prinzips muß der Ausdruck

dx)x(w)x(q)x(w2

EIw

0x

2yy∫=

⎥⎦

⎤⎢⎣

⎡−′′>=<Π

l

Gl. 8-3

extremal gemacht werden. Das erfordert

[ ] 0dx)x(w)x(q)x(w)x(wEI0x

yy =δ−′′δ′′=Πδ ∫=

l

Gl. 8-4

Die ansonsten beliebigen Verschiebungsvariationen δw(x) müssen an den Rändern die Bedingungen

0)(w)0(w =δ=δ l Gl. 8-5

erfüllen. Zur Einarbeitung der Randbedingungen wird Gl. 8-4 zweimal partiell integriert

F.U. Mathiak 8-5

[ ] ( ) ( ) 0dx)x(w)x(q)x(wEI)x(w)x(wEI)x(w)x(wEI0x

yy0

yy0yy =δ⎥⎦⎤

⎢⎣⎡ −″′′+⎥⎦

⎤⎢⎣⎡ δ′′′−′δ′′=Πδ ∫

=

lll

Die Beachtung von Gl. 8-5 und der Momentenrandbedingungen verbleibt

( ) 0dx)x(w)x(q)x(wEI0x

yy =δ⎥⎦⎤

⎢⎣⎡ −″′′∫

=

l

Gl. 8-6

Wegen der Beliebigkeit von δw(x) ist Gl. 8-6 nur erfüllt, wenn die Differentialgleichung

( ) )x(q)x(wEIyy =″′′ Gl. 8-7

gilt.

Beispiel 8-1: Es ist die Stabendverschiebung u eines geraden zylindrischen Stabes mit konstanter Stab-querschnittsfläche A gesucht. Der Stab wird am Ende durch eine Einzellast F belastet.

Abb. 8-3 Dehnstab mit Einzelkraft

Lösung: Von den Verzerrungen verbleibt nur die Dehnung .konstu(x)uε xx ==′=l

Formänderungsenergie: l

l 2

0x

2 u2

EAdxuEA21W =′= ∫

=

Äußere Arbeit: FuAa =

Elastisches Potential: Fuu2

EAAWΠ(u)Π2

a −=−==l

Das elastische Potential ist eine quadratische Funktion der Stabendverschiebung u. Damit liegt lediglich ein parametrisch unbestimmtes Problem vor, das mit den klassischen Me-thoden der Extremwertberechnung gelöst werden kann.


Die erste Ableitung nach u liefert: 0δuFuEAδuuΠδΠ =⎟

⎠⎞

⎜⎝⎛ −=

∂∂

=l

. Wegen der Belie-

bigkeit von δu ist somit die gesuchte Stabendverschiebung: EAFu l

= . Abb. 8-4 zeigt das

Potential Π als quadratische Funktion von u.

Abb. 8-4 Das Potential Π für den Dehnstab mit Einzelkraft

An der Stelle EAFu l

= nimmt EA2

FFuu2

EAAWΠ(u)22

al

l−=−=−= einen Extremwert

an. Das Vorzeichen der zweiten Variation ( ) 0uEAuu

22 >δ=δΠδ∂∂

=Πδl

entscheidet

von welcher Art der Extremwert ist. Hier liegt offensichtlich ein Minimum vor.

F. U. Mathiak 9-1

9 Näherungsverfahren Aufgrund der komplexen Struktur der allgemeinen Grundgleichungen der Elastizitätstheorie

sind analytische Lösungen für die Deformationen und Kraftgrößen eines Körpers nur in weni-

gen Spezialfällen möglich. Aber gerade die in der Praxis auftretenden Probleme gehören nicht

zu diesen Spezialfällen. Um hier überhaupt zu Lösungen zu gelangen, ist der Berechnungsin-

genieur auf Näherungslösungen angewiesen. Oft ist es in der Praxis auch gar nicht erforder-

lich, die im Sinne der Theorie exakte Lösung zu kennen. Vielmehr reicht es häufig aus, für

die das Problem beschreibenden Größen befriedigende Näherungslösungen anzugeben.

9.1 Das Verfahren von Ritz

Ein Weg zur Gewinnung solcher Näherungslösungen ist das Verfahren von Ritz1. Dieses Ver-fahren basiert auf dem Prinzip der virtuellen Verrückungen für elastische Körper. Es setzt voraus, dass das betrachtete statische Problem als Lösung einer Variationsaufgabe aufgefasst werden kann und ein Potential existiert, das es zu minimieren gilt. Zur Vorstellung des Ver-fahrens betrachten wir das System nach Abb. 9-1.

Abb. 9-1 Träger auf zwei Stützen mit konstanter Querlast

1 Walter Ritz, schweizer. Mathematiker und Physiker, 1878-1909

9-2 Näherungslösungen

Für das zu ermittelnde Verschiebungsfeld w(x) wird ein Ansatz der Form

(x)wc(x)w k

n

1kk∑

=

=ˆ Gl. 9-1

gemacht. Dabei sind die kw bekannte linear unabhängige Funktionen von x. Sie sollen belie-

big wählbar sein, bis auf folgende Einschränkungen: Die Funktionen (x)wk müssen mit den kinematischen Lagerungsbedingungen verträg-

lich sein, sie müssen also mindestens die geometrischen Randbedingungen erfüllen.

Abb. 9-2 Geometrische und dynamische Randbedingungen

Hinweis: Wir hätten das Verschiebungsfeld w(x) auch mit der Näherung

(x)wc(x)w(x)w k

n

1kk0 ∑

=

+=

ansetzen können. Die Funktion )x(w0 erfüllt dann sämtliche inhomogenen Randwerte in den

Verschiebungen. Die wk(x) dürfen dann die Randwerte, die durch w0(x) gegeben sind, nicht mehr stören, sie müssen also an denjenigen Stellen, an denen die Randbedingungen gegeben sind, den Wert Null haben. Damit erfüllt )x(w für beliebige ck die Randbedingungen. Im An-

satz nach Gl. 9-1 kommt w0(x) nicht vor, darum müssen diese wk(x), je für sich, die Randwer-te erfüllen. Ansatzfunktionen kw , die den geometrischen Randbedingungen genügen, heißen zulässige

Funktionen. Genügen die Ansatzfunktionen neben den zwingend erforderlichen geometri-schen Randbedingungen auch noch den dynamischen1 Randbedingungen, so spricht man von Vergleichsfunktionen. Die geometrischen Randbedingungen werden auch als wesentliche Randbedingungen und die dynamischen Randbedingungen als natürlichen Randbedingungen bezeichnet. 1 griech. δυναµιζ = Kraft, Stärke

F. U. Mathiak 9-3

Zur Unterscheidung zwischen geometrischen und dynamischen Randbedingungen betrachten wir Abb. 9-2. Hinweis: In der Regel werden die Ansatzfunktionen aus der Erfahrung heraus gewählt. All-gemein kann hinsichtlich der Güte der mit Gl. 9-1 erzielten Approximation gesagt werden, dass diese um so besser ist, je vollständiger mit den Ansatzfunktionen wk(x) die Randbedin-gungen des Problems erfüllt werden. Für den Träger mit Querlast q(x) nach Abb. 9-1 lautet das elastische Potential:

∫∫==

−=ll

0x

'2'

0x

yy q(x)w(x)dx(x)dxw2

EIΠ Gl. 9-2

Mit dem Ritz-Ansatz Gl. 9-1 erhalten wir einen Näherungswert

∫∫==

−=ll

0x

'2'

0x

yy (x)dxwq(x)(x)dxw2

EIΠ Gl. 9-3

und unter Beachtung von ∑=

=n

1k

''kk

'' (x)wc(x)w folgt:

∫ ∑∑= == ⎪⎭

⎪⎬⎫

⎪⎩

⎪⎨⎧

−⎥⎦

⎤⎢⎣

⎡=

l

0x

n

1kkk

2n

1k

''kk

yy dx(x)wcq(x)(x)wc2

EIΠ = Extremum

Da die Ansatzfunktionen bekannte Funktionen sein sollen, ist

n),,1(k)(cΠΠ k K==

nur noch eine Funktion der unbekannten Koeffizienten kc . Damit ist das funktionalkinema-

tisch unbestimmte Problem zur Ermittlung der Durchbiegung w(x) auf ein algebraisches Problem zur Ermittlung der n Unbekannten Koeffizienten kc zurückgeführt. Das Problem

nach Gl. 9-3 ist damit wesentlich einfacher zu lösen als Gl. 9-2, in der die gesuchte Funktion

w(x) unter dem Integral steht. Da $Π ein Extremum des Variationsintegrals ist, muss gelten:

0δccΠΠδ j

n

1j j

=∂∂

= ∑=

und wegen der Beliebigkeit der jcδ folgt

,n)1(j0cΠ

j

==∂∂ Gl. 9-4


Beachten wir in Gl. 9-4 den Wert für Π nach Gl. 9-3, dann erhalten wir

∫∫∑∫ ∑==== =

−⎥⎦

⎤⎢⎣

⎡=⎥

⎦

⎤⎢⎣

⎡−⎟

⎠

⎞⎜⎝

⎛==

∂∂ lll

0xj

''k

''j

0xyy

n

1kk

0xj

''j

n

1k

''kkyy

j

dxq(x)wdxwwEIcdxq(x)wwwcEI0cΠ

Mit den Abkürzungen

∫

∫

=

=

=ϕ

=ϕ

l

l

0xj0j

''k

0x

''jyyjk

dxq(x)w

dxwwEI Gl. 9-5

können wir auch

n),,1(j0c 0jjk

n

1kk K==ϕ−ϕ∑

=

Gl. 9-6

schreiben. Das sind n lineare Gleichungen zur Bestimmung der n unbekannten Koeffizienten

kc . Fassen wir die Werte jkϕ in der symmetrischen Matrix Φ, die unbekannten Ritzparameter

ck im Vektor c und die Lastanteile 0jϕ im Vektor der rechten Seite b zusammen, also

⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢

⎣

⎡

ϕϕ

ϕϕϕϕϕϕ

==

nn1n

n22221

n11211

LL

MMMM

L

L

TΦΦ ;

⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢

⎣

⎡

=

n

2

1

c

cc

Mc ;

⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢

⎣

⎡

ϕ

ϕϕ

=

0n

20

10

Mb

dann kann das Gleichungssystem Gl. 9-6 in der Form

bcΦ =⋅ Gl. 9-7

geschrieben werden.

Beispiel 9-1: Für die Verschiebung des System in Abb. 9-2 machen wir nun konkret den eingliedrigen Ritz-Ansatz

(x)wc(x)wc(x)w 11

n

1kkk == ∑

=

ˆ Gl. 9-8

wobei

F. U. Mathiak 9-5

l

πxsin(x)w1 = Gl. 9-9

gewählt wird. Dieser Ansatz erfüllt mit 0)0(w1 = und 0)(w1 =l die geometrischen und unter

Beachtung von ( )l

lπxsinπ(x)w 2''

1 −= dann mit My(x = 0) = My(x = l ) = 0 auch die dynami-

schen Randbedingungen. Der Ansatz entspricht damit einer Vergleichsfunktion, von der eine

gute Näherungslösung für die Durchbiegung erwartet werden kann. Das Potential Π ist dann

dxwcq(x)wc2

EIdxwq(x)w

2EI

Π0

112''

121

yy

0

2''yy ∫∫ ⎥⎦

⎤⎢⎣

⎡−=⎥

⎦

⎤⎢⎣

⎡−=

ll

= Extremum

Die Extremaleigenschaften dieses Funktionals sind gleichbedeutend mit

1011110

0'2'

10

yy11

cdxwqdxwEIc0cΠ

ϕ−ϕ=−==∂∂

∫∫ll

Gl. 9-10

mit

⎟⎟⎠

⎞⎜⎜⎝

⎛===ϕ

⎟⎟⎠

⎞⎜⎜⎝

⎛=⎟⎟

⎠

⎞⎜⎜⎝

⎛==ϕ

∫∫

∫∫

πq2dxπxsinqdxwq

π2πEIdxπxsinπEIdxwEI

0o

010

010

3

yy0

24

yy2''

10

yy11

l

l

lllll

ll

Aufgelöst nach c1 erhalten wir

yy

40

yy

40

53

yy

0

11

101 EI

q01307,0

EIq

π4

π2πEI

πq2

cll

l

l

==

⎟⎠⎞

⎜⎝⎛

=ϕϕ

= .

Durch Differentiationsprozesse an Gl. 9-8 erhalten wir die Biegemomente

ll

l

πxsinq12901,0πxsinlqπ4wEIM 2

02

03''

yyy ==−=

und die Querkräfte

ll

ll

πxcosq40528,0πxcosqπ4

dxMd

Q 002y

z ===

Ein Vergleich dieser Ergebnisse mit der exakten Lösung und einer weiteren Lösung, die als Ritz-Ansatz lediglich auf einer zulässigen Funktion basiert, zeigt mit l/x=ξ die folgende

Tabelle:


Vergleich Exakte Lösung Vergleichsfunktion Zulässige Funktion

[ ]32

yy

40 21EI24

q)(w ξ+ξ−ξ=ξl πξ=ξ sinc)(w 1 )1(c)(w 1 ξ−ξ=ξ

l0

max,z

qQ 5,0

21= 0,40528 (-18,9%) 0 (-100%)

20

maxy

qM

l 125,0

81= 0,12901 (3,2%) 0,0833

= konst. (-33,3%)

40

maxyy

qwEIl

01302,0384

5= 0,01307 (0,38%) 0,01043 (-19,9%)

Tabelle 9-1 Vergleich zwischen exakter Lösung und Näherungslösung nach Ritz

Es ist ersichtlich, dass die Approximation derjenigen Größe, für die der Ritz-Ansatz gewählt wurde, hier also für die Durchbiegung w(x), besser ist als die hieraus durch Differentiati-onsprozesse abgeleiteten Größen Biegemomente My und Querkräfte Qz. Diese Aussage gilt ganz allgemein. Hinweis: Die im elastischen Potential auftretenden Integrale können recht kompliziert wer-den, so dass in diesen Fällen auf ein numerisches Integrationsverfahren zurückgegriffen wird. Um die Leistungsfähigkeit des Ritz-Verfahrens erneut zu dokumentieren, betrachten wir als weiteres Beispiel einen Biegestab unter Längskraft nach Abb. 9-3. Um hier zu einer Anwen-dung des Verfahrens zu kommen, muss die Arbeit der Kraft P in das Gleichungssystem ein-gebaut werden. Dazu betrachten wir das System in ausgelenkter Lage.

Abb. 9-3 Träger mit Längs- und Querkraftbeanspruchung, Instabilitätsproblem

Eine erste Näherung für die Axialverschiebung des Lastangriffspunktes der Last P erhalten wir unter der Voraussetzung dehnungsfreier Stabachsverformung mit

[ ] ∫∫∫∫===

≈−+=−=ζ=lll

0x

2'

0x

2'

0x

dxw21dx1w1dx)(dsdζ

F. U. Mathiak 9-7

Die Arbeit der äußeren Belastung im elastischen Potential ist um den Anteil

∫=

=l

0x

2 dxw'2PPζ

zu erweitern, also

0δΠExtremumq(x)w(x)dx(x)dxw'2Pdx(x)w

2EI

Π0x0x

22''

0x

yy =⇒=−−= ∫∫∫===

lll

Gl. 9-11

Die Behandlung von Gl. 9-11 mit einem Ritzansatz (x)wc(x)ww(x)n

1kkk∑

=

=≈ ˆ entsprechend

Gl. 9-1 führt zu

∫ ∑∑∑= === ⎪⎭

⎪⎬⎫

⎪⎩

⎪⎨⎧

−⎥⎦

⎤⎢⎣

⎡−⎥

⎦

⎤⎢⎣

⎡=><≈><

l

L0x

n

1kkk

2n

1k

'kk

2n

1k

''kk

yyn1 dx(x)wcq(x)(x)wc

2P(x)wc

2EI

ccΠw(x)Π

Da >< n1 ccΠ L Extremwertcharakter hat, ist 0cc

ˆccΠ j

jn1 =δ

∂Π∂

>=<δ L . Da im Allge-

meinen 0c j ≠δ gilt, muss 0c

ˆ

j

=∂Π∂ (j = 1,n) erfüllt sein. Wir benötigen also die Ableitungen

jc

ˆ

∂Π∂ . Die Differentiation kann unter dem Integral vorgenommen werden:

44 344 2144 344 21444 3444 21

lll

l

l

0jjkjk

0xj

n

1k 0x

'j

'kk

n

1k 0x

''j

''kyyk

0xj

n

1k

'j

'kk

n

1k

''j

''kyyk

0xj

'j

n

1k

'kk

''j

n

1k

''kkyy

j

(x)dxq(x)w(x)dx(x)wwcPdx(x)(x)wwEIc

dx(x)q(x)w(x)(x)wwcP(x)(x)wwEIc

dx(x)q(x)w(x)w(x)wcP(x)w(x)wcEIcΠ

ϕ

==

ψ

==

ϕ=

=

= ==

= ==

∫∑ ∫∑ ∫

∫ ∑∑

∫ ∑∑

−−=

⎪⎭

⎪⎬⎫

⎪⎩

⎪⎨⎧

−⎥⎦

⎤⎢⎣

⎡−⎥

⎦

⎤⎢⎣

⎡=

⎪⎭

⎪⎬⎫

⎪⎩

⎪⎨⎧

−⎥⎦

⎤⎢⎣

⎡−⎥

⎦

⎤⎢⎣

⎡=

∂∂

Damit erhalten wir das lineare Gleichungssystem

n),1(j0)Pψ(c0cΠ

0jjkjk

n

1kk

j

==ϕ−−ϕ==∂∂ ∑

=

Gl. 9-12

In Gl. 9-12 wurde zur Abkürzung


∫∫∫===

=ϕ′′=ψ′′′′=ϕlll

0xj0jk

0xjyyjkk

0xjyyjk dxq(x)w;dxwwEI;dxwwEI Gl. 9-13

gesetzt. Bei einem zweigliedrigen Ritzansatz, also )x(wc)x(wc(x)wc(x)w 2211

2

1kkk +== ∑

=

hätten wir mit Gl. 9-12 das folgende Gleichungssystem

11 12 11 12 1 10

21 22 21 22 2 20

cP c

⎧ ⎫⎡ ⎤ ⎡ ⎤ ⎡ ⎤ ⎡ ⎤ϕ ϕ ψ ψ ϕ⎪ ⎪− =⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎨ ⎬ϕ ϕ ψ ψ ϕ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎪ ⎪⎣ ⎦ ⎣ ⎦ ⎣ ⎦ ⎣ ⎦⎩ ⎭

oder symbolisch

( ) bcΨΦ =⋅− P

zu lösen. Gl. 9-12 geht bei fehlender Querbelastung ( 00j =ϕ ) in ein homogenes Gleichungs-

system über, das für

( ) 0Pdet =− ΨΦ Gl. 9-14

auch von Null verschiedene ck und damit auch von Null verschiedene Transversalverschie-bungen ermöglicht. Das ist das Kennzeichen für das Eintreten des Stabilitätsfalles. Mecha-nisch bedeutet dies, dass es neben der gestreckten Lage des Stabes eine infinitesimal benach-barte ausgelenkte Lage gibt, die ebenfalls eine Gleichgewichtslage ist. Die Eigenwertglei-chung Gl. 9-14 enthält als freien Parameter nur noch die Axialbelastung P. Der kleinste erre-

chenbare Wert dieser Axialbelastung ist ein Näherungswert kritP für die kritische Knicklast

des Systems. Für eine erste Abschätzung der kleinsten kritischen Last kritP setzen wir

(x)wc(x)wc(x)w 11

n

1kkk == ∑

=

ˆ

Bei der Auswahl von w1(x) ist darauf zu achten, dass die Ansatzfunktion möglichst wenig gekrümmt ist, denn mit

11

11krit11krit11 ψ

P0ψP ϕ=⇒=−ϕ Gl. 9-15

wird die kleinste kritische Last dann berechnet, wenn der Zähler in Gl. 9-15, also die Formän-derungsenergie 11ϕ des Stabes, möglichst klein gehalten wird. Das wird erreicht, wenn wir

F. U. Mathiak 9-9

mit w1(x) eine Stabauslenkung wählen, die zwischen 0 < x < l keine weitere Nullstelle be-sitzt, also etwa

)x(1x(x)w1ll

−=

Es handelt sich bei diesem Ansatz um eine zulässige Funktion, da die Verschiebungsrandbe-dingungen an den Stabenden erfüllt werden. Mit Gl. 9-15 errechnen wir

2yy3

yy

0x

2'1

0x

2''1yy

11

11krit

EI12

31

EI4

dxw

dxwEI

ψP

ll

ll

l

===ϕ

=

∫

∫

=

=

Einen wesentlich besseren Wert für die kleinste Knicklast erhalten wir, wenn wir anstatt einer zulässigen Funktion mit

⎥⎥⎦

⎤

⎢⎢⎣

⎡⎟⎠⎞

⎜⎝⎛+⎟

⎠⎞

⎜⎝⎛−=

32

1xx21x(x)wlll

Gl. 9-16

eine Vergleichsfunktion wählen, die, wie durch Ausrechnen leicht bestätigt wird, wegen

41)x(x12)x(w

l

l−=′′ auch die dynamischen Randbedingungen My(x = 0) = My(x = l ) = 0

erfüllt. Mit Gl. 9-16 errechnen wir

2yy

2yy3

yy

0x

2'1

0x

2''1yy

11

11krit

EI882,9

EI17168

35175EI24

dxw

dxwEI

ψP

lll

ll

l

====ϕ

=

∫

∫

=

=

und damit einen Wert, der der exakten Lösung

2yy

2yy2

krit

EI870,9

EIπP

ll==

schon recht nahe kommt. Der Nachteil des Ritzschen Verfahrens liegt auf der Hand. Die Auswahl zulässiger Funktio-nen oder das Auffinden von Vergleichsfunktionen erfordert vom Anwender ein hohes Maß an Erfahrung. Insbesondere bei mehrdimensionalen oder ungleichmäßig berandeten Gebieten stößt das Auffinden geeigneter Ansatzfunktionen oft auf unüberwindliche Schwierigkeiten. Wie wir im Folgenden sehen werden, entfällt dieser Nachteil bei der FE-Methode.

10-1

10 Finite Elemente bei eindimensio-nalen Randwertproblemen Wir haben im vorigen Kapitel die näherungsweise Lösung eines Randwertproblems mit Hilfe des Ritzschen Verfahrens kennen gelernt. Die Herleitung der Methode der finiten Elemente (FEM) kann nun auf einfache Weise erfolgen. Die grundlegende Idee besteht darin, anstelle eines Näherungsansatzes für das Gesamtgebiet Ansätze zu wählen, die sich lediglich auf Teil-bereiche der Struktur beziehen und auch nur dort von Null verschieden sind1. Wir unterstel-len, dass ein zu minimierendes Funktional Π vorliegt. Zur Erläuterung der weiteren Vorge-hensweise beziehen wir uns auf das Beispiel nach Abb. 10-1.

Abb. 10-1 Dehnstab mit linear veränderlichem Querschnitt A(x)

Zum Vergleich mit der noch bereitzustellenden FE-Lösung beschaffen wir uns zunächst die analytische Lösung des Problems. Normalkraft: )x(nFN 0 −+= l

Materialgesetz: [ ])x(nF)x(A

1)x(A

NuEE 0xxxx −+==′=ε=σ l

mit ( )ξ+=⎟⎟⎠

⎞⎜⎜⎝

⎛ −+= A~1Ax

AAA

1AA(x) rl

l

ll

l ( ) ll A/AAA~ r −= , lξ=x

1 Richard Courant, amerikan. Mathematiker dt. Herkunft, 1888-1972

10-2 Finite Elemente bei eindimensionalen Randwertproblemen

Abb. 10-2 Verschiebung u(x), exakte Lösung

Abb. 10-3 Spannung σxx, exakte Lösung

10-3

DGL: )ξA~1(EA)1(nF

)x(EA)x(nF

)x(u 00

+ξ−+

=−+

=′l

ll

Integration der Differentialgleichung und Anpassung an die Randbedingung u(0) = 0 liefert

1. Verschiebung: ( ) ( )⎭⎬⎫

⎩⎨⎧

⎥⎦

⎤⎢⎣

⎡ξ−+⎟

⎠⎞

⎜⎝⎛ +++= ξA~1ln

A~11nξA~1lnF

A~EA1u 2

0lll

2. Dehnung: )ξA~1(EA)1(nF 0

xx +ξ−+

=εl

l

3. Spannung: )ξA~1(A)1(nF

E 0xxxx +

ξ−+=ε=σ

l

l

Mit den Zahlenwerten aus Abb. 10-1 ergibt die numerische Auswertung:

( )2

xxxx2

xx cm/kN20)(max]cm/kN[910525

cm1843,0)(uumax]cm[9,01ln072,001852,0u

=σ=σ→ξ−ξ−

=σ

==→ξ−−ξ=

l

l

Für den Vergleich mit der noch zu beschaffenden FE-Lösung werten wir die Querschnittsflä-chen, die Verschiebungen sowie die Spannungen an diskreten Stellen aus.

l/x=ξ ]cm[)(A 2ξ u(ξ) [cm] σxx(ξ) [kN/cm2] 0,000 10,000 0,0000 2,5000 0,125 8,875 0,0109 2,7465 0,250 7,750 0,0230 3,0645 0,375 6,625 0,0366 3,4906 0,500 5,500 0,0523 4,0909 0,625 4,375 0,0711 5,0000 0,750 3,250 0,0948 6,5385 0,875 2,125 0,1277 9,7059 1,000 1,000 0,1843 20,0000

Tabelle 10-1 Querschnittswerte, Verschiebungen und Spannungen an diskreten Stellen

10.1 Vorgehensweise nach der FE-Methode

Das zuvor analytisch behandelte Beispiel soll nun mit der FE-Methode gelöst werden. Dazu gehen wir in Schritten vor. 1. Schritt Idealisierung des Tragwerkes


Diese Phase der Problemlösung gestaltet sich genauso wie bei der analytischen Vorgehens-weise. Wir stellen zunächst fest, dass es sich bei der vorliegenden Aufgabe um ein statisches Problem mit konstanten Lasten handelt, das nach den Grundlagen der linearen Elastizitätsthe-orie als Stab berechnet werden kann. Die einzige Schnittlast ist die Normalkraft N(x). Wir unterstellen bei der Modellbildung

Konstante Lasten Hookesches Materialgesetz Kleine Verformungen Stabtheorie

2. Schritt Formulierung der Variationsaufgabe Die Formänderungsenergie des Stabes ist

dx)x(u)x(EA21W

0x

2∫=

′=l

Die äußere Arbeit resultiert aus der konstanten Normalkraftschüttung n0 längs der Stabachse und der Einzelkraft F am Stabende

)(Fudx)x(unA0

0a ll

+= ∫

Ausgehend vom elastischen Potential

Extremum)(Fudx)x(undx)x(u)x(EA21AW)u(

0x0

0x

2a =−−′=−=Π ∫∫

==

lll

Gl. 10-1

erhalten wir aus dem Prinzip der virtuellen Verrückung die Extremwertaufgabe

0)(uFdx)x(undx)x(u)x(u)x(EA)u(0x

00x

=δ−δ−′δ′=Πδ ∫∫==

lll

Gl. 10-2

3. Schritt Diskretisierung1 des Tragwerkes Die vorgegebene Aufgabe wird in dem Sinne diskretisiert, dass der Stab in Teilgebiete, die finiten Elemente, zerlegt wird. Bei unserem Einführungsbeispiel (Fachwerk) und der vorlie-genden Stabaufgabe ist die Elementierung bereits vorgegeben. In beiden Fällen entspricht das finite Element einem Fachwerkstab oder Teilen davon. Ähnlich ist übrigens die Situation bei Rahmentragwerken. Dort stellen die Balken oder Balkenstücke die finiten Elemente dar.

1 zu lat. discernere ›absondern‹, ›unterscheiden‹

10-5

Abb. 10-4 Elementierung eines Dehnstabes, 2 Elemente gleicher Länge

Abb. 10-4 zeigt eine mögliche Elementierung unseres Dehnstabes mit zwei gleichlangen E-lementen und drei globalen Knoten. Die globalen Knotenkoordinaten sind in der Knotendatei abgelegt (Tabelle 10-2). In der Elementdatei (Tabelle 10-3) werden den Elementen die glo-balen Anfangs- und Endknoten zugeordnet.

Knotennummer x-Koordinate [cm] 1 0 2 50 3 100

Tabelle 10-2 Knotendatei

Elementnummer Anfangsknoten Endknoten

1 1 2 2 2 3

Tabelle 10-3 Elementdatei

4. Schritt Auswahl des Elementtyps, Ansatzfunktionen Wir betrachten in einem ersten Schritt das 2-Knotenelement entsprechend Abb. 10-5. Jeder Knoten besitzt nur einen kinematischen Freiheitsgrad, das ist die Verschiebung in Stablängs-richtung. Wir sprechen deshalb von einem Element mit zwei Freiheitsgraden.


Abb. 10-5 Dehnstab, 2-Knotenelement

Der Stabanfangsknoten wird mit 1 und der Stabendknoten mit 2 bezeichnet. Hinweis: Anfangs- und Endknoten des Stabes dürfen nicht mit der globalen Knotennumerie-rung nach Abb. 10-4 verwechselt werden. Damit besitzt das Element eine Orientierung. Für die Verschiebung u verwenden wir inner-halb des Elementes einen Verschiebungsansatz, der für alle Elemente identisch ist. Damit entfällt das umständliche Suchen nach geeigneten globalen Ansatzfunktionen, wie das beim Ritz-Verfahren erforderlich ist. Für Elemente, die sich am Rande des Lösungsgebietes befin-den, sind später die vom Problem vorgegebenen Randwerte einzuhalten. An den Element-grenzen müssen die Ansatzfunktionen gewissen Stetigkeitsanforderungen genügen, die vom mechanischen Problem abhängen. Beim Dehnstab müssen an den Elementübergängen die Verschiebungen stetig sein. In der FEM spricht man in diesem Fall von C0-Stetigkeit der Verschiebung u. Hinweis: C0-Stetigkeit reicht z.B. für die Durchbiegung beim Biegebalken nicht aus. Hier muss zusätzlich noch Stetigkeit in der 1. Ableitung von w(x) gefordert werden, die als C1-Stetigkeit bezeichnet wird (Abb. 10-6).

Abb. 10-6 Verletzung der C1-Stetigkeit bei einem Balkenelement

Zur Beschreibung des Verschiebungszustandes innerhalb des Elementes verwenden wir die lokale Koordinate )e()e(x l=ξ ( )10 ≤ξ≤ . Wir benötigen zur Lösung des Problems die glo-

balen Verschiebungen. Dazu stellen wir zunächst den Verschiebungszustand innerhalb eines Elementes auf und sorgen durch entsprechende Wahl der Verschiebungsfunktionen für Ste-tigkeit an den Elementgrenzen.

10-7

Zur Darstellung der Verschiebung innerhalb des Elementes eignen sich speziell Polynome, da diese leicht aufgebaut und auch leicht zu differenzieren sind. Ein Polynom vom Grade n be-sitzt genau n + 1 freie Parameter. Es hat die Darstellung

∑=

ξ=ξ+ξ+ξ+=ξn

0i

ii

nn

2210n aaaaa)(P K Gl. 10-3

Da von der Zustandsgröße u(x) unseres Stabproblems lediglich C0-Stetigkeit gefordert wird, genügt es bei einem 2-Knotenelement einen linearen Ansatz für die Verschiebungen zu wäh-len, also

ξ+=ξ 10 aa)(u ( )10 ≤ξ≤ Gl. 10-4

Dabei entspricht der Wert a0 im Ansatz Gl. 10-4 einer Starrkörperverschiebung des Elemen-tes. Wir drücken nun die Verschiebungen innerhalb des Elementes durch die Stabendver-schiebungen u1 und u2 aus. Das erreichen wir, indem wir die Verschiebungsfunktion Gl. 10-4 an den Rändern ξ = 0 und ξ = 1 auswerten.

121102

1001

uuaaau)1(uuaau)0(u

−=+===ξ====ξ

Einsetzen der Konstanten a0 und a1 in Gl. 10-4 liefert: ξ−+=ξ )uu(u)(u 121 . Wir sortieren

noch etwas um und erhalten

221121 u)(Nu)(Nuu)1()(u ξ+ξ=ξ+ξ−=ξ Gl. 10-5

Damit haben wir die Verschiebungen innerhalb des Elementes zunächst durch die Stabend-verschiebungen ausgedrückt. Die Funktionen vor den Stabendverschiebungen in Gl. 10-5

1 1

2 2

N ( ) 1 L ( )N ( ) L ( )

ξ = −ξ = ξξ = ξ = ξ

Gl. 10-6

sind die Lagrangeschen1 Interpolationspolynome Li, die in der FEM auch Formfunktio-nen2 genannt werden (Abb. 10-7). Die Lagrangeschen Interpolationspolynome nehmen an den Knoten i gerade den Wert 1 und an den übrigen Knoten den Wert 0 an

1 Joseph Louis de Lagrange, frz. Mathematiker und Physiker, 1736-1813 2 In der angelsächsischen Literatur werden diese Funktionen als shape functions bezeichnet


ii i

k

1 fürN ( ) L ( )

0 für (k i)⎧ ξ⎪ξ = ξ = ⎨ ξ = ξ ≠⎪⎩

Abb. 10-7 Formfunktionen, linearer Ansatz

Weiterhin gilt: n

ii 1

L 1=

=∑ .

Zur Herleitung der Elementsteifigkeitsmatrix gehen wir aus vom Funktional Gl. 10-1, das hier noch einmal notiert werden soll

Extremum)(Fudx)x(undx)x(u)x(EA21AW)u(

0x0

0x

2a =−−′=−=Π ∫∫

==

lll

Gl. 10-7

Da sämtliche Energieausdrücke in Gl. 10-7 additiv sind, dürfen wir diese elementweise berech-nen und letztlich durch Summation zum Gesamtpotential zusammenführen. Für die Formän-derungsenergie erhalten wir

)e(n

1e 0x

2)e()e(n

1e

)e(

0x

2 dxuA2

EWdx)x(u)x(EA21W

)e(

)e(∑ ∫∑∫= ===

′==′=ll

Gl. 10-8

und für die Arbeit der äußeren Kräfte folgt

)(Fudx)x(un)(Fudx)x(unAn

1e 0xee0

0x0a

e

e

llll

+=+= ∑ ∫∫= ==

Gl. 10-9

10-9

Die Variation des Funktionals Gl. 10-7 ist

0AW a =δ−δ=Πδ Gl. 10-10

mit

∑ ∫∑= ==

′δ′=δ=δn

1e

)e(

0x

)e(n

1e

)e( dxAuuEWW)e(

)e(

l

Gl. 10-11

und

)(uFdxunAn

1e 0x

)e(0a

)e(

)e(

ll

δ+δ=δ ∑ ∫= =

Gl. 10-12

Wir beschaffen uns jetzt eine Näherungswert für das Potential Π, indem wir für die Verschie-bungen auf Elementebene den linearen Verschiebungsansatz

( ) 221121 u)(Nu)(Nuu1)(u ξ+ξ=ξ+ξ−=ξ Gl. 10-13

wählen. Damit geht das Potential Π(u) über in die Näherung )u(Π . Gl. 10-13 können wir

auch kompakter in Matrizenschreibweise notieren

(e) (e)1 1(e) (e) (e) (e) (e)

1 2(e) (e)2 2

u ux xu(x ) 1 , N (x ), N (x )u u

N u⎡ ⎤ ⎡ ⎤ ⎡ ⎤

⎡ ⎤= − = =⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎣ ⎦⎣ ⎦ ⎣ ⎦⎣ ⎦l l

Gl. 10-14

Die virtuelle Verrückung ist

(e)(e) δu uN=δ Gl. 10-15

Wir benötigen noch die Ableitung der Verschiebungsfunktion u nach der Variablen x, die der Dehnung xxε entspricht. Es gilt

(e)(e)

2

1)e()e(

2

1)e(

2)e(

1)e((e)xx u

u1,1uu

dxdN,

dxdN)x(uε uB=⎥

⎦

⎤⎢⎣

⎡⎥⎦

⎤⎢⎣

⎡−=⎥

⎦

⎤⎢⎣

⎡⎥⎦

⎤⎢⎣

⎡=′=

ll Gl. 10-16

Die Matrix


[ ] ⎥⎦

⎤⎢⎣

⎡−=⎥

⎦

⎤⎢⎣

⎡== )e()e()e(

2)e(

121

(e) 1,1dxdN

,dxdN

B,Bll

B Gl. 10-17

enthält die Ableitungen der Ansatzfunktionen. Die Dehnungen auf Elementebene lassen sich dann, ausgedrückt durch die Knotenverschiebungen, wie folgt schreiben

(e)(e))e((e)xx )x(uε uB=′= Gl. 10-18

Beachten wir mit Gl. 10-16

(e)(e)u uB δ=′δ

dann liefert das Einsetzen in Gl. 10-11

(e)(e)T(e)(e)

0

)e((e))e(T(e)T(e)

0

)e((e)(e))e((e)(e))e(

)e(

)e(

dxAE

dxAE)u(W

ukuuBBu

uBuB

(e)k

δ=δ=

δ=δ

=

∫

∫

444 3444 21

l

l

Gl. 10-19

Bezeichnen A1 und A2 die Querschnittsflächen am Stabanfang (Knoten 1) bzw. Stabende

(Knoten 2) dann gilt: )e(

)e(

121)e( x)AA(A)x(A

l−+= . Ausmultiplizieren der Matrizen und

anschließende Integration liefert die Elementsteifigkeitsmatrix

( )⎥⎦

⎤⎢⎣

⎡−

−=⎥

⎦

⎤⎢⎣

⎡−

−+=

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

−+−−−

−−−−+== ∫∫

1111AE

1111AAE

21

dxx)AA(Ax)AA(A

x)AA(Ax)AA(AEdx)x(AE

)e(

)e(m

)e(

)e(21

)e(

0

)e(

)e(

)e(

121)e(

)e(

121

)e(

)e(

121)e(

)e(

121

2)e(

)e(

0

)e()e((e))e(T(e)(e)

)e()e(

ll

ll

lll

ll

BBk

Gl. 10-20

die offensichtlich nur von der Länge )e(l und der Dehnsteifigkeit )e(

m)e( AE des betrachteten

Elementes abhängt. In Gl. 10-20 ist

10-11

( )21)e(

m AA21A += Gl. 10-21

der Mittelwert der Querschnittsflächen am Stabanfang und Stabende. Ist die Querschnittsflä-che A(e) = A1 =A2 innerhalb des Elementes konstant, dann geht Gl. 10-20 über in

⎥⎦

⎤⎢⎣

⎡−

−=

1111AE

)e(

)e()e((e)

lk Gl. 10-22

Die virtuelle Arbeit der äußeren Kräfte setzt sich aus zwei Anteilen zusammen:

1. Arbeit der Linienlast n0 an der virtuellen Verschiebung δu 2. Arbeit der Einzellast F an der virtuelle Verschiebung δu( l )

Wir betrachten zuerst die Linienlast n0. Die Einzellast F am Stabende wird später berücksich-tigt. Wir erhalten durch Summation über sämtliche Elemente

∑∑ ∫==

δ=δ=δn

1e

)e(a

n

1e 0

)e()e(0a Adx)x(unA

)e(l

Gl. 10-23

mit dem auf das Element bezogenen Anteil

(e)(e)T

0x

)e(0

(e)T(e)T

0

)e(0

(e)(e)

0

)e()e(0

)e(a

(e)

)e(

e

)e()e(

dxndxndx)x(unA puNuuN

p

δ=δ=δ=δ=δ

=

=∫∫∫

44 344 21

lll

Gl. 10-24

Die Integration liefert für n = n0 = konst. den Elementlastvektor

(e )

e

(e)(e)T (e) 0

0x 0

1nn dx12

(e)p N=

⎡ ⎤= = ⎢ ⎥

⎣ ⎦∫l l Gl. 10-25

5. Schritt Aufbau des globalen Gleichungssystems Die auf Elementebene berechneten Matrizen werden jetzt zum globalen Gleichungssystem zusammengebaut. Der Zusammenbau muss so erfolgen, dass die Kompatibilität in den Ver-formungen über die Elementgrenzen hinweg gewährleistet ist. Dazu benötigen wir den Zu-sammenhang zwischen den lokalen Elementknotenverschiebungen (u1, u2) und den globalen Systemknotenverschiebungen


⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

3

2

1

vvv

v Gl. 10-26

Hierzu führen wir auf Elementebene die Zuordnungsmatrix A(e) ein, die sich aus der Element-datei ermitteln lässt. Es gilt dann

vAu (e)(e) = Gl. 10-27

Die Matrix A(e) ist eine Boolesche Matrix, die nur die Informationen 0 oder 1 enthält (s.h. Beispiel Fachwerk). Für das Element 2 erhalten wir zum Beispiel

vAu )2(

3

2

1)2(

2

1)2(

vvv

100010

uu

=⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

⎥⎦

⎤⎢⎣

⎡=⎥

⎦

⎤⎢⎣

⎡= Gl. 10-28

Das Potential des Elementes ist dann

(e)(e)T(e)(e)(e)T

0

)e(0

)e(

0

2)e()e(

21dxundxuAE

21

)e()e(

puuku −=−′=Π ∫∫ll

Gl. 10-29

Am Gesamtpotential fehlen jetzt noch die Arbeiten der äußeren Kräfte, die aus der Kraft F am rechten Rand und der zunächst unbekannten Auflagerkraft R am linken Rand bestehen. Um aus der Reaktionskraft R eine äußere Kraft zu machen, muss das kinematische Modell durch das entsprechende statische Modell ersetzt werden (Abb. 10-8).

Abb. 10-8 Statisches Modell zur Darstellung der Reaktionskraft R

Die Reaktionskraft R wird am Systemknoten "1" und die Kraft F am Systemknoten "3" in die Konstruktion eingeleitet. Führen wir mit

10-13

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

F0R

FFF

3

2

1

F Gl. 10-30

den Knotenlastvektor ein, dann ist die Arbeit der eingeprägten Kräfte

vFTF,aA =

Die Gesamtenergie des Stabes ist dann endgültig

[ ]

[ ]

Extremum21

221

221)u(ˆ

TT

Tn

1e

(e)(e)T(e)(e)(e)TT

n

1e

T(e)(e)T(e)(e)(e)Tn

1e

)e(

=−=

−−=

−−=Π=Π

∑

∑∑

=

==

PvvKv

FvpAvAkAv

vFpuuku

Gl. 10-31

mit der globalen Systemsteifigkeitsmatrix

∑=

=n

1e

(e)(e)(e)T AkAK Gl. 10-32

und dem globalen Knotenlastvektor

FpAP += ∑=

n

1e

(e)(e)T Gl. 10-33

Das Prinzip der virtuellen Verrückung fordert

[ ] 021

21ˆ TTTT =−δ=δ−δ+δ=Πδ PvKvPvvKvvKv Gl. 10-34

In Gl. 10-34 wurde die Symmetrie der Systemsteifigkeitsmatrix berücksichtigt1. Wegen der Beliebigkeit von vδ muss dann

1 Dann gilt: vKvvKv TT δ=δ


PvK = Gl. 10-35

erfüllt sein. 6. Schritt Aufbau des modifizierten Gleichungssystems, Einbau der Randbedingun-gen Mit Hilfe der Zuordnungsmatrizen A(e) lässt sich nun das Gleichungssystem für das noch nicht gefesselte System aufbauen. Wir beginnen mit der Systemmatrix K

)2()2()T2()1()1()T1()2()1(2

1e

(e)(e)(e)T AkAAkAKKAkAK +=+== ∑=

Gl. 10-36

Zum Aufbau der Elementsteifigkeitsmatrix k(e) benötigen wir gemäß Gl. 10-21 neben den Elastizitätsmoduli die Mittelwerte der Elementquerschnittsflächen )e(

mA

Element Elastizitätsmodul [kN/cm2] (e)

mA [cm2] Elementlänge [cm]

1 3000,0 0,5 (10 + 5,5) = 7,75 50,0 2 3000,0 0,5 (5.5 + 1,0) = 3,25 50,0

Tabelle 10-4 Elementgrößen, zwei Elemente gleicher Länge

Damit erhalten wir nach Gl. 10-20 die Elementsteifigkeitsmatrizen

⎥⎦

⎤⎢⎣

⎡−

−=⎥

⎦

⎤⎢⎣

⎡−

−=⎥

⎦

⎤⎢⎣

⎡−

−=⎥

⎦

⎤⎢⎣

⎡−

−=

1111

1951111AE

1111

4651111AE

)2(

)2(m

)2()2(

)1(

)1(m

)1()1(

llkk

sowie unter Beachtung von Gl. 10-32 die globalen Elementmatrizen

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡−

−=⎥

⎦

⎤⎢⎣

⎡⎥⎦

⎤⎢⎣

⎡−

−

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡==

00004654650465465

010001

465465465465

001001

)1()1()T1()1( AkAK

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−−=⎥

⎦

⎤⎢⎣

⎡⎥⎦

⎤⎢⎣

⎡−

−

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡==

19519501951950

000

100010

195195195195

100100

)2()2()T2()2( AkAK

Der Zusammenbau liefert die globale Systemsteifigkeitsmatrix

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−−−

−=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−−+−

−=+=

1951950195660465

0465465

1951950195195465465

0465465)2()1( KKK

10-15

Nun wird die rechte Seite aufgebaut. Es gilt

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡+

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡+

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=++=+= ∑

= 25,2150,225,1

2000

110

25,1011

25,1)2()T2()1()T1(n

1e

(e)(e)T FpApAFpAP

Die Systemgleichung Gl. 10-35 lautet dann

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−−−

−

25,2150,225,1

vvv

1951950195660465

0465465

3

2

1

Wie man leicht zeigen kann, ist die Systemmatrix K singulär (det K = 0). Das System ist of-fensichtlich kinematisch, denn wir haben noch nicht die Lagerungsbedingung des Stabes am linken Rand (x = 0) berücksichtigt. Wegen v1 = 0 und damit auch 0v1 =δ können die 1. Spal-

te und die 1. Zeile des Gleichungssystems gestrichen werden. Es verbleibt dann die reduzierte Systemgleichung

⎥⎦

⎤⎢⎣

⎡=⎥

⎦

⎤⎢⎣

⎡⎥⎦

⎤⎢⎣

⎡−

−25,2150,2

vv

195195195660

3

2 → ⎥⎦

⎤⎢⎣

⎡=⎥

⎦

⎤⎢⎣

⎡1600,00511,0

vv

3

2 [cm]

Setzen wir die nun bekannten Knotenverschiebungen v2 und v3 in den Knotenverschiebungs-vektor ein, dann erhalten wir die Systemgleichung

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡ +=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−−−

−

25,2150,2

R25,1

1600,00511,00

1951950195660465

0465465 Gl. 10-37

in der in der 1. Zeile nur noch die unbekannten Reaktionskraft R steht. Ausmultiplizieren lie-fert

R25,1v465 2 +=−

und damit kN0,2525,10511,0465R −=−⋅−=

7. Schritt Rückrechnung Mit den Knotenverschiebungen liegen dann auch die verbleibenden Zustandsgrößen fest. Da wir einen linearen Verschiebungsansatz gewählt haben, sind die Verzerrungen (und damit auch die Spannungen) innerhalb des Elementes konstant1. Mit Gl. 10-16 gilt

1 In der angelsächsischen Literatur wird deshalb ein solches Element als constant strain element bezeichnet.


vABuBε (e)(e)(e)(e)(e) == Gl. 10-38

Element 1:

3e022,1500511,0

1600,00511,00

010001

501,

501

vvv

0100011,1

3

2

1

)1()1()1()1()1(

−==

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

⎥⎦

⎤⎢⎣

⎡⎥⎦⎤

⎢⎣⎡−=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

⎥⎦

⎤⎢⎣

⎡⎥⎦⎤

⎢⎣⎡−==

llvABε

Element 2:

( ) 3e178.20511,016,0501

1600.00511.00

100010

501,

501

vvv

1000101,1

3

2

1

)2()2()2()2()2(

−=−=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

⎥⎦

⎤⎢⎣

⎡⎥⎦⎤

⎢⎣⎡−=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

⎥⎦

⎤⎢⎣

⎡⎥⎦⎤

⎢⎣⎡−==

llvABε

Für die Spannungen gilt

(e))e((e) E εσ = Gl. 10-39

Element 1: 2)1()1()1( cm/kN07,33E022,13E0,3E =−⋅+== εσ

Element 2: 2)2()2()2( cm/kN53,63E178,23E0,3E =−⋅+== εσ

Die grafischen Darstellungen der Verschiebungen und Spannungen in Abb. 10-9 und Abb. 10-10 zeigen deutlich ein starkes Anwachsen der Zustandsgrößen in der Nähe des rechten Randes. Das trifft besonders für die Spannungen zu. Die Elementierung des Stabes mit nur zwei Elementen unter Verwendung des Zweiknotenelementes mit linearem Verschiebungsan-satz ist offensichtlich nicht in der Lage, insbesondere die Spannungen befriedigend wieder-zugeben. Die mit der FEM ermittelten Spannungen entsprechen jedoch exakt den theoreti-schen Werten in Elementmitte. Diese Mittelung im energetischen Sinne ist charakteristisch für die FE-Methode. Die größte Spannung tritt bei l=x auf. Der relative Fehler beträgt im-merhin

%4,6720

53,620

an

FEan =−

=σ

σ−σ=σ∆

ein für die sinnvolle Auslegung des Systems zu hoher Wert.

10-17

Abb. 10-9 Verschiebung u, Ergebnis für zwei Elemente gleicher Länge

Abb. 10-10 Spannung σxx, Ergebnis für zwei Elemente gleicher Länge


Mit dem vorgestellten 2-Knoten-Element lassen sich die Ergebnisse durch folgende Modifika-tionen wesentlich verbessern: 1. Erhöhung der Anzahl der Elemente, wobei alle Elemente dieselbe Elementlänge besitzen

h-Adaptivität (h = Elementlänge) 2. Feinere Elementierung im Bereich starker Änderung der Zustandsgrößen

r-Adaptivität (r = Knotenabstand) Hinweis: Die Konvergenz einer Näherungslösung gegen die im Sinne der Theorie exakte Lö-sung kann dadurch beobachtet werden, dass in eine grobe Vernetzung ein feineres Netz gelegt wird, wobei sich die Lage der Knoten und der Elementgrenzen der gröberen Vernetzung nicht ändert. Nähert sich dann die Lösung einem Grenzwert, dann kann auf Konvergenz der Nähe-rungslösung gegen die exakte Lösung geschlossen werden. Wir erhöhen in einem ersten Schritt die Anzahl der Elemente von zwei auf vier. Den vier E-lementen gleicher Länge sind 5 Systemknoten mit den entsprechenden Freiheitsgraden zuge-ordnet.

Abb. 10-11 Elementierung eines Dehnstabes, 4 Elemente gleicher Länge

Knotennummer x-Koordinate [cm]

1 0 2 25 3 50 4 75 5 100

Tabelle 10-5 Knotendatei, 4 Elemente gleicher Länge

10-19


1 1 2 2 2 3 3 3 4 4 4 5

Tabelle 10-6 Elementdatei, 4 Elemente gleicher Länge

Zur Berechnung der Elementsteifigkeitsmatrizen benötigen wir wieder die gemittelten Quer-schnittsflächen )e(

mA

Element Elastizitätsmodul [kN/cm2] (e)

mA [cm2] Elementlänge [cm]

1 3000,0 0.5 (10,0 + 7,75) = 8,875 25 2 3000,0 0.5 (7,75 + 5,50) = 6,625 25 3 3000,0 0.5 (5,50 + 3,25) = 4,375 25 4 3000,0 0.5 (3,25 + 1.00) = 2,125 25

Tabelle 10-7 Elementgrößen, 4 Elemente gleicher Länge

Elementlastvektoren

⎥⎦

⎤⎢⎣

⎡=⎥

⎦

⎤⎢⎣

⎡⋅=⎥

⎦

⎤⎢⎣

⎡=

11

625,011

225,05

11

2n )e(

0(e) lp

Elementsteifigkeitsmatrizen

⎥⎦

⎤⎢⎣

⎡−

−=⎥

⎦

⎤⎢⎣

⎡−

−=⎥

⎦

⎤⎢⎣

⎡−

−=⎥

⎦

⎤⎢⎣

⎡−

−=

⎥⎦

⎤⎢⎣

⎡−

−=⎥

⎦

⎤⎢⎣

⎡−

−=⎥

⎦

⎤⎢⎣

⎡−

−=⎥

⎦

⎤⎢⎣

⎡−

−=

255255255255

1111AE

525525525525

1111AE

795795795795

1111AE

0651065106510651

1111AE

)4(

)4()4()4(

)3(

)3()3()3(

)2(

)2()2()2(

)1(

)1()1()1(

ll

ll

kk

kk

Vektor der rechten Seite

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

=

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

+

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

+

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

+

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

+

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

=+= ∑=

625,20250,1250,1250,1625,0

200000

11000

625,0

01100

625,0

00110

625,0

00011

625,04

1e

(e)(e)T FpAP


Systemgleichung

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡ +

=

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

−−+−

−+−−+−

−

625,20250,1250,1250,1

R625,0

vvvv0

2552550002552555255250005255257957950007957951065106500010651065

5

4

3

2

Knotenverschiebungen

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

=

1745,00936,00520,00229,00

v

Reaktionskraft

kN00,25625,002289,01065R =−⋅−=

Verzerrungen

(1) (1) (1) (1) (1)

00,0229

1 0 0 0 01 1 0,0229, 9,160e 40,05200 1 0 0 025 25 25

0,09360,1745

ε B u B A v

⎡ ⎤⎢ ⎥⎢ ⎥⎡ ⎤⎡ ⎤ ⎢ ⎥= = = − = = −⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦⎢ ⎥⎢ ⎥⎣ ⎦

(2) (2) (2) (2) (2)

00,0229

0 1 0 0 01 1 0,052 0,0229, 1,164 30,05200 0 1 0 025 25 25

0,09360,1745

ε B u B A v

⎡ ⎤⎢ ⎥⎢ ⎥⎡ ⎤⎡ ⎤ −⎢ ⎥= = = − = = −⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦⎢ ⎥⎢ ⎥⎣ ⎦

(3) (3) (3) (3) (3)

00,0229

0 0 1 0 01 1 0,0936 0,052, 1,664e 30,05200 0 0 1 025 25 25

0,09360,1745

ε B u B A v

⎡ ⎤⎢ ⎥⎢ ⎥⎡ ⎤⎡ ⎤ −⎢ ⎥= = = − = = −⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦⎢ ⎥⎢ ⎥⎣ ⎦

(4) (4) (4) (4) (4)

00,0229

0 0 0 1 01 1 0,1745 0,0936, 3, 236e 30,05200 0 0 0 125 25 25

0,09360,1745

ε B u B A v

⎡ ⎤⎢ ⎥⎢ ⎥⎡ ⎤⎡ ⎤ −⎢ ⎥= = = − = = −⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦⎢ ⎥⎢ ⎥⎣ ⎦

10-21

Abb. 10-12 Verschiebung u, Ergebnis für vier Elemente gleicher Länge

Abb. 10-13 Spannung σxx, Ergebnis für vier Elemente gleicher Länge


Spannungen 2)1()1()1( cm/kN75,24E160,93E0,3E =−⋅+== εσ

2)2()2()2( cm/kN49,33E164,13E0,3E =−⋅+== εσ 2)3()3()3( cm/kN99,43E664,13E0,3E =−⋅+== εσ 2)4()4()4( cm/kN71,93E236,33E0,3E =−⋅+== εσ

10.2 Stabelement mit quadratischem Verschiebungsansatz Wird der Grad des Polynoms für die Verschiebungsfunktion )(u ξ erhöht, so sind weitere

Zwischenknoten erforderlich. Die automatische Anpassung der Polynomordnung in der An-satzfunktion wird als p-Adaptivität bezeichnet. Wählen wir z.B. ein Polynom 2. Ordnung, also

2210 aaa)(u ξ+ξ+=ξ Gl. 10-40

dann können wir den zusätzlichen Freiwert (hier a2) einem weiteren Knoten zuordnen, z.B. dem Mittelpunkt des Elementes in Abb. 10-14.

Abb. 10-14 Dehnstab, 3-Knoten-Element

Die drei Konstanten a0, a1 und a2 werden aus dem linearen Gleichungssystem

2103

2102

01

aaau)1(u4a2aau)21(u

au)0(u

++==++==

==

zu

3212

3211

10

u2u4u2auu4u3a

ua

+−=−+−=

= Gl. 10-41

ermittelt. Einsetzen dieser Konstanten in Gl. 10-40 liefert

10-23

( ) ( ) ( ) 32

22

12 u2u4u231)(u ξ+ξ−+ξ−ξ+ξ+ξ−=ξ Gl. 10-42

Mit den Formfunktionen

( )( )12)(N

14)(N)1)(12()(N

3

2

1

−ξξ=ξξ−ξ=ξ

−ξ−ξ=ξ Gl. 10-43

lauten dann die Elementverschiebungen

332211 u)(Nu)(Nu)(N)(u ξ+ξ+ξ=ξ Gl. 10-44

Abb. 10-15 Formfunktionen, quadratischer Ansatz

In Abb. 10-15 sind die Formfunktionen für das 3-Knoten-Element dargestellt. Hinweis: Wir hätten selbstverständlich den Knoten mit 10 2 <ξ< an jede beliebige Stelle

zwischen beide Endknoten legen können. Die Formfunktionen errechnen sich für diesen all-gemeinen Fall dann zu

( ) ( )1

)(N;)1(

1)(N;))(1(

)(N2

23

222

2

21 −ξ

ξ−ξξ=ξ

−ξξξ−ξ

=ξξ

ξ−ξ−ξ=ξ Gl. 10-45

Für 2/12 =ξ gehen die obigen Funktionen in die einfache Form Gl. 10-43 über, was die Wahl

des Knotens in der Elementmitte nachträglich als sinnvoll erscheinen lässt.


Die Aufstellung des linearen Gleichungssystems Gl. 10-41 zur Berechnung der Polynomkon-stanten kann vermieden werden, wenn man bei einem Polynom p-ten Grades das Bildungsge-setz der Lagrangeschen Interpolationspolynome

)()(

)()(

)()(

)()(

)()(

)()(

)()()(L

1pk

1p

pk

p

1kk

1k

1kk

1k

2k

2

1k

11p

)km(1m mk

mk

+

+

+

+

−

−+

≠= ξ−ξ

ξ−ξξ−ξξ−ξ

ξ−ξξ−ξ

ξ−ξξ−ξ

ξ−ξξ−ξ

ξ−ξξ−ξ

=ξ−ξξ−ξ

=ξ ∏ LL

beachtet. In unserem Fall des quadratischen Ansatzes Gl. 10-40 ist p = 2 und wir erhalten:

321

1 2 1 3

( )( ) ( 1 2) ( 1)L ( ) (2 1)( 1)( ) ( ) (0 1 2) (0 1)

ξ − ξξ −ξ ξ − ξ −ξ = = = ξ − ξ −

ξ − ξ ξ −ξ − −

312

2 1 2 3

( )( ) ( 0) ( 1)L ( ) 4 (1 )( ) ( ) (1 2 0) (1 2 1)

ξ − ξξ − ξ ξ − ξ −ξ = = = ξ − ξ

ξ − ξ ξ −ξ − −

1 23

3 1 3 2

( ) ( ) ( 0) ( 1 2)L ( ) (2 1)( ) ( ) (1 0) (1 1 2)ξ − ξ ξ −ξ ξ − ξ −

ξ = = = ξ ξ −ξ − ξ ξ −ξ − −

Die Lagrangeschen Interpolationspolynome sind identisch mit den Formfunktionen, die wir im Vektor der Formfunktionen N(e) zusammenfassen.

[ ] (e)(e)

3

2

1)e(

3)e(

2)e(

1)e(

uuu

)x(N),x(N),x(N)x(u uN=⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡= Gl. 10-46

Mit diesem Ansatz ist C0-Stetigkeit über die Elementgrenzen hinweg gesichert.

(e)(e)u uN δ=δ Gl. 10-47

Wir benötigen wieder die Ableitung der Verschiebungsfunktion u . Es gilt

[ ] (e)(e)

3

2

1

)e(

3

2

1

)e(3

)e(2

)e(1)e(

uuu

14,84,341

uuu

dxdN

,dxdN

,dxdN

)x(u uB=⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡−ξξ−−ξ=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

⎥⎦⎤

⎢⎣⎡=′

l Gl. 10-48

Die Querschnittsfläche des Stabes ist linear veränderlich. Bezeichnen A1 und A3 die Quer-schnittsflächen an den entsprechenden Elementknoten, dann gilt:

( ) 31)e(

)e(

131)e( AA)1(xAAA)x(A ξ+ξ−=−+=

l

und für die Elementsteifigkeitsmatrix erhalten wir ( 13 A/A=α )

10-25

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

α+α+−α+α+−α+α+−α+α+−α+

== ∫113)31(41

)31(4)1(16)3(41)3(4311

6AEdxAE )e(

1)e(

0

)e((e))e((e)T(e))e(

l

l

BBk Gl. 10-49

Unter der Voraussetzung elementweise konstanter Querschnittsfläche A1 = A3 = )e(0A = konst.

( 1=α ) folgt aus Gl. 10-49

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−−−

−=

7818168187

3AE

)e(

)e(0

)e()e(

lk Gl. 10-50

Entsprechend erhalten wir mit n0 = konst. den Elementlastvektor

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=ξ

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−ξξξ−ξξ+ξ−

== ∫∫=ξ 1

41

6n

d)12()1(4

231ndxn

)e(0

1

0

2

)e(0

)e(

00

(e)T(e)

)e(

ll

l

Np Gl. 10-51

Abb. 10-16 Dehnstab, 2 Elemente gleicher Länge (quadratischer Verschiebungsansatz)

Zur Berechnung der Elementsteifigkeitsmatrizen nach Gl. 10-49 benötigen wir noch die Quer-schnittsflächen an den Elementrändern.

Element E(e) [kN/cm2] (e)1A [cm2] (e)

3A [cm2] α = A3/A1 (e)l [cm]

1 3000 10 5.50 0.55 50 2 3000 5.50 1.00 0.1818 50

Tabelle 10-8 Elementgrößen


⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−−−

−=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−−−

−=

275340653401040700

65700635

9051060155106024801420

15514201265)2()1( kk

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡⋅

=⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

417,0667,1417,0

141

65,05

141

6n )e(

0(e) lp

Systemgleichung des gefesselten Systems

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡ +

=

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

−−−

−+−−−

−

417,20667,1834,0667,1

R417,0

vvvv0

2753406500340104070000657006359051060155001060248014200015514201265

5

4

3

2

Lösung 1.) Knotenverschiebungen

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

=

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

=

18022,009571,005228,002302,0

0

vvvvv

5

4

3

2

1

v [cm]

2.) Reaktionskraft

kN25417,005228,015502302,01420R −=−⋅+⋅−= 3.) Elementverschiebungen Element 1:

[ ]

]cm[10248,110979,318022,009571,005228,002302,0

0

001000001000001

)12(),1(4,231)(u

222

2)1()1()1()1(

ξ⋅+ξ⋅=

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡−ξξξ−ξξ+ξ−===ξ

−−

vANuN

10-27

Element 2:

[ ]

]cm[10216,810577,410228,518022,009571,005228,002302,0

0

100000100000100

)12(),1(4,231)(u

2222

2)2()2()2()2(

ξ⋅+ξ⋅+⋅=

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡−ξξξ−ξξ+ξ−===ξ

−−−

vANuN

4.) Spannungen

vABuB )e()e()e()e()e()e()e()e()e( EEE ==ε=σ Element 1:

[ ]

]cm/kN[499,1388,2

18022,009571,005228,002302,0

0

001000001000001

14,84,3450

3000E

2

)1()1()1()1(

ξ+=

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡−ξξ−−ξ==σ vAB

Element 2:

[ ]

]cm/kN[859,9746,218022,009571,005228,002302,0

0

100000100000100

14,84,3450

3000E

2

)2()2()2()2(

ξ+=

⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢

⎣

⎡

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡−ξξ−−ξ==σ vAB


Abb. 10-17 Verschiebung u, Ergebnis für zwei Elemente gleicher Länge (quadratischer Ansatz)

Abb. 10-18 Spannung σxx, Ergebnis für zwei Elemente gleicher Länge (quadratischer Ansatz)

10-29

Abb. 10-19 Verschiebung am Stabende in Abhängigkeit von der Anzahl der Elemente

Abb. 10-20 Spannung am Stabende in Abhängigkeit von der Anzahl der Elemente


Ein Blick auf das Verschiebungsfeld in Abb. 10-17 zeigt, dass, obwohl nur zwei Elemente verwendet wurden, bereits eine sehr gute Übereinstimmung mit den analytischen Werten be-steht. Für die Spannungen gilt das nur im ersten Element (Abb. 10-18). Dort ist die Span-nungsänderung jedoch nicht so stark, wie in der zweiten Hälfte. Am Elementübergang tritt ein Spannungssprung auf. Mit feiner werdender Elementierung lassen sich die Ergebnisse für die Verschiebungen und Spannungen noch erheblich verbessern. Abb. 10-19 und Abb. 10-20 zei-gen die Entwicklung der Zustandsgrößen bei Verwendung gleicher Elementlängen und zu-nehmender Anzahl der Elemente. Das Element mit quadratischem Verschiebungsansatz zeigt für beide Zustandsgrößen die besseren Ergebnisse. Allerdings ist beim Element mit quadrati-schem Verschiebungsansatz infolge des zusätzlichen Mittenknotens, und damit einer zusätzli-chen Unbekannten je Element, die Rechenzeit größer als beim linearen Element. Den analytischen Lösungen für die Verschiebungen und Spannungen ist zu entnehmen, dass die Änderungen der Zustandsgrößen (ihre Gradienten) mit Annäherung an den rechten Rand zunehmen. Es liegt daher nahe, auch die Elementierung zum rechten Rand hin zu verdichten. Wir wählen abschließend ein Elementnetz entsprechend Abb. 10-21. Die rechte Stabhälfte wurde nochmals in 2 Elemente gleicher Länge (25 cm) unterteilt. Die analytische Lösung des linken Bereichs zeigt eine nahezu lineare Veränderlichkeit der Zustandsgrößen, hier kann also recht grob elementiert werden.

Abb. 10-21 Elementierung eines Dehnstabes, 3 Elemente ungleicher Länge, quadratischer Ansatz

10-31

Knotennummer x-Koordinate [cm]

1 0 2 25 3 50 4 62,5 5 75 6 87,5 7 100

Tabelle 10-9 Knotendatei, 3 Elemente ungleicher Länge, quadratischer Verschiebungsansatz

Elementnummer Anfangsknoten Mittenknoten Endknoten

1 1 2 3 2 3 4 5 2 5 6 7

Tabelle 10-10 Elementdatei, 3 Elemente ungleicher Länge, quadratischer Verschiebungsansatz

Zur Berechnung der Elementsteifigkeitsmatrizen nach Gl. 10-49 benötigen wir wieder die Querschnittsflächen an den Elementrändern.

Element E(e) [kN/cm2] (e)

1A [cm2] (e)3A [cm2] (e)l [cm]

1 3000 10 5.50 50 2 3000 5,50 3,25 25 3 3000 3,25 1.00 25

Tabelle 10-11 Elementgrößen

Systemgleichung

⎥⎥⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢⎢⎢

⎣

⎡ +

=

⎥⎥⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

⎥⎥⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−−

−−

−−

2083,208333,04166,08333,06250,06667,1

R4167,0

vvvvvv0

415500136085860182000122028000017515802310000010602480000015514201265

7

6

5

4

3

2

sym.


Abb. 10-22 Verschiebungen, 3 Elemente ungleicher Länge, quadratischer Ansatz

Abb. 10-23 Spannungen, 3 Elemente ungleicher Länge, quadratischer Ansatz

10-33

Knotenverschiebungen:

[ ]18353,012802,009479,007110,005228,002302,00T =v [cm]

Elementverschiebungen: Element 1: ξ(0,03979 + 0,01249ξ) [cm] Element 2: 0,05228 + 0,032763ξ + 0,009738 ξ2 [cm] Element 3: 0,09478 + 0,044211ξ + 0,044529 ξ2 [cm] Elementspannungen: Element 1: 2,3876 + 1,4983ξ [kN/cm2] Element 2: 3,9316 + 2,3372ξ [kN/cm2] Element 3: 5,3053 + 10,6870ξ [kN/cm2] Die errechneten Verschiebungen sind praktisch deckungsgleich mit den analytischen Werten. Bei den Spannungen hat sich im rechten Bereich eine wesentliche Verbesserung im Vergleich zu einer äquidistanten Elementierung mit 3 Elementen ergeben. Der Maximalwert am rechten

Rand ( 2FE cm/kN00,16=σ ) liegt näher an der analytischen Lösung, der relative Fehler be-

trägt aber immer noch etwa 20%. Eine Verbesserung der Spannungsergebnisse kann durch feinere Elementierung in der rechten Stabhälfte erzielt werden.

10.3 Statische Kondensation

In den vorangegangenen Untersuchungen wurde gezeigt, dass bei Verwendung eines Poly-nomansatzes n-ter Ordnung für die Verschiebungen genau n + 1 Freiwerte anfallen. Bei einem linearen Verschiebungsansatz waren das 2 Freiwerte, die wir den Knotenverschiebungen u1 und u2 zuordneten. Bei Verwendung eines Polynoms 2. Ordnung (quadratischer Verschie-bungsansatz, 3 Freiwerte) wurde zur Abdeckung des dritten Freiwertes ein zusätzlicher Kno-ten in Elementmitte erforderlich. Dieser Mittenknoten stellt einen inneren Knotenpunkt des Elementes dar verglichen mit den äußeren Knoten, die am Rand des Elementes liegen. Dieser Knotenfreiwert ist nur mit den äußeren Knotenwerten des Elementes verknüpft und zeigt kei-ne Wirkung auf die Nachbarelemente, was eine Folge der Ansatzfunktionen ist, die als lokale Träger nur auf der Elementebene definiert sind. Es wird deshalb im Folgenden versucht, die-sen inneren Knotenfreiwert bereits auf Elementebene durch die äußeren Knotenwerte zu er-


setzen. Dieser Vorgang wird statische Kondensation1 genannt, weil der Elimination der in-neren Knotenvariablen die Extremalbedingung des elastischen Potentials und damit dem stati-schen Gleichgewicht zugrunde gelegt ist. Ausgangspunkt unserer Untersuchungen ist der auf das Element entfallende Anteil des elasti-schen Potentials

(e)(e)T(e)(e)T(e))e(

21ˆ puuku −=Π Gl. 10-52

Durch Summation aller Elementbeiträge erhalten wir das vollständige Potential

Extremum21ˆˆ

n

1e

(e)(e)T(e)(e)(e)Tn

1e

)e( =⎟⎠⎞

⎜⎝⎛ −=Π=Π ∑∑

==

puuku Gl. 10-53

Im Elementlastvektor p(e) sind sämtliche Elementbeiträge zum Lastvektor der rechten Seite zusammengefasst. Die Variation des elastischen Potentials liefert

( ) 0ˆˆn

1e

(e)(e)(e)(e)Tn

1e

)e( =−δ=Πδ=Πδ ∑∑==

puku

Für jedes Element ist also 0][ (e)(e)(e)(e)T =−δ puku oder

(e)(e)(e) puk = Gl. 10-54

sicherzustellen. Es wird nun eine Umsortierung derart vorgenommen, dass im Elementkno-tenverschiebungsvektor u(e) die Knotenvariablen zu den äußeren Knotenpunkten im Vektor

(e)au und die Variablen zu den inneren Knotenpunkten im Vektor (e)

iu zusammengefasst sind.

Das erfordert eine Umsortierung2 der Elementsteifigkeitsmatrix k(e) und des Vektors der rech-ten Seite p(e) in

⎥⎦

⎤⎢⎣

⎡=⎥

⎦

⎤⎢⎣

⎡⎥⎦

⎤⎢⎣

⎡→= (e)

i

(e)a

(e)i

(e)a

(e)ii

(e)ia

(e)ai

(e)aa

pp

uu

kkkk

puk (e)(e)(e) ~~~ Gl. 10-55

Damit zerfällt das obige Gleichungssystem in die beiden Matrizengleichungen

1 spätl. ›Verdichtung‹ 2 was durch einfache Zeilen- und Spaltentausche immer möglich ist

10-35

)e(i

)e(i

)e(ii

)e(a

)e(ia

)e(a

)e(i

)e(ai

)e(a

)e(aa

pukuk

pukuk

=+

=+ Gl. 10-56

Aus der zweiten Gleichung kann )e(iu sofort ermittelt werden, wenn unterstellt wird, dass die

Submatrix (e)iik regulär ist

][ )e(a

)e(ia

)e(i

1)e(ii

)e(i ukpku −=

− Gl. 10-57

Setzen wir dieses Ergebnis in die erste Gleichung ein, dann erhalten wir zunächst

[ ] )e(a

)e(a

)e(ia

)e(i

1)e(ii

)e(ai

)e(a

)e(aa pukpkkuk =−+

−

und zusammengefasst

)e(i

1)e(ii

)e(ai

)e(a

)e(a

)e(ia

1)e(ii

)e(ai

)e(aa pkkpu]kkk[k −−

−=− Gl. 10-58

Der Gl. 10-58 entnehmen wir die kondensierte Elementsteifigkeitsmatrix

)e(ia

1)e(ii

)e(ai

)e(aa

(e)ˆ kkkkk a−

−= Gl. 10-59

und den kondensierten Elementlastvektor

)e(i

1)e(ii

)e(ai

)e(a

(e)ˆ pkkppa−

−= Gl. 10-60

Mit Gl. 10-59 und Gl. 10-60 können wir dann Gl. 10-58 kompakter notieren

(e))e(a

(e) ˆˆaa puk = Gl. 10-61

Wir wenden die obigen Gleichungen auf das Stabelement mit quadratischem Verschiebungs-ansatz an. Es soll zunächst die kondensierte Elementsteifigkeitsmatrix berechnet werden. Es galt mit Gl. 10-49

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

α+α+−α+α+−α+α+−α+α+−α+

=113)31(41

)31(4)1(16)3(41)3(4311

6AE

)e(1

)e((e)

lk


Das Element besitzt neben den beiden Außenknoten einen zusätzlichen Mittenknoten, den wir durch statische Kondensation auf Elementebene eliminieren wollen. Dadurch reduziert sich die Elementsteifigkeitsmatrix der Dimension ]33[ × auf eine Matrix der Dimension ]22[ × .

Der herauszukondensierende Knoten ist der Mittenknoten 2. Wir haben also in der Steifig-keitsmatrix die 2. und 3. Spalte und Zeile zu tauschen. Nach der Umsortierung erhalten wir

⎥⎦

⎤⎢⎣

⎡=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

α+α+−α+−α+−α+α+α+−α+α+

= )e(ii

)e(ia

)e(ai

)e(aa

)e(1

)e((e)

)1(16)31(4)3(4)31(41131)3(41311

6AE

kkkk

kl

Gl. 10-62

Entsprechende Umsortierungen am Knotenverschiebungsvektor und Elementlastvektor

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

3

2

1(e)

uuu

u ; ⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

141

6n e0l(e)p

liefern

⎥⎦

⎤⎢⎣

⎡=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡= )e(

i

)e(a

2

3

1(e)

uuu

uu

u ; ⎥⎦

⎤⎢⎣

⎡=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡= (e)

i

(e)a(e)

pp

p411

6n )e(

0l Gl. 10-63

Gl. 10-62 entnehmen wir die Submatrizen

⎥⎦

⎤⎢⎣

⎡α+α+α+α+

=1131

13116

AE)e(

1)e(

(e)aa

lk (e)T

ia)e(1

)e((e)ai 31

33

AE2 kk =⎥⎦

⎤⎢⎣

⎡α+α+

−=l

[ ]α+= 13

AE8)e(

1)e(

(e)ii

lk [ ])1(1

AE83

1)e(

)e(1)e(

ii α+=− lk

Gl. 10-64

und aus Gl. 10-63 folgen

[ ]13

n2;

11

6n )e(

0)e(

0 ll=⎥

⎦

⎤⎢⎣

⎡= (e)

i(e)a pp Gl. 10-65

Die kondensierte Elementsteifigkeitsmatrix ist dann

⎥⎦

⎤⎢⎣

⎡

α+α+α+α+−α+α+−α+α+

α+=−=

−

22

22

)e(1

)e()e(

ia1)e(

ii)e(

ai)e(

aa(e)

41)41()41(41

)1(3AEˆl

kkkkk a Gl. 10-66

10-37

und für den kondensierten Elementlastvektor erhalten wir

⎥⎦

⎤⎢⎣

⎡α+α+

α+=−=

−

212

)1(3nˆ

)e(0)e(

i1)e(

ii)e(

ai)e(

a(e) l

pkkppa Gl. 10-67

Die Verschiebung des Innenknotens ist

⎥⎥⎦

⎤

⎢⎢⎣

⎡α++α++

α+=−=

−

311

)e(

)e(0)e(

a)e(

ia)e(

i1)e(

ii)e(

i u)31(u)3(AE

n)1(4

1][2

lukpku Gl. 10-68

Wir elementieren den Dehnstab entsprechend Abb. 10-16 mit zwei Elementen gleicher Länge. Entsprechend Gl. 10-66 und Gl. 10-67 erhalten wir dann die Elementgrößen

Element 1: cm50;55.0105.5;cm5.5A;cm10A;cm/kN3000E 24

21

2 ===α=== l

Element 2: cm50;1818.05.51;cm1A;cm5.5A;cm/kN3000E 24

21

2 ===α=== l

⎥⎦

⎤⎢⎣

⎡=⎥

⎦

⎤⎢⎣

⎡=

⎥⎦

⎤⎢⎣

⎡=⎥

⎦

⎤⎢⎣

⎡=

961538462,0538461539,1

ˆ;8461539,1638461539,163-8461539,163-8461539,163ˆ

129032258,1370967742,1

ˆ;9354838,4519354838,451-9354838,451-9354838,451ˆ

)2(a

)2(a

)1(a

)1(a

pk

pk Gl. 10-69

Die Systemgleichung des gefesselten Systems kann dann unter Berücksichtigung der Einzel-kraft am Knoten 5 (20kN) leicht aufgebaut werden

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡ +=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

961538,20667494,2

R370968,1

vv0

8461539,1638461539,163-08461539,163-7816377,6159354838,451-

09354838,451-9354838,451

5

3

1.) Knotenverschiebungen

]cm[18022,005228,0

0

vvv

v

5

3

1

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡=

2.) Reaktionskraft

kN25370968,105228,09354838,451R −=−⋅−=


Diese Lösungen sind bereits bekannt. Die Verschiebungen der Innenknoten werden auf Ele-mentebene mit Gl. 10-68 berechnet.

cm05228,0][ )1(a

)1(ia

)1(i

1)1(ii

)1(i =−=

− ukpku

cm09571,0][ )2(a

)2(ia

)2(i

1)2(ii

)2(i =−=

− ukpku

Für den Sonderfall elementweise konstanter Querschnittsfläche A1 = A3 = A0 = konst. erhalten wir

⎥⎦

⎤⎢⎣

⎡−

−=

1111AEˆ

(e)

)e(0

(e)(e)a

lk Gl. 10-70

was Gl. 10-22 entspricht. Der Elementlastvektor ergibt sich zu

⎥⎦

⎤⎢⎣

⎡=

11

2nˆ

)e(0(e) l

ap Gl. 10-71

Hinweis: Die zur Bildung der Elementsteifigkeitsmatrizen und der Elementlastvektoren erfor-derlichen Integrationen werden in kommerziellen FE-Programmen numerisch durchgeführt (s.h. Anhang A, Numerische Integration).

F. U. Mathiak 11-1

11 Balkenelemente 11.1 Die gerade oder einachsige Biegung

Abb. 11-1 Träger mit Querbelastung q(x), EIyy = konst. (Koordinatensystem entsprechend Stabstatik) Ist der Balken statisch bestimmt gelagert, dann lassen sich die Biegemomente My(x) aus den Gleichgewichtsbedingungen allein ermitteln. Die Biegelinie kann dann aus der Differential-gleichung 2. Ordnung

yy yEI w (x) M (x)′′ = − Gl. 11-1

durch zweimaliges Integrieren gewonnen werden. Die beiden dabei anfallenden Integrations-konstanten werden aus den Randbedingungen bestimmt. Ist das Tragwerk statisch unbestimmt gelagert, so ist unter Beachtung der lokalen Gleichgewichtsbedingungen,

yz yy

zyy

dM (x)Q (x) EI (x) w(x)

dxdQ (x) q(x) EI (x) w(x)

dx

′⎡ ⎤′′= = − ⎣ ⎦

″⎡ ⎤′′= − = − ⎣ ⎦

Gl. 11-2

auf die Differentialgleichung 4. Ordnung

11-2 Balkenelement

yyEI (x) w (x) q(x)″

⎡ ⎤′′ =⎣ ⎦ Gl. 11-3

überzugehen. Ist die Biegesteifigkeit EIyy mindestens abschnittsweise konstant, dann folgt aus Gl. 11-3

IVyyEI w (x) q(x)= Gl. 11-4

11.2 Ein Balkenelement mit kubischem Verschiebungsan-satz Von der Biegelinie eines Balkens verlangen wir einen stetigen Verlauf der Verschiebung w(x). Soll die Biegelinie keine Knicke aufweisen, dann muss auch Stetigkeit der Tangenten-neigung )x(w′ gefordert werden. Um diesen Forderungen zu genügen, benötigen wir für das

Verschiebungsfeld w(x) mindestens C1-stetige Ansatzfunktionen. Wenn wir uns an die Herleitung der Steifigkeitsmatrix für ein Stabelement erinnern, dann waren den Koeffizienten des C0–stetigen Verschiebungsansatzes u(x) den Knotenwerten zu-geordnet. Bei Verwendung eines Polynoms 3. Grades mit vier Konstanten waren neben den beiden Randknoten zusätzlich zwei Mittenknoten erforderlich. Wir können bei einem 2-Knotenelement mit kubischem Verschiebungsansatz die beiden zusätzlichen Parameter jedoch auch dazu benutzen, die Stetigkeitsanforderungen an den Randknoten zu erhöhen. Wir wählen die lokalen Koordinaten entsprechend Abb. 11-2. Das Balkenelement liegt damit in der x-z-Ebene. Die Verschiebungsgrößen zeigen in positive Koordinatenrichtungen und posi-tive Drehungen drehen im Sinne der Rechtsschraubenregel um die senkrecht zur x-z-Ebene stehende y-Achse.

Abb. 11-2 Ein 2-Knotenelement

Wir wählen dazu die Knotenvariablen w und w'. Neben der Verschiebung w verfügt der Knoten eines Balkenelementes damit über einen zusätzlichen Freiheitsgrad, den Drehfrei-heitsgrad w'. Jeder Knoten besitzt also zwei und das 2-Knoten-Element damit insgesamt 4

F. U. Mathiak 11-3

Freiheitsgrade. Wir wählen als Verschiebungsansatz ein Polynom 3. Grades mit vier Freiwer-ten

2 30 1 2 3w( ) a a a aξ = + ξ + ξ + ξ Gl. 11-5

Damit werden innerhalb des Elementes die Biegemomente linear und die Querkräfte konstant approximiert. Differentiation unter Beachtung der Vorschrift

ξ

=→=ξ→=ξdd1

dxddx1dx

)e()e()e(

)e()e(

)e(

lll Gl. 11-6

liefert die Ableitungen

21 2 3(e) (e) (e)

2 2

2 3(e)2 (e)2 2 (e)2

3 3

3(e)3 (e)3 3 (e)3

IV

dw 1 dw 1w (a 2a 3a )dx dd w 1 d w 1w (2a 6a )

dx dd w 1 d w 1w 6a

dx dw 0

′ = = = + ξ + ξξ

′′ = = = + ξξ

′′′ = = =ξ

=

l l

l l

l l

Gl. 11-7

Hinweis: Wegen IVw (x) 0= ist w(x) nach Gl. 11-5 Lösung der homogen Differentialgleichung

der Balkenbiegung IVyyEI w (x) 0= .

Die vier freien Konstanten a0, a1, a2, a3 in Gl. 11-5 sind nun mit den Knotenfreiwerten zu ver-knüpfen. An den Knoten muss die Ansatzfunktion den folgenden Bedingungen genügen

11 0 1 (e)

2 0 1 2 3 2 1 2 3(e)

aw(0) w a w (0) w

1w(1) w a a a a w (1) w (a 2a 3a )

′ ′= = = =

′ ′= = + + + = = + +

l

l

Die Auflösung des obigen Gleichungssystems ergibt

2)e(

21)e(

13

2)e(

21)e(

12

1)e(

1

10

ww2ww2aww3w2w3a

wawa

′+−′+=

′−+′−−=

′=

=

ll

ll

l Gl. 11-8

Setzen wir diese Werte in Gl. 11-5 ein und sortieren nach den Knotenvariablen, dann erhalten wir

11-4 Balkenelement

2 3 (e) 2 3 2 3 (e) 3 21 1 2 2w (1 3 2 ) w ( 2 )w (3 2 )w ( )w′ ′= − ξ + ξ + ξ − ξ + ξ + ξ − ξ + ξ −ξl l Gl. 11-9

Die Funktionen vor den Knotenfreiwerten

2 3 21

(e) 2 3 (e) 22

2 3 23

(e) 2 3 (e) 24

N ( ) 1 3 2 (2 1)( 1)

N ( ) ( 2 ) ( 1)

N ( ) 3 2 (3 2 )

N ( ) ( ) ( 1)

ξ = − ξ + ξ = ξ + ξ −

ξ = ξ − ξ + ξ = ξ ξ −

ξ = ξ − ξ = ξ − ξ

ξ = −ξ + ξ = ξ ξ −

l l

l l

Gl. 11-10

sind die Formfunktionen 3. Grades. Mit Gl. 11-10 können wir die Elementverschiebungen auch in folgende Darstellung bringen

2

12

)e(2

021

11

)e(1

01

24231211

w)(Hw)(Hw)(Hw)(Hw)(Nw)(Nw)(Nw)(N)(w

′ξ+ξ+′ξ+ξ=

′ξ+ξ+′ξ+ξ=ξ

ll Gl. 11-11

Die Funktionen

0 2 1 21 10 2 1 22 2

H ( ) (2 1)( 1) H ( ) ( 1)

H ( ) (3 2 ) H ( ) ( 1)

ξ = ξ + ξ − ξ = ξ ξ −

ξ = ξ − ξ ξ = ξ ξ − Gl. 11-12

heißen Hermitesche1 Interpolationspolynome 3. Grades.

Abb. 11-3 Hermite-Polynome 3. Grades

1 Charles Hermite, frz. Mathematiker, 1822-1901

F. U. Mathiak 11-5

Sie besitzen die Eigenschaft, dass entweder ihre Ordinate oder ihre Tangentenneigung an ei-nem der Stützpunkte den Wert 1 besitzt. Alle anderen Stützwerte sind Null (Abb. 11-3). Die Herleitung der Elementsteifigkeitsmatrix und des Elementlastvektors erfolgt wieder mit Hilfe des Prinzips der virtuellen Verrückung. Dazu schreiben wir die Elementverschiebungen noch etwas um. Wir führen zunächst wieder den Vektor der Formfunktionen

[ ])x(N),x(N),x(N),x(N )e(4

)e(3

)e(2

)e(1

)e( =N Gl. 11-13

und den Vektor der Element-Knotenvariablen

(e)T1 1 2 2w w w wz ′ ′⎡ ⎤= ⎣ ⎦ Gl. 11-14

ein. Damit geht Gl. 11-11 über in die kompaktere Form

)e()e()(w zN=ξ Gl. 11-15

Ausgehend vom elastischen Potential1 für den Biegebalken

Extremumdx)x(w)x(qdx)x(wEI21AW

0x0x

2yya =−′′=−=Π ∫∫

==

ll

Gl. 11-16

liefert die Variation

0dx)x(w)x(qdx)x(w)x(wEI0x0x

yy =δ−′′δ′′=Πδ ∫∫==

ll

Aufgrund der Additivität der Energieausdrücke dürfen wir die Integration über die Balken-länge durch die Summe der Integrale über die einzelnen Elemente ersetzten, also

Extremumdx)x(w)x(qdx)x(wEI21AW

n

1e 0x

)e(

0x

)e(2yya

)e(

)e(

)e(

)e(

=⎥⎥⎦

⎤

⎢⎢⎣

⎡−′′=−=Π ∑ ∫∫

= ==

ll

Gl. 11-17

Die Variation des Funktionals ergibt

1 Im Falle zusätzlicher Belastungen in Form von Einzelkräften und Einzelmomenten ist die äußere Arbeit ent-sprechend zu ergänzen.

11-6 Balkenelement

0dxwqdxwwEIn

1e 0x

)e(

0x

)e(yy

)e(

)e(

)e(

)e(

=⎥⎥⎦

⎤

⎢⎢⎣

⎡δ−′′δ′′=Πδ ∑ ∫∫

= ==

ll

Gl. 11-18

Substituieren wir den Verschiebungsansatz Gl. 11-15 in das variierte Potential Gl. 11-18, dann erhalten wir einen Näherungswert

0dxwqdxwwEIˆn

1e 0x

)e(

0x

)e(yy

)e(

)e(

)e(

)e(

=⎥⎥⎦

⎤

⎢⎢⎣

⎡δ−′′δ′′=Πδ ∑ ∫∫

= ==

ll

Gl. 11-19

Wir benötigen in Gl. 11-19 die zweite Ableitung von w . Mit Gl. 11-15 folgt

)e()e(

2

2

1

1

2)e(

42

2)e(

32

2)e(

22

2)e(

12

)e(

wwww

dx

Nd,

dx

Nd,

dx

Nd,

dx

Nd)x(w zB=

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

′

′⎥⎦

⎤⎢⎣

⎡=′′ Gl. 11-20

mit der Ableitungsmatrix

( ) ( ) ( ) ( )[ ]ξ+−ξ+−−ξ+−ξ+−= 3122163222161 )e()e(2)e(

)e( lll

B Gl. 11-21

Einsetzen in Gl. 11-19 liefert

[ ] 0

dxqdxEI

dxqdxEIˆ

n

1e

)e()e()e(T)e(

n

1e 0

)e(T)e(

0

)e()e()e(yy

T)e(T)e(

n

1e 0

)e()e()e(

0

)e()e()e(yy

)e()e(

)e()e(

)e()e(

=−δ=

⎥⎥⎦

⎤

⎢⎢⎣

⎡−δ=

⎥⎥⎦

⎤

⎢⎢⎣

⎡δ−δ=Πδ

∑

∑ ∫∫

∑ ∫∫

=

=

=

pzkz

NzBBz

zNzBzB

ll

ll

Gl. 11-22

Wegen der Beliebigkeit der Verschiebungsvariationen (e)zδ muss für jedes Element

(e) (e) (e) 0k z - p = Gl. 11-23

erfüllt sein. In Gl. 11-22 bezeichnet

F. U. Mathiak 11-7

∫=)e(

0

)e()e(yy

T)e()e( dxEIl

BBk Gl. 11-24

die symmetrische Elementsteifigkeitsmatrix und

(e)

(e) (e)T (e)

0

q dxp N= ∫l

Gl. 11-25

den Elementlastvektor. Zur Ermittlung der Elementsteifigkeitsmatrix ist das Produkt

)e(T)e( BB zu bilden. Die Biegesteifigkeit EIyy sei innerhalb des Elementes konstant und kann deshalb vor das Integral gezogen werden. Im Einzelnen erhalten wir

dξ

ξ)31(4.)ξ6ξ51(12ξ)21(36)ξ9ξ92(4)ξ6ξ72(12ξ)32(4)ξ6ξ52(12ξ)21(36)ξ6ξ72(12ξ)21(36

EI

dxEI

1

0ξ2)e(

2)e(2

2)e(22)e(

2)e(22)e(2

3)e(

yy

0

)e()e(yy

T)e()e(

e

)e(

∫

∫

=

⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢

⎣

⎡

−+−−−+−+−−−+−−−+−−

=

=

l

l

ll

ll

l

l

sym

BBk

Gl. 11-26 Die Integration liefert

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

−−−

=

2)e(

)e(

2)e((e)2(e)

)e((e)

3(e)

yy(e)

4.612

264612612

EI

l

l

lll

ll

l

sym

k Gl. 11-27

Die Werte in den einzelnen Spalten der Steifigkeitsmatrix können wie folgt gedeutet werden. Wird zum Beispiel dem Knoten 1 die Einheitsverschiebung v1 = "1" (LE) eingeprägt (Abb.

11-4), dann filtert dieser Verschiebungsvektor die erste Spalte aus der Steifigkeitsmatrix her-aus (Gl. 11-28), deren Werte den Festhaltekräften eines beidseitig eingespannten Trägers infol-ge dieses Verformungszustandes entsprechen.

11 12 13 14 11(e)

yy21 22 23 24 1(e) (e)(e)3

31 32 33 34 2(e)

41 42 43 44 2

k k k k Vv 1(LE) 12EIk k k k M0 6

1(LE)k k k k V0 12k k k k M0 6

k z =

⎡ ⎤ ⎡ ⎤⎡ ⎤ ⎡ ⎤=⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥= ⋅ =⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥−⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥

⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎣ ⎦ ⎣ ⎦ ⎣ ⎦⎣ ⎦

l

l

l

Gl. 11-28

11-8 Balkenelement

Abb. 11-4 Einheitsverschiebung am Knoten 1, Festhaltekräfte Die zweite Spalte enthält demzufolge die Festhaltekräfte infolge einer Einheitsverdrehung am Knoten 1 usw.. Unterstellen wir auf Elementebene eine linear veränderliche Belastung q, also

( )(e)

1 2 1 1 2(e)

xq q q q (1 )q q= + − = −ξ + ξl

Gl. 11-29

dann gilt für den Elementlastvektor zunächst (e) 1

(e) (e)T (e) (e) (e)T

0 0

q dx q dp N N= = ξ∫ ∫l

l und nach

Einsetzen des Vektors der Formfunktionen N(e) und der Belastung q

{ }( )

( )

2 31 2 L1

(e)1 (e) 2 3 (e)1 2 L1(e) (e)

1 2 2 3L21 20

(e)(e) 2 3L21 2

21q 9q V1 3 23q 2q M( 2 )

(1 )q q dV9q 21q603 2

M2q 3q( )

p

⎡ ⎤+⎡ ⎤ ⎡ ⎤− ξ + ξ⎢ ⎥⎢ ⎥ ⎢ ⎥+ξ − ξ + ξ ⎢ ⎥⎢ ⎥ ⎢ ⎥= − ξ + ξ ξ = =⎢ ⎥⎢ ⎥ ⎢ ⎥+ξ − ξ ⎢ ⎥⎢ ⎥ ⎢ ⎥− +−ξ + ξ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎣ ⎦⎣ ⎦ ⎣ ⎦

∫ll l

l

ll

Gl. 11-30

Abb. 11-5 Beidseitig eingespannter Träger, äquivalente Knotenlasten Für den Sonderfall konstanter Elementbelastung ( 021 qqq == ) geht Gl. 11-30 über in

F. U. Mathiak 11-9

L1(e)(e)

L1(e) 0

L2(e)

L2

V6MqV612

M

p

⎡ ⎤⎡ ⎤⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥= =⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥

−⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦

ll

l

Gl. 11-31

Hinweis: Der Elementlastvektor p(e) enthält die statisch äquivalenten Knotenlasten (Einzel-kräfte VL1, VL2 und Einzelmomente ML1, ML2), die sich als Auflagerkräfte und Vollein-spannmomente1 eines beidseitig eingespannten Trägers unter Querbelastung q(x) ergeben.

11.2.1 Beispiel 10.1 Wir testen das soeben entwickelte Element am Beispiel der Abb. 11-6. Die analytische Lösung

4 32 20 0

yy yy

20 0

q qw(x) (3 2 )(1 ) w (x) (8 15 6)48EI 48EI

q qM(x) (1 )(1 4 ) Q(x) (5 8 )8 8

′= ξ − ξ −ξ = ξ ξ − ξ +

= − −ξ − ξ = − ξ

l l

l l Gl. 11-32

ist bekannt.

Abb. 11-6 Elementierung eines Balkens, 2 Elemente gleicher Länge

Die Verschiebungen entsprechen einem Polynom 4. Grades, sie lassen sich also mit unserem kubischen Verschiebungsansatz nicht exakt wiedergeben. Das gilt dann auch für die Biege-momente und die Querkräfte.

1 hierbei ist die Änderung der Vorzeichenregel gegenüber der aus der Statik bekannten Festsetzung zu beachten

11-10 Balkenelement

Die Diskretisierung des Balkens erfolgt durch zwei Elemente gleicher Länge. Wir benötigen wieder eine Knoten- und eine Elementdatei.

Knotennummer x-Koordinate [m]

1 0 2 5,00 3 10,0



1 1 2 2 2 3


Systemwerte: m5)e( =l , 2)e( m/MN30000E = , Balken mit Rechteckquerschnitt b/h =

40/80cm ( 4yyI 0,0171m= ), 2

yy MNm512EI = , m/MN2,0qq 0 == .

Da die Systemwerte für beide Elemente gleich sind, sind auch beide Elementsteifigkeitsmatri-zen und Elementlastvektoren gleich.

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

−−−

==

600,409.880,122152,49800,204880,122600,409880,122152,49880,122152,49

)2()1(

sym

kk

⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢

⎣

⎡

−

==

4167,05,0

4167,05,0

)2()1( pp

Zum Aufbau der globalen Steifigkeitsmatrix und des globalen Lastvektors benötigen wir den Zusammenhang zwischen den lokalen Elementknotenfreiheitsgraden und den globalen Sys-temknotenfreiheitsgraden, die wir im Vektor

[ ]332211 vvvvvv ′′′=Tv Gl. 11-33

zusammenfassen. Die Einordnung der Elementsteifigkeitsmatrizen und der Elementlastvekto-ren in die Systemmatrix bzw. die rechte Seite erfolgt mit Hilfe der Knoten- und Elementdatei. Für die Steifigkeitsmatrix ergibt sich folgende Parkettierung

F. U. Mathiak 11-11

Abb. 11-7 Bildung der Gesamtsteifigkeitsmatrix

⎥⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−

+−+

=

⎥⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

′

′

′

⎥⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−−+−+−+

−−

4167,05,0

4167,04167,05,05,0

4167,05,0

vvvvvv

600,409800,122152,49

800,204880,122600,409600,409880,122152,49880,122880,122152,49152,49

00800,204880,122600,40900880,122152,49880,122152,49

3

3

2

2

1

1

sym. Fassen wir die Komponenten in der obigen Gleichung zusammen, dann erhalten wir

⎥⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−

=

⎥⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

′

′

′

⎥⎥⎥⎥⎥⎥⎥⎥

⎦

⎤

⎢⎢⎢⎢⎢⎢⎢⎢

⎣

⎡

−−−

−−

4167,05000,00

0000,14167,05000,0

vvvvvv

600,409.800,122152,49

800,204880,122200,819880,122152,490304,98

00800,204880,122600,40900880,122152,49880,122152,49

3

3

2

2

1

1

sym

Es fehlt noch der Einbau der Randbedingungen. Am linken Rand müssen 0vv 11 =′= erfüllt

sein, und am rechten Rand ist 0v3 = zu fordern. Die aus diesen kinematischen Zwängen re-

sultierenden Reaktionslasten V1, M1 und V3 erscheinen dann auf der rechten Seite im Vektor der Unbekannten.

2

2

3

49,152 122,880 49,152 122,880 0 0 0 0,5000409,600 122,880 204,800 0 0 0 0,4167

98,304 0 49,152 122,880 v 1,0000819,200 122,880 204,800 v 0

49,152 122,800 0 0,5000. 409,600 v 0sym

⎡ ⎤− ⎡ ⎤⎢ ⎥ ⎢ ⎥−⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥−⎢ ⎥ =⎢ ⎥− ′⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥−⎢ ⎥ ⎢ ⎥

′ −⎢ ⎥ ⎢ ⎥⎣ ⎦⎣ ⎦

1

1

3

VM00

V,4167 0

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥

+⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦

Wegen 0vvv 311 ==′= und 0δvvδδv 311 ==′= können die 1., 2., 5. Zeile und Spalte aus

dem Gleichungssystem gestrichen werden. Es verbleibt

11-12 Balkenelement

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

−=

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

′′

⎥⎥⎥

⎦

⎤

⎢⎢⎢

⎣

⎡

4167,00000,00000,1

vvv

600,409.800,204200,819880,1220304,98

3

2

2

sym

mit der Lösung: 00813818,0v00203454,0vm02034526,0v 322 −=′=′=

Damit sind auch die Reaktionslasten bekannt:

MN25,15,0v880,122v125,49V 221 −=−′⋅+⋅−=

MNm5,24167,0v800,204v880,122M 221 −=−′⋅+⋅−=

MN75,05000,0v800,122v880,122v125,49V 3223 −=−′⋅−′⋅−⋅−=

11.2.1.1 Rückrechnung Mit dem Vektor v der globalen Knotenwerte liegen sämtliche Zustandsgrößen fest. Für beide Elemente gilt:

]55, 23, 5105, 231[NN 32323232)2()1( ξ+ξ−ξ−ξξ+ξ−ξξ+ξ−==

]20.140.0, 48.024.0, 20.180.0, 48.024.0[-BB )2()1( ξ+−ξ−ξ+−ξ+==

Element 1 (1)

1

1(1) (1) (1)1 2 3 4 1 2 3 4

2

2

2 23 4

2 3

w 0w 0

w N , N , N , N N , N , N , Nw 0,0203453w 0,0020345

0.0203453N 0.0020345N 0.0203453 (2 3) 0.0020345 5 (1 )

0.0508634 - 0.0305181 [m]

N z

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥′⎢ ⎥ ⎢ ⎥⎡ ⎤ ⎡ ⎤= = =⎣ ⎦ ⎣ ⎦⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥′ ⎢ ⎥⎢ ⎥ ⎣ ⎦⎣ ⎦

= + = − ξ ξ − − ⋅ ξ − ξ

= ξ ξ

( ) ( )

(1)1

1(1) (1) (1)y yy yy 1 2 3 4 yy 1 2 3 4

2

2

yy

w 0w 0

M EI B z EI B , B , B , B EI B , B , B , Bw 0,0203453w 0,0020345

EI 0.0203453 0.24 1 0.48 0.0020345 0.40 1.20

2.08334 3.75006 [MNm]

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥′⎢ ⎥ ⎢ ⎥⎡ ⎤ ⎡ ⎤= − = − = −⎣ ⎦ ⎣ ⎦⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥′ ⎢ ⎥⎢ ⎥ ⎣ ⎦⎣ ⎦

⎡ ⎤= − ⋅ − ξ + − + ξ⎣ ⎦= − + ξ

(1)(1)

zdM 3.75006Q 0.75 [MN]

dx 5.00= = =

Element 2

F. U. Mathiak 11-13

(2)1

1(2) (2) (2)1 2 3 4 1 2 3 4

2

2

1 2 42 3 2 3

w 0,0203453w 0,0020345

w N , N , N , N N , N , N , Nw 0w 0,0081382

0.0203453N 0.0020345N 0,0081382N

0.0203453(1 3 2 ) 0.0020345(5 10 5 ) 0,00

N z

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥′⎢ ⎥ ⎢ ⎥⎡ ⎤ ⎡ ⎤= = =⎣ ⎦ ⎣ ⎦⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥′ −⎢ ⎥⎢ ⎥ ⎣ ⎦⎣ ⎦

= + −

= − ξ + ξ + ξ − ξ + ξ − 2 3

2 3

81382( 5 5 )0.0203453 0.0101725 0.04069000 0.01017220

− ξ + ξ

= + ξ − ξ + ξ

Abb. 11-8 Verschiebung w, 2 Elemente gleicher Länge

Abb. 11-9 Biegemomente M, 2 Elemente gleicher Länge

11-14 Balkenelement

Abb. 11-10 Querkräfte Q, 2 Elemente gleicher Länge

(2)1

1(2) (2) (2)y yy yy 1 2 3 4 yy 1 2 3 4

2

2

yy 1 2 4

w 0,0203453w 0,0020345

M EI B z EI B ,B ,B ,B EI B ,B ,B ,Bw 0w 0,0081382

EI 0,0203453 B 0.0020345 B 0,0081382 B

1,666662400 1,249959

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥′⎢ ⎥ ⎢ ⎥⎡ ⎤ ⎡ ⎤= − = − = −⎣ ⎦ ⎣ ⎦⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥′ −⎢ ⎥⎢ ⎥ ⎣ ⎦⎣ ⎦

⎡ ⎤= − ⋅ + ⋅ − ⋅⎣ ⎦= − 936 [MNm]ξ

(2)(2)

zdM 1,249959936Q 0.25 [MN]

dx 5.00−

= = = −

Offensichtlich werden trotz der groben Diskretisierung die Verschiebungen recht gut wieder-gegeben. Das gilt nicht für die Biegemomente My und schon gar nicht für die Querkräfte Qz. Für eine Bemessung sind die Ergebnisse der Schnittlasten unbrauchbar. Um hier zu akzeptab-leren Lösungen zu kommen, müsste, bei Beibehaltung des Elementes, wesentlich feiner ele-mentiert werden.

11.3 Ein Balkenelement mit quintischem Verschiebungsan-satz Neben der feineren Elementierung besteht eine weitere Möglichkeit der Ergebnisverbesserung in der Erhöhung des Polynomgrades der Verschiebungsfunktion w. Bei der Wahl eines Ver-schiebungsansatzes mit sechs Freiwerten

2 3 4 50 1 2 3 4 5w( ) a a a a a aξ = + ξ + ξ + ξ + ξ + ξ Gl. 11-34

F. U. Mathiak 11-15

könnten wir die beiden zusätzlichen Freiwerte dazu benutzen, die Stetigkeitsforderungen an beiden Randknoten zu erhöhen, also neben w und w ′ zusätzlich Stetigkeit in w ′′ zu fordern.

Abb. 11-11 Stetigkeiten am Übergang zweier benachbarter Elemente

Damit hätten wir als zusätzlichen Knotenfreiwert die linearisierte Krümmung w ′′ , was aber wegen y yyM EI w′′= − nur dann auch Stetigkeit in My bedeutet, wenn die Biegesteifigkeiten

EIyy zweier benachbarter Elemente gleich sind (Abb. 11-11). Es macht also mechanisch keinen großen Sinn, ein solches Element weiter zu betrachten, besser ist es, die beiden hinzukom-menden Freiwerte einem zusätzlichen Knoten zuordnen, etwa dem Elementmittelpunkt (Abb.

11-12) und jedem der Knoten die Variablen w und w' zuzuordnen.

Abb. 11-12 Balkenelement mit drei Knoten

Mit dem Ansatz Gl. 11-34 werden auf Elementebene die Biegemomente kubisch und die Quer-kräfte quadratisch approximiert. Das ist eine wesentliche Verbesserung gegenüber dem kubi-schen Ansatz. Für den weiteren Rechengang benötigen wir die Ableitungen der Verschie-bungsfunktion

11-16 Balkenelement

2 3 41 2 3 4 5(e) (e) (e)

2 22 3

2 3 4 5(e)2 (e)2 2 (e)2

3 32

3 4 5(e)3 (e)3 3 (e)3

4 4IV

4(e)4 (e)4 4 (e)4

dw 1 dw 1w (a 2a 3a 4a 5a )ddx

d w 1 d w 1w (2a 6a 12a 20a )dx dd w 1 d w 1w (6a 24a 60a )

dx dd w 1 d w 1w (24a

dx d

′ = = = + ξ + ξ + ξ + ξξ

′′ = = = + ξ + ξ + ξξ

′′′ = = = + ξ + ξξ

= = =ξ

l l

l l

l l

l l5120a )+ ξ

Gl. 11-35

Die 6 freien Konstanten a0, a1, a2, a3, a4 und a5 sind mit den Knotenfreiwerten zu verknüpfen. An den 3 Knoten müssen die Ansatzfunktion sowie deren 1. Ableitung den Knotenvariablen entsprechen. Das führt auf das lineare Gleichungssystem

1 0

11 (e)

2 0 1 2 3 4 5

2 1 2 3 4 5(e)

3 0 1 2 3 4 5

3 1 2 3 4 5(e)

w(0) w aaw (0) w

1 1 1 1 1w(1/ 2) w a a a a a a2 4 8 16 32

1 3 1 5w (1/ 2) w (a a a a a )4 2 16

w(1) w a a a a a a1w (1) w (a 2a 3a 4a 5a )

= =

′ ′= =

= = + + + + +

′ ′= = + + + +

= = + + + + +

′ ′= = + + + +

l

l

l

Gl. 11-36

aus dem die Konstanten berechnet werden. Einsetzen der Konstanten in Gl. 11-34 und Koeffi-zientenvergleich in den Knotenvariablen liefert die Formfunktionen 5. Grades für unser Balkenelement

2 2 (e) 2 21 2

2 2 (e) 2 23 4

2 2 (e) 2 25 6

N ( ) (6 1)(2 1) ( 1) N ( ) (2 1) ( 1)

N ( ) 16 ( 1) N ( ) 8 (2 1)( 1)

N ( ) (6 7)(2 1) N ( ) (2 1) ( 1)

ξ = ξ + ξ − ξ − ξ = ξ ξ − ξ −

ξ = ξ ξ − ξ = ξ ξ − ξ −

ξ = −ξ ξ − ξ − ξ = ξ ξ − ξ −

l

l

l

Gl. 11-37

Die Elementverschiebungen gehen dann mit Gl. 11-34 über in

1 1 2 1 3 2 4 2 5 3 6 30 (e) 1 0 (e) 1 0 (e) 11 1 1 1 2 2 2 2 3 3 3 3

w( ) N ( )w N ( )w N ( )w N ( )w N ( )w N ( )w

H ( )w H ( )w H ( )w H ( )w H ( )w H ( )w

′ ′ ′ξ = ξ + ξ + ξ + ξ + ξ + ξ

′ ′ ′= ξ + ξ + ξ + ξ + ξ + ξl l l Gl. 11-38

Die Funktionen

F. U. Mathiak 11-17

0 2 2 0 2 2 0 2 21 2 31 2 2 1 2 2 1 2 21 2 3

H ( ) (6 1)(2 1) ( 1) H ( ) 16 ( 1) H ( ) (6 7)(2 1)

H ( ) (2 1) ( 1) H ( ) 8 (2 1)( 1) H ( ) (2 1) ( 1)

ξ = ξ + ξ − ξ − ξ = ξ ξ − ξ = −ξ ξ − ξ −

ξ = ξ ξ − ξ − ξ = ξ ξ − ξ − ξ = ξ ξ − ξ − sind die Hermiteschen Interpolationspolynome 5. Grades (Abb. 11-13 und Abb. 11-14).

Abb. 11-13 Hermite-Polynome 5. Grades

Abb. 11-14 Hermite-Polynome 5. Grades Mit Einführung des Vektors der Formfunktionen

1 2 3 4 5 6N , N , N , N , N , N(e)N ⎡ ⎤= ⎣ ⎦

und des Vektors der Knotenvariablen

11-18 Balkenelement

1 1 2 2 3 3w w w w w w(e) Tz ⎡ ⎤′ ′ ′= ⎢ ⎥⎣ ⎦ Gl. 11-39

lassen sich die Elementverschiebungen wie folgt in Matrizenschreibweise notieren

w (e) (e)N z= Gl. 11-40

Der weitere Rechengang entspricht dem des 2-Knotenelementes. Unter Beachtung von Gl.

11-24 und

1

12 2 22 2 2

2 (e) (e)3 5 61 2 4(e)2 (e)2 (e)2 (e)2 (e)2 (e)2

2

3

3

wwwd N d N d Nd N d N d N

w , , , , ,wdx dx dx dx dx dxww

B z

⎡ ⎤⎢ ⎥′⎢ ⎥⎢ ⎥⎡ ⎤

′′ = =⎢ ⎥⎢ ⎥ ′⎢ ⎥ ⎢ ⎥⎣ ⎦⎢ ⎥⎢ ⎥

′⎢ ⎥⎣ ⎦

Gl. 11-41

und der Ableitungsmatrix 2 2 22 2 2

(e) 3 5 61 2 4(e)2 (e)2 (e)2 (e)2 (e)2 (e)2

d N d N d Nd N d N d N, , , , ,

dx dx dx dx dx dxB

⎡ ⎤= ⎢ ⎥⎢ ⎥⎣ ⎦

mit

(e) 2 31 (e)2

2B (-23 198 - 408 240 )= + ξ ξ + ξl

(e) 2 32 (e)

2B (-6 39 - 72 40 )= + ξ ξ + ξl

(e) 23 (e)2

32B (1 6 6 )= − ξ + ξl

(e) 2 34 (e)

16B (-1 12 - 30 20 )= + ξ ξ + ξl

(e) 2 35 (e)2

2B (7 102 312 240 )= − ξ + ξ − ξl

(e) 2 36 (e)

2B (-1 15 - 48 40 )= + ξ ξ + ξl

Gl. 11-42

erhalten wir unter Beachtung der Formfunktionen nach Gl. 11-37 die symmetrische Element-steifigkeitsmatrix, die jetzt die Größe [ 6 6× ] besitzt.

(e) (e) (e)

(e)2 (e) (e)2 (e) (e)2

(e)yy(e)3 (e)2 (e) (e)2

(e)

(e)2

5092 1138 3584 1920 1508 242332 896 320 242 38

EI 7168 0 3584 89635 1280 1920 320

5092 1138332

(e)k

sym.

⎡ ⎤− −⎢ ⎥

− −⎢ ⎥⎢ ⎥−

= ⎢ ⎥−⎢ ⎥

⎢ ⎥−⎢ ⎥⎢ ⎥⎣ ⎦

l l l

l l l l l

l

l l l l

l

l

Gl. 11-43

Zur Ermittlung des Elementlastvektors

F. U. Mathiak 11-19

(e )

(e) (e)T (e) (e)

0

q dxp N= ∫l

Gl. 11-44

interpolieren wir die Linienlast q durch eine quadratische Verteilung mittels der Lagrange-schen Interpolationspolynome

(e)1 2 3q (2 1)( 1)q 4 ( 1)q (2 1)q= ξ − ξ − − ξ ξ − + ξ ξ − Gl. 11-45

die für 2 1 31q (q q )2

= + in die lineare Verteilung

(e)

1 3q (1 )q q= − ξ + ξ Gl. 11-46

übergeht.

Abb. 11-15 Elementlasten Die Auswertung der Gl. 11-44 liefert die Elementlastvektoren für eine quadratische Lastvertei-lung

L1(e) (e)

L11(e)

L2(e) (e) 2

L23

L3(e) (e)

L3

57 44 3 V3 4 0 M

q12 192 16 Vq8 0 8 M420q3 44 57 V

0 4 4 M

(e)p

⎡ ⎤− ⎡ ⎤⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎡ ⎤⎢ ⎥ ⎢ ⎥⎢ ⎥⎢ ⎥= = ⎢ ⎥⎢ ⎥−⎢ ⎥ ⎢ ⎥⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦−⎢ ⎥ ⎢ ⎥

− −⎢ ⎥ ⎢ ⎥⎣ ⎦⎣ ⎦

l l

l

l l

l l

Gl. 11-47

sowie für eine lineare Lastverteilung.

11-20 Balkenelement

( )

( )

( )

1 3 L1(e)

1 3 L1(e)

1 3 L2(e)

1 3 L2

L31 3(e)

L31 3

79q 19q V5q 2q M

112(q q ) V8 q q M420

V19q 79qM2q 5q

(e)p

⎡ ⎤+ ⎡ ⎤⎢ ⎥ ⎢ ⎥+⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥+⎢ ⎥= = ⎢ ⎥− −⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥+⎢ ⎥ ⎢ ⎥⎢ ⎥− + ⎢ ⎥⎣ ⎦⎣ ⎦

l

l

l

l

Gl. 11-48

Die Elementlasten VLi und MLi (i = 1,3) sind die statisch äquivalenten Knotenlasten (Kräfte V und Momente M).

Abb. 11-16 Beidseitig eingespannter Träger, äquivalente Knotenlasten Für den praktisch wichtigen Fall konstanter Elementbelastung ( 021 qqq == ) geht Gl. 11-48

über in

L1(e)

L1(e) (e)

L20

L2

L3(e)

L3

V14MV32q

M060V14

M

(e)p

⎡ ⎤⎡ ⎤⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥

= = ⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥

− ⎢ ⎥⎢ ⎥⎣ ⎦ ⎣ ⎦

l

l

l

Gl. 11-49

11.3.1 Beispiel 10.2

F. U. Mathiak 11-21

Abb. 11-17 Elementierung eines Balkens, 2 Elemente gleicher Länge

Wir testen das Element aus Kap.11.3 am Beispiel der Abb. 11-17. Der Grad des Polynoms Gl.

11-34 reicht aus, die Zustandsgrößen eine linear veränderliche Belastung exakt wiederzugeben. Bei einer konstanten Belastung entspricht die analytische Verschiebungsfunktion nämlich einem Polynom 4. Grades. Wir benötigen deshalb nur ein Element.

Knotennummer x-Koordinate [m]

1 0 2 5,00 3 10,00


Elementnummer Anfangsknoten Mittenknoten Endknoten

1 1 2 3


Systemwerte: (e) 10m=l , 2)e( m/MN30000E = , Balken mit Rechteckquerschnitt b/h =

40/80cm ( 4yyI 0,0171m= ), 2

yy MNm512EI = , m/MN2,0qq 0 == .

Die Knotenfreiheitsgrade fassen wir im Vektor der globalen Systemfreiheitsgrade

1 1 2 2 3 3v v v v v vTv ′ ′ ′⎡ ⎤= ⎣ ⎦ Gl. 11-50

zusammen. Mit den obigen Werten ergibt sich die Systemgleichung des ungefesselten Sys-tems

11-22 Balkenelement

(e) (e) (e)

1(e)2 (e) (e)2 (e) (e)2

1(e)

yy 2(e)3 (e)2 (e) (e)2

2(e)

3(e)2

3

v5092 1138 3584 1920 1508 242v332 896 320 242 38

EI v7168 0 3584 896v35 1280 1920 320v5092 1138v332sym.

⎡ ⎤ ⎡− −⎢ ⎥ ⎢ ′− −⎢ ⎥ ⎢⎢ ⎥ ⎢−⎢ ⎥ ⎢

′−⎢ ⎥ ⎢⎢ ⎥ ⎢−⎢ ⎥ ⎢

′⎢ ⎥ ⎢⎣⎣ ⎦

l l l

l l l l l

l

l l l l

l

l

(e)

(e)0

(e)

14

32q060

14

⎤ ⎡ ⎤⎥ ⎢ ⎥⎥ ⎢ ⎥⎥ ⎢ ⎥=⎥ ⎢ ⎥

⎥ ⎢ ⎥⎥ ⎢ ⎥⎥ ⎢ ⎥

−⎥ ⎢ ⎥⎣ ⎦⎦

l

l

l

Gl. 11-51 Am linken Rand müssen 0vv 11 =′= erfüllt sein, und am rechten Rand ist 3v 0= zu fordern.

Die aus diesen kinematischen Zwängen resultierenden Reaktionslasten V1, M1 und V3 er-scheinen dann auf der rechten Seite im Vektor der Unbekannten.

(e) (e) (e)

(e)2 (e) (e)2 (e) (e)2

(e)yy 2(e)3 (e)2 (e) (e)2

2(e)

(e)23

05092 1138 3584 1920 1508 2420332 896 320 242 38

EI v7168 0 3584 896v35 1280 1920 32005092 1138v332sym.

⎡ ⎤ ⎡ ⎤− −⎢ ⎥ ⎢ ⎥

− −⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢−⎢ ⎥ ⎢

′−⎢ ⎥ ⎢⎢ ⎥ ⎢−⎢ ⎥ ⎢

′⎢ ⎥ ⎢⎣ ⎦⎣ ⎦

l l l

l l l l l

l

l l l l

l

l

1(e)

1(e)

0

3(e)

14 VM

32 0q0 060

14 V0

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥

⎥ ⎢ ⎥ ⎢ ⎥= +⎥ ⎢ ⎥ ⎢ ⎥

⎥ ⎢ ⎥ ⎢ ⎥⎥ ⎢ ⎥ ⎢ ⎥⎥ ⎢ ⎥ ⎢ ⎥

−⎥ ⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦

l

l

l

Gl. 11-52

Wegen 1 1 3v v v 0′= = = und 1 1 3δv δv δv 0′= = = können die 1., 2., 5. Zeile und Spalte aus dem

Gleichungssystem Gl. 11-51 gestrichen werden. Es verbleibt das reduzierte Gleichungssystem

(e)2 (e)

yy (e)2 (e)2 02(e)3

(e)2 (e)3

7168 0 896 v 32EI q

1280 320 v 06035

. 332 vsym

⎡ ⎤ ⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥ ⎢ ⎥′ =⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥′ −⎣ ⎦ ⎣ ⎦⎣ ⎦

ll

l ll

l l

mit der Lösung: 4 3 3

0 0 02 2 3

yy yy yy

q q qv v v192EI 192EI 48EI

′ ′= = = −l l l

F. U. Mathiak 11-23

Der System-Knotenverschiebungsvektor ist dann

1

13

2 0

2 yy

3

3

v 0v 0v qv 1192EIv 0v 4

v

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥′⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥

= =⎢ ⎥ ⎢ ⎥′⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥′ −⎢ ⎥ ⎢ ⎥⎣ ⎦⎣ ⎦

ll

Aus Gl. 11-52 ermitteln wir die Auflagerreaktionsgrößen

21 0 1 0 3 0

5 1 3V q M q V q8 8 8

= − = − = −l l l

11.3.1.1 Rückrechnung

Die Verschiebungen ermitteln wir mittels Gl. 11-40. Mit (e) 10m=l erhalten wir zunächst

2 2

2 2

2 2

2 2

2 2

2 2

(1 6 )(2 -1) ( -1)10 (2 -1) ( -1)

16 ( -1)80 (2 -1)( -1)

(6 - 7)(2 -1)10 ( -1)(2 -1)

TN

⎡ ⎤+ ξ ξ ξ⎢ ⎥

ξ ξ ξ⎢ ⎥⎢ ⎥ξ ξ

= ⎢ ⎥ξ ξ ξ⎢ ⎥

⎢ ⎥−ξ ξ ξ⎢ ⎥ξ ξ ξ⎢ ⎥⎣ ⎦

und damit 4

2 20

yy

qw (3 5 2 )

48EI= ξ − ξ + ξ

l

Das ist die theoretisch exakte Lösung. Die Schnittlasten folgen in bekannter Weise durch Dif-ferentiationsoperationen an der Verschiebungsfunktion.

11.3.2 Statische Kondensation Durch statische Kondensation lässt sich der Mittenknoten eliminieren. Dazu wird in Ana-

logie zu Kap. 9.3 eine Portionierung der Gleichung (e) (e) (e)k z p= vorgenommen:

11-24 Balkenelement

(e) (e) (e)aa ai a(e) (e) (e)ia ii i

(e)a(e)i

k k z pk k z p

⎡ ⎤ ⎡ ⎤ ⎡ ⎤=⎢ ⎥ ⎢ ⎥ ⎢ ⎥

⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦ ⎣ ⎦ Gl. 11-53

Dabei sind im Vektor (e)az die Verschiebungskomponenten der Außenknoten und im Vektor

(e)iz die Verschiebungskomponenten des Innenknotens zusammengefasst. Eine entsprechende

Umsortierung ist in der Elementsteifigkeitsmatrix k(e) und dem Elementlastvektor p(e) vorzu-nehmen. Wir erhalten im Einzelnen, wenn wir uns auf eine linear veränderliche Querlastver-teilung q beschränken

1

1

3

3

2

2

w

ww

ww

w

(e)a(e)(e)i

zz

z

⎡ ⎤⎢ ⎥

′⎢ ⎥⎢ ⎥⎡ ⎤ ⎢ ⎥= =⎢ ⎥ ⎢ ⎥′⎢ ⎥⎣ ⎦ ⎢ ⎥⎢ ⎥⎢ ⎥

′⎢ ⎥⎣ ⎦

1 3 1(e)

1 3 1(e) (e)a 1 3 3(e) (e)i 1 3 3

1 3 2(e)

1 3 2

79q 19q V(5q 2q ) M

19q 79q V(2q 5q ) M420

112(q q ) V8 (q q ) M

(e) pp

p

⎡ ⎤ ⎡ ⎤+⎢ ⎥ ⎢ ⎥+⎢ ⎥ ⎢ ⎥

⎡ ⎤ ⎢ ⎥ ⎢ ⎥+= = =⎢ ⎥ ⎢ ⎥ ⎢ ⎥

− +⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎣ ⎦⎢ ⎥ ⎢ ⎥+⎢ ⎥ ⎢ ⎥− −⎢ ⎥ ⎢ ⎥⎣ ⎦⎣ ⎦

l

l

l

l

(e) (e)

(e)2 (e) (e)2yy(e)3 (e)

(e)2

5092 1138 1508 242EI 332 242 38

35 5092 1138. 323

(e)aak

sym

⎡ ⎤−⎢ ⎥

−⎢ ⎥= ⎢ ⎥−⎢ ⎥⎢ ⎥⎣ ⎦

l l

l l l

l l

l

(e)

(e) (e)2yy(e)3 (e)

(e) (e)2

3605 1920EI 896 320

35 3605 1920896 320

(e) (e)Tai iak k

⎡ ⎤−⎢ ⎥−⎢ ⎥= =⎢ ⎥− −⎢ ⎥⎢ ⎥⎣ ⎦

l

l l

l l

l l

yy(e)2(e)3

EI 7168 00 128035

(e)iik

⎡ ⎤= ⎢ ⎥

⎣ ⎦ll

(e)31

yy(e)2

1 035 7168

1EI 01280

(e)iik −

⎡ ⎤⎢ ⎥⎢ ⎥=⎢ ⎥⎢ ⎥⎣ ⎦

l

l

Hinweis: Zur Berechnung der kondensierten Elementsteifigkeitsmatrix und des kondensierten

Elementlastvektors wird die Inverse von (e)iik benötigt, was offensichtlich immer möglich ist,

da die Determinante von (e)iik nicht singulär ist.

Entsprechend Kap. 9 folgt dann die kondensierte Elementsteifigkeitsmatrix

F. U. Mathiak 11-25

(e) (e)

(e)2 (e) (e)2yy(e) (e) (e) (e) 1 (e)

aa ai ii ia (e)(e)3

(e)2

12 6 12 64 6 2EIˆ

12 6. 4

ak k k k k

sym

−

⎡ ⎤−⎢ ⎥

−⎢ ⎥= − = ⎢ ⎥−⎢ ⎥⎢ ⎥⎣ ⎦

l l

l l l

ll

l

Gl. 11-54

die mit der Steifigkeitsmatrix des 2-Knotenelementes übereinstimmt, und der kondensierte Elementlastvektor

1 3(e)(e)

1 3(e) (e) (e) (e) 1 (e)a ai ii i

1 3(e)

1 3

3(7q 3q )(3q 2q )

ˆ3(3q 7q )60(2q 3q )

ap p k k p−

⎡ ⎤+⎢ ⎥+⎢ ⎥= − = ⎢ ⎥+⎢ ⎥− +⎢ ⎥⎣ ⎦

ll

l

Gl. 11-55

der im Falle konstanter Belastung 1 3 0q q q= = übergeht in

(e)(e)(e) 0

(e)

6qˆ

612ap

⎡ ⎤⎢ ⎥⎢ ⎥=⎢ ⎥⎢ ⎥−⎢ ⎥⎣ ⎦

ll

l

Gl. 11-56

Mit Gl. 11-54 und Gl. 11-55 können wir dann kompakter notieren

(e) (e) (e)a

ˆ ˆa ak z p= Gl. 11-57

Die Verformung des Innenknotens (e) (e) 1 (e) (e) (e)i ii i ia a[ ]z k p k u−= − errechnet sich zu

1(e)(e) (e)

(e) 3 1 32 1(e)i

3yy2(e) (e)1 3

3

w4 4(q q )w w1768

12 12 wEI 81w 2 2(q q )1920 w

z

⎡ ⎤⎡ ⎤⎡ ⎤ ⎢ ⎥−+ ⎢ ⎥⎢ ⎥ ′⎡ ⎤ ⎢ ⎥⎢ ⎥⎢ ⎥= = +⎢ ⎥ ⎢ ⎥′ ⎢ ⎥⎢ ⎥⎢ ⎥ ⎢ ⎥− − −⎣ ⎦ − − ⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ′⎣ ⎦ ⎣ ⎦

l l ll

l l

Gl. 11-58

und für den Sonderfall 1 3 0q q q= =

11-26 Balkenelement

1(e) (e)

(e) 42 1(e) 0

i3yy2

(e) (e)

3

w4 4w w1q 1

12 12 w0384EI 8w 2 2w

z

⎡ ⎤⎡ ⎤ ⎢ ⎥−⎢ ⎥ ′⎡ ⎤ ⎡ ⎤ ⎢ ⎥⎢ ⎥= = +⎢ ⎥ ⎢ ⎥ ⎢ ⎥′ ⎢ ⎥⎢ ⎥ ⎣ ⎦ ⎢ ⎥− − −⎣ ⎦ ⎢ ⎥ ⎢ ⎥′⎣ ⎦ ⎣ ⎦

l ll

l l

Gl. 11-59

11.3.3 Beispiel 10.3 Wir testen das soeben entwickelte Element am Beispiel des Trägers nach Abb. 11-17 und ver-wenden wieder nur ein Element über die gesamte Trägerlänge. Da wir nur ein Element ver-wenden, stimmt die kondensierte Elementsteifigkeitsmatrix mit der Systemsteifigkeitsmatrix überein. Das zu lösende Gleichungssystem ist dann

1

2 2yy 1 03

32

3

v12 6 12 6 64 6 2EI v q

12 6 v 612. 4 vsym

⎡ ⎤⎡ ⎤− ⎡ ⎤⎢ ⎥⎢ ⎥ ⎢ ⎥′− ⎢ ⎥⎢ ⎥ ⎢ ⎥=⎢ ⎥⎢ ⎥ ⎢ ⎥− ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ −′ ⎢ ⎥⎢ ⎥ ⎣ ⎦⎣ ⎦ ⎣ ⎦

l l

l l l ll

ll

l l

Der Einbau der homogenen Randbedingungen ergibt

12

yy 103

32

3

012 6 12 6 6 V04 6 2EI Mq012 6 6 V12

. 4 0vsym

⎡ ⎤⎡ ⎤− ⎡ ⎤ ⎡ ⎤⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥− ⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥= +⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥−⎢ ⎥⎢ ⎥ ⎢ ⎥ ⎢ ⎥

′ −⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎣ ⎦ ⎣ ⎦⎣ ⎦ ⎣ ⎦

l l

l l l ll

ll

l l

Aus der letzten Gleichung folgt 3

03

yy

qv

48EI′ = −

l, womit dann aus den drei ersten Gleichungen

die Auflagerreaktionsgrößen ermittelt werden können:

yy 0 0 01 3 02

2 2 2yy 20 0 0

1 3 02

yy 0 0 03 3 02

6EI q q q 5V v q2 8 2 8

2EI q q q 1M v q12 24 12 8

6EI q q q 3V v q2 8 2 8

′= − = − − = −

′= − = − = −

′= − − = − = −

l l ll

l

l l ll

l

l l ll

l

Die Verformungen des Innenknotens sind

F. U. Mathiak 11-27

4 42 0 0 3

yy yy2

3

04 4 0v 1 1q q v1

012 120 0 2384EI 8 384EI 8v 2 2v

iz

⎡ ⎤⎡ ⎤ ⎢ ⎥−⎢ ⎥⎡ ⎤ ′⎡ ⎤ ⎡ ⎤ ⎡ ⎤⎢ ⎥⎢ ⎥= = + = −⎢ ⎥ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥⎢ ⎥′ ⎢ ⎥⎢ ⎥ ⎣ ⎦ ⎣ ⎦ ⎣ ⎦− − −⎣ ⎦ ⎢ ⎥⎢ ⎥ ′⎢ ⎥⎣ ⎦ ⎣ ⎦

l lll l

l l

und damit 4

02

yy

qv192EI

=l sowie

30

2 3yy

q1v v4 192EI

′ ′= − =l in Übereinstimmung mit des Ergebnis-

sen des Beispiels 10.2.

11-28 Balkenelement

11.4 Der elastisch gebettete Balken

Abb. 11-18 Balken auf nachgiebiger Unterlage Wir betrachten einen Balken, der vollständig auf einer elastischen Unterlage liegt1. Der Bal-ken sei durch Streckenlasten und Einzellasten in z-Richtung belastet (Abb. 11-18). Nach Wink-ler wird angenommen, dass der Bodendruck qB(x) proportional zur lokalen Eindringtiefe w(x) ist:

Bq (x) c w(x)= Gl. 11-60

Die Konstante c heißt Bettungszahl.

1 22 2

Masse Nc , Einheit kgm sLänge (Zeit) m

− −⎡ ⎤ = =⎣ ⎦ ⋅

Material Bettungszahl c [MN/m2]

Sand, locker, rund 10...15 Sand, mitteldicht, rund 50...100 Sand, dicht, eckig 150...250 Geschiebemergel, fest 30...100 Lehm, halbfest 20...50 Torf 0,4...1

Tabelle 11-1 Rechenwerte von Bettungszahlen c einiger ausgewählter Böden für Vorentwürfe2 Die Differentialgleichung der elastischen Linie ermitteln wir aus der linearisierten Momenten-Krümmungsbeziehung yy yEI w (x) M (x)′′ = − sowie unter Beachtung der lokalen Gleichge-

wichtsbedingungen y BM (x) q(x) q (x)′′ ⎡ ⎤= − −⎣ ⎦ zu yyEI w (x) q(x) c w(x)″

⎡ ⎤′′ = −⎣ ⎦ .

1 Solche Systeme treten z.B. im Bauwesen bei Flachgründungen auf. 2 Nach Empfehlungen des Arbeitsausschusses Ufereinfassungen - EAU, 8. Aufl. 1990

F. U. Mathiak 11-29

Für einen Balken mit konstanter Biegesteifigkeit EIyy folgt daraus IV

yyEI w (x) c w(x) q(x)+ = und nach Division durch die Biegesteifigkeit EIyy

IV

yy yy

c q(x)w (x) w(x)EI EI

+ = Gl. 11-61

Substituieren wir noch mit

yy yy4 44 EI 4EI

L Lc c

= → = Gl. 11-62

die als charakteristische Länge bezeichnete Konstante L, so erhalten wir abschließend die inhomogene lineare Differentialgleichung 4. Ordnung mit konstanten Koeffizienten für die Biegung des elastisch gebetteten Balkens:

IV4

yy

4 q(x)w (x) w(x)EIL

+ = Gl. 11-63

Die allgemeine Lösung der Gl. 11-63 setzt sich additiv zusammen aus der Lösung wh der zuge-

hörigen homogenen Differentialgleichung IVh h4

4w (x) w (x) 0L

+ = mit der Lösung ( x / Lζ = )

( ) ( )h 1 2 3 4w (x) e C cos C sin e C cos C sinζ −ζ= ζ + ζ + ζ + ζ Gl. 11-64

und einer beliebigen partikulären Lösung wp der inhomogenen Differentialgleichung

IVp p4

yy

4 q(x)w (x) w (x)EIL

+ = Gl. 11-65

so dass für die Gesamtlösung h pw(x) w (x) w (x)= + gilt. Ein für die Praxis wichtiger Belas-

tungsfall ist q(x) = q0 = konst., für den durch Einsetzen in Gl. 11-65 die partikuläre Lösung

p 0w (x) q / c= nachgewiesen werden kann.

11.4.1 Die Steifigkeitsmatrix für die elastischer Bettung Um hier zu einer Lösung im Sinne der Methode der Finiten Elemente zu kommen, denken wir uns den Boden ersetzt durch beliebig dicht gelagerte linear elastische Federn mit der Feder-

11-30 Balkenelement

steifigkeit c. Auf die infinitesimale Balkenlänge dx entfällt die Formänderungsenergie

2F

1dW c w (x)dx2

= , und die gesamte Formänderungsenergie ist dann

2F F

x 0

1W dW c w (x)dx2 =

= =∫ ∫l

Gl. 11-66

und deren Variation

Fx 0

W c w(x) w(x) dx=

δ = δ∫l

Gl. 11-67

Unter Beachtung von w (e) (e)N z= und w (e) (e)N zδ = δ folgt aus Gl. 11-67

(e) (e)T (e) (e)F F

x 0

W c N N dx(e)T (e) (e)T (e) (e)cz z z k z

=

δ = δ = δ∫l

In der obigen Beziehung ist

(e )

(e) (e) (e)F

0

c dx(e)T (e)ck N N= ∫

l

Gl. 11-68

der Anteil der Elementsteifigkeitsmatrix der elastischen Bettung, zu der der Anteil aus der Balkenkrümmung hinzuzurechnen ist

ˆ (e) (e) (e)ck = k + k Gl. 11-69

Bei einem kubischen Verschiebungsansatz nach Gl. 11-5 erhalten wir nach kurzer Rechnung

(e) (e)

(e)2 (e) (e)2(e) (e)(e) F

(e)

(e)2

156 22 54 134 13 3c

156 22420. 4

ck

sym

⎡ ⎤−⎢ ⎥

−⎢ ⎥= ⎢ ⎥−⎢ ⎥⎢ ⎥⎣ ⎦

l l

l l ll

l

l

Gl. 11-70

und entsprechend für den quintischen Verschiebungsansatz nach Gl. 11-34

F. U. Mathiak 11-31

(e) (e) (e)

(e)2 (e) (e)2 (e) (e)2

(e)(e) (e)(e) F

(e)2 (e) (e)2

(e)

(e)2

2092 114 880 160 262 298 88 12 29 3

5632 0 880 88c128 160 1213860

2092 114. 8

ck

sym

⎡ ⎤− −⎢ ⎥

− −⎢ ⎥⎢ ⎥−⎢ ⎥=

−⎢ ⎥⎢ ⎥−⎢ ⎥⎢ ⎥⎣ ⎦

l l l

l l l l l

ll

l l l

l

l

Gl. 11-71

Hinweis: Numerische Untersuchungen zeigen, dass aufgrund der hochgradig nichtlinearen Lösungsfunktionen des elastisch gebetteten Balkens (s.h. Gl. 11-64) der kubische Verschie-bungsansatz zu unbefriedigenden Ergebnissen führt. Hier ist es besser, das höherwertige 3-Knoten-Element zu verwenden.

11.4.2 Beispiel 10.4 Für den links eingespannten und rechts frei gelagerten Balken auf elastischer Unterlage (Bet-tungszahl c) mit konstanter Belastung q = q0 sind sämtliche Zustandsgrößen zu ermitteln.

Geg.: l = 10m, E = 30000 MN/m2, Iyy = 5,10 10-3 m4, q0 = 1MN/m, c = 50 MN/m2

Abb. 11-19 Elastisch gebetteter Träger mit konstanter Belastung q0

Mit den Abkürzungen yy44EI

L 1.87c

= = ; 0.535L

λ = =l

l ; xξ =

l; ζ = ξλ lautet die voll-

ständige analytische Lösung: ( ) ( ) 01 2 3 4

qw(x) e C cos C sin e C cos C sin

cζ −ζ= ζ + ζ + ζ + ζ + .

Die Konstanten

( )( )

2

1 2 2

q sin cos cos e coshC

2k cos cosh

−λλ λ − λ − λ=

λ + λ

l

( )( )

2

2 2 2


2k cos cosh

−λλ λ + λ − λ= −

λ + λ

l

( )( )

2

3 2 2


2k cos cosh

λλ λ + λ + λ= −

λ + λ

l

( )( )

2

4 2 2


2k cos cosh

λλ λ − λ + λ= −

λ + λ

l

ergeben sich aus den Randbedingungen am linken Rand: w(x 0) 0; w (x 0) 0′= = = = und

den Kraftrandbedingungen am rechten Rand: y zM (x ) 0; Q (x ) 0= = = =l l .

11-32 Balkenelement

Für die FE-Lösung verwenden wir das 3-Knotenelement. Wir diskretisieren den Balken durch ein Element der Länge 10m=l . Damit ist die Elementsteifigkeitsmatrix identisch mit der Systemsteifigkeitsmatrix. Die Steifigkeitsmatrix für den elastisch gebetteten Balken setzt sich

zusammen aus der Steifigkeitsmatrix (e)k nach Gl. 11-43 und der Elementsteifigkeitsmatrix für

die elastische Bettung (e)ck nach Gl. 11-71. Einsetzen der Systemwerte liefert

97.728 90.872 16.079 26.211 2.860 .117173.991 -7.422 96.596 -.117 5.789

234.509 0 16.079 7.422ˆ1021.303 -26.211 96.596

97.728 -90.872. 173.991

(e) (e) (e)ck k + k

sym

⎡ ⎤⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥= =⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎣ ⎦

Nach Gl. 11-49 erhalten wir den Elementlastvektor, der bei unserem Beispiel auch dem Sys-temlastvektor entspricht:

(e)

(e)0

(e)

14 2.3331.667

32 5.333q0 060

14 2.333-1.667

(e)p

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥

= =⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥−⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦

l

l

l

Damit ist folgendes Gleichungssystem zu lösen

1

1

2

2

3

3

97.728 90.872 16.079 26.211 2.860 .117 v 2.333173.991 -7.422 96.596 -.117 5.789 v 1.667

234.509 0 16.079 7.422 v 5.3331021.303 -26.211 96.596 v 0

97.728 -90.872 v 2.. 173.991 vsym

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥′⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ =⎢ ⎥

′⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥

′⎢ ⎥ ⎢ ⎥⎣ ⎦⎣ ⎦

333-1.667

⎡ ⎤⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎢ ⎥⎣ ⎦

Wir stellen zunächst fest, dass sich die Randbedingungen 1 1v 0; v 0′= = am linken Rand prob-

lemlos in das obige Gleichungssystem einbauen lassen. Für den rechten freien Rand, an dem Biegemoment und Querkraft verschwinden müssten, stehen im Vektor der Systemfreiheits-grade nur die Verschiebung 3v und die Tangentenneigung 3v′ zur Randwertvorgabe zur Ver-

fügung. Damit ist aber die Anpassung der Lösung an die Kraftrandbedingungen des rechten

F. U. Mathiak 11-33

Randes nicht möglich. Wir können deshalb nur die homogenen Verschiebungsrandbedingun-gen am linken Rand berücksichtigen und erhalten

2

2

3

3

97.728 90.872 16.079 26.211 2.860 .117 0 2.333173.991 -7.422 96.596 -.117 5.789 0 1.667

234.509 0 16.079 7.422 v 5.3331021.303 -26.211 96.596 v 0

97.728 -90.872 v 2.333. 173.991 vsym

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ =⎢ ⎥

′⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥

′⎢ ⎥ ⎢ ⎥⎣ ⎦⎣ ⎦

1

1

VM000

-1.667 0

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥

+⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦

Das obige Gleichungssystem besitzt die Lösungen

1

1

2

2

3

3

v 0v 0v .02135089v .00053696v .02036687v .00015066

⎡ ⎤ ⎡ ⎤⎢ ⎥ ⎢ ⎥′⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥

=⎢ ⎥ ⎢ ⎥′⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥′ −⎢ ⎥ ⎢ ⎥⎣ ⎦⎣ ⎦

1

1

V 1.9177 MNM 1.7765 MNm0 00 00 00 0

⎡ ⎤ ⎡ ⎤−⎢ ⎥ ⎢ ⎥−⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥

=⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎢ ⎥ ⎢ ⎥⎣ ⎦ ⎣ ⎦

2 3 4 5w(x) .4427323269 1.211408646 1.197961273 .4089180835= ξ − ξ + ξ − ξ [m]

2 3 4w (x) .08854646538 -.3634225938 .4791845092 .2044590418′ = ξ ξ + ξ − ξ [/] 2 3w (x) .008854646538 -.07268451876 .1437553528 -.08178361672 ′′ = ξ + ξ ξ [1/m]

2w (x) .007268451876 .02875107056 .02453508502′′′ = − + ξ − ξ [1/m2]

Schnittlasten: 2 3

y yyM (x) EI w (x) 1.3548 11.1207 21.9946 12.5129′′= − = − + ξ − ξ + ξ [MNm] 2

z yyQ (x) EI w (x) 1.1121 4.3989 3.7539′′′= − = − ξ + ξ [MN]

Abb. 11-20 Resultierende Kraftgrößen

11-34 Balkenelement

Ein Vergleich der Querkraft zQ (x 0) 1.1121MN= = mit der berechneten Auflagerkraft

1V 1.9177 MN= − zeigt, dass die Auflagerkraft nicht als Grenzwert der Querkraft Qz(x) am

linken Rand ermittelt werden kann. Vielmehr steht V1 im energetischen Mittel mit der Resul-

tierenden 1

0 B0

R [q q (x)]d 1.9177 MNξ=

= − ξ =∫l aus der äußeren Belastung q0 und dem entge-

gengesetzt wirkenden Bodendruck Bq (x) c w(x)= im globalen Gleichgewicht: 1V R 0+ = .

Der Angriffspunkt der Resultierenden liegt bei

12

0 B0

[q q ]da 0.926m

Rξ=

ξ − ξ

= =∫l

. Momen-

tengleichgewicht bezogen auf den linken Rand erfordert: 1M aR 0+ = , womit auch

1M 1.7765 MNm= − nachgewiesen ist.

Die folgenden Darstellungen zeigen den Vergleich der analytischen Lösung mit derjenigen nach der FE-Methode. Am rechten Rand, an dem keine Randwerte vorgegeben wurden, weicht die FE-Lösung deutlich von der analytischen Lösung ab. Weiterhin ist auch hier wie-der zu beobachten, dass die durch Differentiationen an der Verschiebungsfunktion gewonne-nen Zustandsgrößen (Biegemomente und Querkräfte) doch erheblich von der theoretisch ex-akten Lösung abweichen.

Abb. 11-21 Durchbiegung des elastisch gebetteten Balkens [m]

F. U. Mathiak 11-35

Abb. 11-22 Biegemoment des elastisch gebetteten Balkens [MNm]

Abb. 11-23 Querkraft des elastisch gebetteten Balkens [MN]

7 Der Arbeits- und Energiebegriff in der Elastostatikuser · 2 q q 6 4 d f 3 f A 1 0 2 2 a l l l l...

Documents

Transcript of 7 Der Arbeits- und Energiebegriff in der Elastostatikuser · 2 q q 6 4 d f 3 f A 1 0 2 2 a l l l l...