BZM Heft 1 09

BEITRÄGE
ZUR
MATERIALTHEORIE
Herausgeber: Univ.-Prof. Dr.-Ing. habil. Alexander Lion

___________________________________________________________________
Dipl.-Ing. Nico Diercks
Validierung und ANSYS-

Implementierung des Schlimmer-
Mahnken Materialmodells zur
Festigkeitsberechnung von
Klebeverbindungen
___________________________________________________________________
Heft 01/09 ISSN 1862-5703

Validierung und ANSYS-Implementierung des
Schlimmer-Mahnken Materialmodells zur
Klebeverbindungen
von
Dipl.-Ing. Nico Diercks
Universität der Bundeswehr München

Fakultät für Luft- und Raumfahrttechnik
Institut für Mechanik (LRT 4)
Werner-Heisenberg-Weg 39
85579 Neubiberg
Deutschland
Heft 01/09 der Reihe
Beiträge zur Materialtheorie
ISSN 1862-5703
Herausgeber:
Univ.-Prof. Dr.-Ing. habil. Alexander Lion

Fakultät für Luft- und Raumfahrttechnik
Institut für Mechanik (LRT 4)
Werner-Heisenberg-Weg 39
85579 Neubiberg
Deutschland
___________________________________________________________________
Dipl.-Ing. Nico Diercks:

Validierung und ANSYS-Implementierung
des Schlimmer-Mahnken Materialmodells zur
Klebeverbindungen
Beiträge zur Materialtheorie 01/09

Neubiberg 2009
ISSN 1862-5703
Postanschrift des Verfassers:
Nico Diercks
Hansaring 111
50670 Köln
Dieses Werk ist urheberrechtlich geschützt. Jede Verwertung außerhalb der Grenzen
des Urheberrechtsgesetzes ist unzulässig und strafbar. Insbesondere gilt dies für die
Übersetzung, den Nachdruck sowie die Speicherung auf Mikrofilm, mit vergleichba-
ren Verfahren oder in Datenverarbeitungsanlagen.
Schlimmer-Mahnken Materialmodells
zur Festigkeitsberechnung von Klebeverbindungen
Inhaltsverzeichnis
Inhaltsverzeichnis
Abbildungsverzeichnis IV
1. Einleitung 1
1.1. Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.2. Ziel der Arbeit . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.3. Aufbau der Arbeit . . . . . . . . . . . . . . . . . . . . . . . . . 2
2. Grundlagen 3
2.1. Plastizität . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
2.1.1. Allgemeines . . . . . . . . . . . . . . . . . . . . . . . . . 3
2.1.2. eindimensionale Plastizitätstheorie . . . . . . . . . . . . 3
2.1.2.1. Einführung . . . . . . . . . . . . . . . . . . . . 3
2.1.2.2. Beispiel: „Prandtl-Reuss-Plastizität", eindimen-
sional . . . . . . . . . . . . . . . . . . . . . . . 6
2.1.2.3. Elasto-Plastizität mit linearer Verfestigung . . . 7
2.1.2.4. Zusammenhang zwischen den Moduln . . . . . 8
2.1.3. Verallgemeinerung auf 3 Dimensionen am Beispiel der
Prandtl-Reuss-Plastizität . . . . . . . . . . . . . . . . . . 10
2.1.3.1. Vorbemerkungen . . . . . . . . . . . . . . . . . 10
2.1.3.2. Zerlegung in deviatorische und volumetrische
Anteile . . . . . . . . . . . . . . . . . . . . . . 11
2.1.3.3. Fließfunktion (nach von-Mises) . . . . . . . . . 12
2.1.3.4. Fließbedingung . . . . . . . . . . . . . . . . . . 12
2.1.3.5. Fließregel . . . . . . . . . . . . . . . . . . . . . 12
2.1.3.6. Konsistenzbedingung . . . . . . . . . . . . . . . 13
2.1.3.7. Belastungsbedingung . . . . . . . . . . . . . . . 14
2.1.3.8. Materialgleichung . . . . . . . . . . . . . . . . . 14
2.1.3.9. Plastizität mit Verfestigung . . . . . . . . . . . 15
© Nico Diercks I
Inhaltsverzeichnis
2.2. Grundlagen der Finite-Elemente-Methode . . . . . . . . . . . . 17

2.2.1. Allgemeines . . . . . . . . . . . . . . . . . . . . . . . . . 17
2.2.2. Die schwache Form der Bewegungsdifferentialgleichung . 18
2.2.2.1. Das Prinzip von d’Alembert in der Lagrange’schen
Fassung . . . . . . . . . . . . . . . . . . . . . . 18
2.2.2.2. Beispiel: Längselastischer Stab . . . . . . . . . 19
2.2.3. Nichtlineare FEM . . . . . . . . . . . . . . . . . . . . . . 23
2.2.3.1. Formen der Nichtlinearität . . . . . . . . . . . . 23
2.2.3.2. Das Newton-Raphson Verfahren . . . . . . . . . 24
2.3. Das Schlimmer-Mahnken Modell . . . . . . . . . . . . . . . . . . 27
2.3.1. Allgemeines . . . . . . . . . . . . . . . . . . . . . . . . . 27
2.3.2. Kurzvorstellung und Besonderheiten des Modells . . . . 27
2.3.2.1. Notationen . . . . . . . . . . . . . . . . . . . . 28
2.3.3. Konstitutive Gleichungen . . . . . . . . . . . . . . . . . . 28
2.3.3.1. Zerlegung der Verzerrungen . . . . . . . . . . . 28
2.3.3.2. Fließfunktion und plastisches Potential . . . . . 29
2.3.3.3. Evolutionsgleichungen . . . . . . . . . . . . . . 30
2.3.3.4. Thermodynamische Konsistenz . . . . . . . . . 31
3. Implementierung 33
3.1. Numerische Vorbereitung . . . . . . . . . . . . . . . . . . . . . . 33
3.1.1. Allgemeines . . . . . . . . . . . . . . . . . . . . . . . . . 33
3.1.1.1. Strukturierung des Lösungsweges . . . . . . . . 33
3.1.1.2. Aufgabe der Subroutine . . . . . . . . . . . . . 33
3.1.2. numerische Formulierung des Modells . . . . . . . . . . . 34
3.1.2.1. Evolutionsgleichungen . . . . . . . . . . . . . . 34
3.1.2.2. Ermittlung der Trial-Spannungen . . . . . . . . 35
3.1.2.3. Berechnung des plastischen Multiplikators . . . 36
3.1.2.4. Spannungsupdate und Bestimmung der plasti-
schen Dehnungen . . . . . . . . . . . . . . . . . 39
3.1.2.5. Ermittlung der Tangentensteifigkeit . . . . . . . 40
3.2. ANSYS-spezifische Implementierung . . . . . . . . . . . . . . . 42
3.2.1. Transformation auf ANSYS-Größen . . . . . . . . . . . . 42
3.2.2. Die Routine usermat.3d . . . . . . . . . . . . . . . . . . 44
3.2.3. Verwendung der Subroutine durch den Anwender . . . . 45
© Nico Diercks II
Inhaltsverzeichnis
3.2.3.1. Aufruf der Routine . . . . . . . . . . . . . . . . 45

3.2.3.2. Eingabe der Materialparameter . . . . . . . . . 46
4. Modellvalidierung 48
4.1. Experimente mit einfacher Beanspruchung . . . . . . . . . . . . 49
4.1.1. Zugversuch . . . . . . . . . . . . . . . . . . . . . . . . . 49
4.1.2. Kompressionsversuch . . . . . . . . . . . . . . . . . . . . 50
4.1.3. Zugscherversuch . . . . . . . . . . . . . . . . . . . . . . . 51
4.2. Kombinierte Zug-Schubbelastung . . . . . . . . . . . . . . . . . 54
4.2.1. Versuchsaufbau . . . . . . . . . . . . . . . . . . . . . . . 54
4.2.2. Versuchsdurchführung . . . . . . . . . . . . . . . . . . . 55
4.2.3. Auswerung . . . . . . . . . . . . . . . . . . . . . . . . . 55
5. Fazit und kritische Bewertung 58
Literaturverzeichnis 59
Eidesstattliche Erklärung 61
A. Mathematische Grundlagen i
A.1. Tensoralgebra . . . . . . . . . . . . . . . . . . . . . . . . . . . . i
A.1.1. Einige Konventionen . . . . . . . . . . . . . . . . . . . . i
A.1.2. Tensorbegriff . . . . . . . . . . . . . . . . . . . . . . . . ii
A.1.3. Produkte von Tensoren . . . . . . . . . . . . . . . . . . . iii
A.1.4. Tensoren höherer Stufe . . . . . . . . . . . . . . . . . . . iv
A.2. Tensoranalysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . vi
A.2.1. Einführung des Funktionsbegriffes . . . . . . . . . . . . . vi
A.2.2. Tensorielle Ableitungen von Funktionen . . . . . . . . . . vi
B. Quellcode viii
B.1. Verzeichnis der verwendeten Variablen . . . . . . . . . . . . . . viii
B.1.1. Lokale Variablen . . . . . . . . . . . . . . . . . . . . . . viii
B.1.2. Durch ANSYS definierte Variablen . . . . . . . . . . . . x
B.1.3. ANSYS-Steifigkeitsmatrix . . . . . . . . . . . . . . . . . xii
B.2. Quellcode der usermat.3d . . . . . . . . . . . . . . . . . . . . . xiii
B.3. ANSYS-Macro zur Validierungssimulation . . . . . . . . . . . . xxv
© Nico Diercks III

Abbildungsverzeichnis
Abbildungsverzeichnis
2.1. σ(ε) eines Zugversuches . . . . . . . . . . . . . . . . . . . . . . 4

2.2. σ(ε) nach Prandtl-Reuss . . . . . . . . . . . . . . . . . . . . . . 6
2.3. σ(ε) mit Verfestigung . . . . . . . . . . . . . . . . . . . . . . . . 7
2.4. infinitesimales Stück der σ(ε) -Kurve . . . . . . . . . . . . . . . 9
2.5. 3-D Fließfläche des Prandtl-Reuss-Modells . . . . . . . . . . . . 13
2.6. Beispiel für FE-Vernetzung . . . . . . . . . . . . . . . . . . . . . 17
2.7. Stabelement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.8. Newton-Raphson Verfahren . . . . . . . . . . . . . . . . . . . . 25
4.1. Zugprobe nach EN ISO 527-2 . . . . . . . . . . . . . . . . . . . 49

4.2. Last-Verformungskurve des Zugversuches . . . . . . . . . . . . . 50
4.3. plastische Dehnung der Zugprobe . . . . . . . . . . . . . . . . . 50
4.4. Druckprobe . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
4.5. Last-Verformungskurve des Kompressionsversuches . . . . . . . 51
4.6. Zugscherprobe . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
4.7. Last-Verformungskurve des Zugscherversuches . . . . . . . . . . 52
4.8. v. Mises-Spannung der Zugscherprobe . . . . . . . . . . . . . . . 53
4.9. Experimenteller Aufbau zur kombinierten Belastung . . . . . . . 54
4.10. Spannungsverläufe . . . . . . . . . . . . . . . . . . . . . . . . . 56
4.11. axiale, plastische Dehnung . . . . . . . . . . . . . . . . . . . . . 57
© Nico Diercks IV
1. Einleitung
1. Einleitung
1.1. Motivation
Bei modernen Fertigungsprozessen, sowohl in der Automobilindustrie als auch

in der Luftfahrttechnik, werden immer öfter Klebeverbindungen zum Fügen
zweier Werkstücke verwendet. Da in diesen Bereichen höchste Anforderungen
an die Zuverlässigkeit der Strukturbauteile gestellt werden, werden genaue Fe-
stigkeitsrechnungen benötigt, wobei üblicherweise die Finite-Element-Methode
(FEM) angewandt wird.
Um diese Zuverlässigkeit zu erzielen, sind bei den Klebeverbindungen, zum

einen aufgrund der komplizierten Verbindungen und zum anderen aufgrund
fehlender numerischer Verfahren, noch viele Experimente zum Festigkeitsnach-
weis erforderlich. Diese Experimente sind meist sehr aufwändig und damit teu-
er, worin der Bedarf nach zuverlässigen numerischen Berechnungsmethoden
begründet liegt.
1.2. Ziel der Arbeit
Im Rahmen dieser Arbeit wird ein Materialmodell der Elastoplastizität in die

kommerzielle FEM-Software ANSYS implementiert, um es anhand von Simu-
lationsrechnungen durch Vergleich mit experimentellen Ergebnisse zu validie-
ren. Dieses Modell wurde von Schlimmer und Mahnken für die Behandlung
von Klebstoffen erarbeitet und im März 2005 veröffentlicht.
Das Hauptziel der Arbeit ist es nach Validierung des Modells, ein Werkzeug
für die numerische Festigkeitsberechnung von verklebten Bauteilen bereit zu
stellen. Damit kann in Zukunft die Zuverlässigkeit von Strukturberechnungen
© Nico Diercks 1
1. Einleitung
erhöht werden und gleichzeitig Kosten gespart werden, die teure Experimente
entstehen würden.
1.3. Aufbau der Arbeit
Der Hauptteil dieser Arbeit gliedert sich in 3 Abschnitte:

Im ersten dieser Abschnitte werden die Grundlagen erläutert und der Leser
an die Thematik herangeführt. Dabei sind es drei wesentliche Säulen, auf die
diese Arbeit aufbaut. Dieses sind
1. Plastizitätstheorie
2. Die Finite-Elemente-Methode
3. Das Schlimmer-Mahnken Materialmodell
Dabei wird eine kurze Einführung in die ersten beiden Gebiete geboten, um
dem Leser aufzuzeigen, wie diese Probleme allgemein bearbeitet werden, um
im Anschluss daran das schließlich zu implementierende Modell vorzustellen.
Der zweite Teil befasst sich mit der Implementierung des Modells an sich.
Dies ist, neben der Validierung durch Experimente, der Hauptteil der Arbeit.
Hierbei werden die Gleichungen ausführlich behandelt und in eine FE-taugliche
Formulierung gebracht. Anschließend wird auf die Besonderheiten der ANSYS-
Implementierung eingegangen.
Im dritten Abschnitt folgt die Validierung des Modells. Dazu werden verschie-
dene Experimente mithilfe der zuvor getätigten FE-Implementierung in AN-
SYS simuliert und anschließend die Ergebnisse mit den Messdaten der Expe-
rimente verglichen, um so Aussagen über die Verlässlichkeit der Simulation
treffen zu können.
Im Anschluss an den Hauptteil folgt die Bewertung der gewonnenen Erkennt-
nisse und ein Ausblick, welche Arbeiten in Zukunft angestellt werden sollten,
um daran anzuknüpfen.
© Nico Diercks 2
2. Grundlagen
2. Grundlagen
2.1. Plastizität
2.1.1. Allgemeines
Plastizität allgemein beschreibt einen Effekt fester Stoffe, bei dem unter Bela-
stung Verformungen auftreten, die nach Rücknahme der Belastungen dauerhaft
bleiben, also irreversibel sind. Üblicherweise treten neben diesen irreversiblen,
plastischen Dehnungen auch reversible, elastische Deformationen auf. In die-
sem Fall spricht man von Elasto-Plastizität. Ein wichtiges Merkmal der Elasto-
Plastizität besteht darin, dass keine Zeitabhängigkeit des Materialverhaltens
auftritt, womit das Verhalten unabhängig von der Belastungsgeschwindigkeit
ist. Sind zeitabhängige Effekte zu beobachten (z.B. beim Kriechen von Metallen
oder bei Kunststoffen oberhalb der Glastemperatur), spricht man von Visko-
Elastizität oder Visko-Plastizität. In der vorliegenden Arbeit sollen Klebstoffe
für Strukturbauteile modelliert werden. Da deren Einsatz üblicherweise un-
terhalb der Glastemperatur stattfindet, ist es nicht notwendig zeitabhängige
Effekte zu berücksichtigen, weshalb im weiteren auf die Elasto-Plastizität ein-
gegangen wird. Weitere Ausführungen zu dem Thema sind unter [5] zu fin-
den.
2.1.2. eindimensionale Plastizitätstheorie
2.1.2.1. Einführung
Um die Effekte und die Herangehensweise an selbige zu verdeutlichen, wird

zunächst ein eindimensionaler Spannungszustand betrachtet. Ein Beispiel hier-
für wäre ein einfacher Zugversuch, bei dem die Spannungs-Dehnungskurve be-
© Nico Diercks 3
2. Grundlagen
B
σ
A
σF0
C
εpl εel ε
ε
Abbildung 2.1.: σ(ε) eines Zugversuches
trachtet wird. In Abbildung 2.1 ist eine solche Spannungs-Dehnungskurve skiz-

ziert. Im Bereich zwischen Ursprung und dem Punkt A besteht ein linearer Zu-
sammenhang zwischen Dehnung und der zugehörigen Spannung (Hooke’sches
Gesetz). Bei Belastungen σ ≤ σF0 und anschließender Entlastung wird wieder
der Ursprung erreicht; es bleiben keine Deformationen zurück. Dies ist der Be-
reich der reinen Elastizität.
Wird der Stab jedoch weiter belastet (σ > σF0 ) bis zum Punkt B und anschlie-
ßend entlastet, so folgt die Entlastung einer Geraden (parallel zur elastischen
Geraden bis Punkt A) zum Punkt C und es bleibt eine plastische Dehnung εpl
zurück, wobei die additive Zerlegung der Dehnung in die Anteile εel und εpl
nur mit der Einschränkung auf kleine Deformationen möglich ist.
Eine erneute Belastung folgt der Graden von Punkt C zurück zu Punkt B, bis
es zu weiterer plastischer Dehnung entlang der vorgegebenen Strecke kommt.
Zur mathematischen Formulierung wird nun der Punkt B betrachtet und eine
Fallunterscheidung vorgenommen.
© Nico Diercks 4
2. Grundlagen
1. Entlastung: dσ < 0, hier gilt das elastische Materialgesetz (Hooke’sches

Gesetz) σ = Eε bzw.:
dσ
dε = dεel = , dεpl = 0 (2.1)
E
2. Belastungssteigerung: dσ > 0, wobei weitere Deformationen auftreten

dε = dεel + dεpl .
Zur mathematischen Formulierung der Fallunterscheidung wird eine Fließfunk-

tion formuliert:
F (σ, ε) = σ − σF (ε)
wobei σ die aktuelle Spannung und σF (ε) die Fließkurve gem. Abbildung 2.1
bezeichnet. Ist F (σ, ε) < 0 verhält sich das Material rein elastisch und die
Dehnungen werden anhand von Glg. 2.1 bestimmt. Plastische Dehnungen ent-
stehen für
F (σ, ε) = 0 . (2.2)
Diese Bedingung wird allgemein Fließbedingung genannt. Ist diese und die
Belastungsbedingung dσ > 0 erfüllt, so gilt
dε = dεel + dεpl .
Nach Bilden des totalen Differentials der Fließbedingung (Glg. 2.2) folgt die
Konsistenzbedingung:
∂F ∂F
dF = dσ + dε = 0
∂σ ∂ε
Das Materialgesetz der Hooke’schen Elastizität (Glg. 2.1) muss wegen der pla-
stischen Dehnung um einen Term erweitert werden zu
dσ dεpl
dε = dεel + dεpl = + dσ .
E dσ
dσ(ε )
An dieser Stelle wird der plastische Tangentenmodul Tp = dεplpl eingeführt,
der die lokale Steigung, abhängig von der plastischen Dehnung, beschreibt und
es folgt das Materialgesetz

1 1
dε = + dσ . (2.3)
E Tp
© Nico Diercks 5
2. Grundlagen
2.1.2.2. Beispiel: „Prandtl-Reuss-Plastizität", eindimensional
Bei diesem sehr einfachen Modell der Plastizität nach Prandtl-Reuss verhält
sich das Material zunächst linear elastisch bis zum Erreichen der Fließspan-
nung σF und anschließend ideal plastisch, d.h. es ist keine weitere Spannungs-
steigerung möglich und das Material beginnt uneingeschränkt zu fließen. Die
σ
σF0
σF0
Abbildung 2.2.: σ(ε) nach Prandtl-Reuss
zugehörige Fließfunktion lautet
1
F (σ, ε) = F (σ) = (σ 2 − σF2 )
2
Für |σ| < σF ist F (σ) < 0; das Material reagiert also elastisch
dσ
dε = dεel = , dεpl = 0
E
Wird σ weiter gesteigert bis |σ| = σF , so folgt F (σ) = 0. Auch wenn in

diesem Modell eine weitere Steigerung der Spannung dσ > 0 nicht möglich
ist, so treten dennoch plastische Dehnungen auf, dεpl 6= 0. Eine anschließende
Entlastung führt entlang der elastischen Geraden dεel = E1 dσ
© Nico Diercks 6
2. Grundlagen
2.1.2.3. Elasto-Plastizität mit linearer Verfestigung
Nach der relativ einfachen „Prandtl-Reuss-Plastizitätfolgt nun die Betrachtung

eines Plastizitätsmodells mit linearer Verfestigung, wie es in Abbildung 2.3 dar-
gestellt ist. Dieses skizzierte Materialverhalten kann durch die Fließfunktion
σ − σF
σF0
ε
εpl εel
Abbildung 2.3.: σ(ε) mit Verfestigung
1
(σ − Bεpl )2 − σF2 ,

F (σ, εpl ) = B = const.
2
beschrieben werden. Für das Auftreten von plastischen Deformation liefert die
Fließbedingung
F (σ, εpl ) = 0 ⇒ (σ − Bεpl )2 = σF2
σ = σF + Bεpl .
Nach Bildung des Differentials und Umstellen nach dεpl
dσ
dεpl =
B
mithilfe des Materialgesetzes der Elastizität folgt
dσ dσ
dε = dεel + dεpl = + .
E B
© Nico Diercks 7
2. Grundlagen
Nach Integration erhält man

1 1
ε= + σ+C
E B
Die Bestimmung der Integrationskonstante C erfolgt durch

σF 1 1 σF
ε(σF ) = = + σ+C ⇒ C=
E E B B
σ σ − σF
⇒ ε= +
E B
E
⇒ σ= (σF + εB) (2.4)
E+B
oder mit Verwendung des Tangentenmoduls T (Steigung der Kurve oberhalb
von σF ):
σ − σF σF
T = ⇒ σ = σF + T ε − (2.5)
ε − σEF E
2.1.2.4. Zusammenhang zwischen den Moduln
In den vorherigen Kapiteln wurde das Materialverhalten mithilfe folgender

Größen beschrieben:
1. dem Elastizitätsmodul
dσ(εel )
E=
dεel
2. dem Tangentenmodul
dσ(ε)
T =
dε
3. dem plastischen Tangentenmodul
dσ(εpl )
Tp =
dεpl
Weiterhin ist der Lösung von Kapitel 2.1.2.3 noch die (relativ unanschauliche)
Konstante B enthalten. Diesem gilt es nun auf den Grund zu gehen und einen
Zusammenhang zwischen diesen Größen herzustellen. Dazu wird in Abb. 2.4
ein infinitesimales Stück der σ(ε) -Kurve betrachtet.
© Nico Diercks 8
2. Grundlagen
dσ
dεpl dεel
dε
Abbildung 2.4.: infinitesimales Stück der σ(ε) -Kurve
Es gilt:
dσ = Tp dεpl = Tp (dε − dεel )

dσ dσ
mit dεel = ⇒ dσ = Tp dε −
E E

Tp
1+ dσ = Tp dε
E
dσ Tp
= =T.
dε 1 + TEp
Somit ist der Zusammenhang gefunden:
1 1 1
= + (2.6)
T E Tp
Zur Bestimmung der Konstanten B wird Glg. 2.4 erneut betrachtet:
E E 1
σ= (σF + εB) = σF + 1 1 ε (2.7)
E+B E+B E
+ B
© Nico Diercks 9
2. Grundlagen
Weiterhin lautet Glg. 2.5

σF
σ = σF + T ε −
E
T
= 1− σF + T ε
E
!
1 1
= 1− 1 σF + T ε
E
+ T1p E

Tp
= 1− σF + T ε
E + Tp
E 1
= σF + 1 ε. (2.8)
E + Tp E
+ T1p
Aus den Glgn. 2.4 und 2.8 folgt
E 1 E 1
σF + 1 1 ε = σF + 1 1 ε,
E+B E
+ B
E + Tp E
+ Tp
wobei man durch Koeffizientenvergleich sofort erkennt
B = Tp
Damit lautet die Beziehung zwischen Dehnungen und Spannungen für ein Ma-
terial mit linearer Verfestigung
σF E
σ = σF + T ε − = σF + εTp (2.9)
E E + Tp
2.1.3. Verallgemeinerung auf 3 Dimensionen am Beispiel

der Prandtl-Reuss-Plastizität
2.1.3.1. Vorbemerkungen
Bei der Formulierung der Gleichungen für die Modellierung der plastischen Ef-
fekte unter dreidimensionalen Spannungs- und Verzerrungszuständen ist das
Vorgehen ganz analog zum eindimensionalen Fall, nur dass Spannungen und
© Nico Diercks 10
2. Grundlagen
Verzerrungen nicht nur aus einer Komponente, sondern aus Tensoren mit je-
weils 9 Komponenten bestehen.
   
σ11 σ12 σ13 ε11 ε12 ε13
σ =  σ21 σ22 σ23  ε =  ε21 ε22 ε23 
   
σ31 σ32 σ33 ε31 ε32 ε33
Wegen der Symmetrie der Tensoren
σ = σT ε = εT
reduziert sich die Zahl der Komponenten auf insgesamt jeweils 6 verschiedene
Komponenten. Davon abgesehen werden entsprechende Gleichungen wie im
vorherigen Kapitel verwendet:
1. Fließfunktion
2. Fließbedingung
3. Fließregel
4. Konsistenzbedingung
5. Belastungsbedingung
6. bei Verfestigung: Evolutionsgleichungen für die Spannungen
Dieses Beispiel ist in etwas anderer Form in [10] zu finden, wo auch weitere
Ausführungen enthalten sind.
2.1.3.2. Zerlegung in deviatorische und volumetrische Anteile
Zunächst werden Spannungs- und Verzerrungstensor in die deviatorischen und

die volumetrischen Anteile zerlegt. Die Deviatorkomponenten beschreiben all-
gemein die Gestaltänderungskomponenten, während die Volumenanteile den
hydrostatischen Druck, bzw. die Volumenänderung ausdrücken.
ε = ε dev + ε vol σ = σ dev + σ vol
© Nico Diercks 11
2. Grundlagen
Für einen beliebigen Tensor A gilt:
1
Avol = Sp(A)1 mit Sp(A) = 1 : A
3
2.1.3.3. Fließfunktion (nach von-Mises)
Die Fließfunktion wird ganz ähnlich wie im 1-D-Fall formuliert, nur mit dem
Unterschied, dass nun alle Komponenten des deviatorischen Spannungstensors
eingehen:
1
F (σ dev ) = σ dev : σ dev − k 2
2
k ist eine Werkstoffkonstante analog zur Fließspannung σF darstellt, die den
Radius der Fließfläche beschreibt. In diesem Fall wird angenommen, dass die
Plastizität nur im Zusammenhang mit Gestaltänderung steht.
2.1.3.4. Fließbedingung
Wie auch bei allen vorherigen Betrachtungen lautet die Fließbedingung
F (σ dev ) = 0 (2.10)
Somit beschreibt die Fließregel zusammen mit der Fließbedingung eine Zylin-
derfläche mit dem Radius k = √ σF
3
im Hauptspannungsraum, wie in Abb. 2.5
gezeigt. Befindet sich der aktuelle Spannungszustand innerhalb dieser Zylin-
derfläche, so sind die Deformationen rein elastisch, während auf der Fläche
plastische Deformationen auftreten. Ein Zustand außerhalb der Fläche kann
nicht erreicht werden, da nach Einsetzen des Fließens keine weitere Spannungs-
zunahme stattfindet.
2.1.3.5. Fließregel
Zur Bestimmung des Zuwachses an plastischen Dehnungen wird folgende Fließ-

regel formuliert:
∂F
εpl = λ̇
dε = λ̇σ dev (2.11)
∂σ
© Nico Diercks 12
2. Grundlagen
σ2
σ1
σ3
Raumdiagonale σ1 = σ2 = σ3
Abbildung 2.5.: 3-D Fließfläche des Prandtl-Reuss-Modells
wobei der plastische Multiplikator λ̇ eingeführt wird. Die Größe λ̇ ist keine
freie an Experimente anzupassende Materialfunktion, sondern ergibt sich aus
der Konsistenzbedingung.
2.1.3.6. Konsistenzbedingung
Wie auch im 1-D-Fall lautet die Konsistenzbedingung
dF = 0
∂F ∂F
dF = : dσ + dk
∂σ ∂k
wegen k = const. ⇒ dk = 0
∂F
⇒ dF = : dσ dev = σ dev : dσ dev = 0
∂σ dev
Im allgemeinen Fall kann nun der plastische Multiplikator λ̇ bestimmt werden.

Bei diesem Modell genügen dafür schon die Fließbedingung (Glg. 2.10) und
die Fließregel (Glg. 2.11):
1
F (σ dev ) = σ dev : σ dev − k 2 = 0
2
© Nico Diercks 13
2. Grundlagen
1
k 2 λ̇ = λ̇σ dev : σ dev
2
mit Glg. 2.11
1
k 2 λ̇ = ε pl : σ dev
2
Da sich die plastischen Dehnungen nur auf den Deviatoranteil auswirken, gilt
schließlich für den plastischen Multiplikator
1
λ̇ = ε dev : σ dev . (2.12)
2k 2
2.1.3.7. Belastungsbedingung
Die Belastungsbedingung aus dem 1-D-Fall dσ > 0 muss nun erweitert werden
zu
ε>0
B = σ dev : dε
Für das Auftreten plastischer Deformationen muss sowohl die Fließbedingung,

als auch die Belastungsbedingung erfüllt sein.
2.1.3.8. Materialgleichung
Aus diesem Satz an Gleichungen kann nun die Materialgleichung unter Be-
rücksichtigung der plastischen Verformungen formuliert werden:
ε = dε
dε εel + dε
εpl
Mit den elastischen Dehnungen aus dem Hooke’schen Gesetz

1+ν ν
εel =
dε σ+ Sp(σ)1
E 1+ν
und den plastischen Dehnungen aus der Fließregel (Glg. 2.11) folgt die Mate-
rialgleichung für die Verzerrungsraten:

1+ν ν 1
ε=
dε dσ + εdev ) σ dev
Sp(dσ)1 + 2 (σ dev : dε
E 1+ν 2k
© Nico Diercks 14
2. Grundlagen
2.1.3.9. Plastizität mit Verfestigung
Um das doch sehr einfache Modell der Prandtl-Reuss-Plastizität zu erweitern

besteht die Möglichkeit, isotrope oder kinematische Verfestigung zu berück-
sichtigen. Eine isotrope, also richtungsunabhängige, Verfestigung bedeutet,
dass die Fließspannung σF wie in Abschnitt 2.1.2.3 mit zunehmender pla-
stischer Verformung ansteigt. Dies entspricht im 3-D-Fall einer Vergrößerung
des Radius der Fließfläche. Kinematische Verfestigung hingegen ist mit einem
Verschieben der Fließfläche zu modellieren. Es treten also je nach Belastungs-
richtung andere Fließspannungen auf. Da das Schlimmer-Mahnken-Modell nur
isotrope Verfestigung modelliert, wird auch nur auf diese näher eingegangen.
Im Vergleich zum vorherigen Beispiel bleibt die Fließfunktion in ihrer Form
also gleich, jedoch ist k keine Konstante mehr, sondern eine Materialfunktion,
die mit Hilfe von Evolutionsgleichungen bestimmt werden muss.
1 k2
F (σ dev , k) = σ dev : σ dev −
2 3
Der Faktor 1/3 wurde dabei aus rechentechnischen Gründen eingeführt. Die
Fließbedingung und die Fließregel bleiben unverändert:
∂F
F (σ dev ) = 0 und εpl = λ̇
dε = λ̇σ dev (2.13)
∂σ
Zur Beschreibung der Verfestigung wird die Evolutionsgleichung für die Fließ-
spannung
dk = k 0 (sp )dsp
als Funktion der plastischen Bogenlänge sp eingeführt, die sich nach

r
2
dsp = εpl : dε
dε εpl (2.14)
3
entwickelt. Da durch die Verfestigung auch bei plastischer Dehnung ein Span-
nungszuwachs möglich ist, lautet die Belastungsbedingung
B = σ dev : dσ > 0
© Nico Diercks 15
2. Grundlagen
Die Konsistenzbedingung mit dF = 0 bleibt zwar gleich, allerdings wird die

Formulierung aufgrund von dk 6= 0 um einen Term umfangreicher:
dF = 0
∂F ∂F
= : dσ + dk
∂σ ∂k
2
= σ dev : dσ − kk 0 (sp )dsp = 0
3
mit Glg. 2.14
r
2 0 2
= σ dev : dσ − kk (sp ) dεεpl : dε
εpl
3 3
mit Glg. 2.13 (Fließregel)
r
2 0 2
= σ dev : dσ − kk (sp )λ̇ σ dev : σ dev
3 3
mit der Fließbedingung
s
2
2 0 2
0 = σ dev : dσ − kk (sp )λ̇ k2 (2.15)
3 3
σ dev : dσ
⇒ λ̇ = 2 2 0
3
k k (sp ) 23
Nach Berechnung des plastischen Multiplikators ist es nun unter Anwendung
der Fließregel (Glg. 2.13) ein Leichtes, den Zuwachs an plastischen Dehnungen
zu berechnen. Die elastischen Dehnungen werden wieder mit dem Hooke’schen
Gesetz bestimmt und es folgt:
ε = dε
dε εel + dε
εpl

1+ν ν σ dev : dσ
ε=
dε dσ + Sp(dσ)1 + 2 σ dev
E 1+ν 2
k k 0 (s )
p
3
© Nico Diercks 16
2. Grundlagen
2.2. Grundlagen der Finite-Elemente-Methode
2.2.1. Allgemeines
Die Finite-Element-Methode ist ein Näherungsverfahren mit dessen Hilfe die

Berechnung komplizierter Strukturen mit hinreichender Genauigkeit möglich
ist und bei denen die klassischen Methoden mit dem Versuch analytische Lö-
sungen zu finden, nicht mehr anwendbar sind. Ursprünglich wurde dieses Ver-
fahren in den fünfziger Jahren von Flugzeugbau-Ingenieuren entwickelt und
machte anschließend, beflügelt durch die Entwicklungen auf dem Gebiet der
Rechner- und Computertechnik, eine derart rasante Entwicklung durch, dass
es mittlerweile nicht nur im Bereich der Mechanik, sondern darüber hinaus in
der Thermodynamik, der Strömungsmechanik und auf vielen anderen Gebie-
ten als praktisch wichtigstes Näherungsverfahren zählt.
reale Geometrie
Finite Elemente
Abbildung 2.6.: Beispiel für FE-Vernetzung
Das Grundprinzip dieses Verfahrens ist relativ simpel: Ein Kontinuum, beste-
hend aus unendlich vielen (infiniten) materiellen Punkten, wird in eine end-
liche (finite) Anzahl von Knoten und Elementen zerlegt, die für die folgende
Rechnung das Kontinuum repräsentieren (siehe Abb. 2.6). Dabei wird dieses
dann (zunächst) nicht mehr durch ein kontinuierliches Verschiebungsfeld, son-
dern lediglich durch die Verschiebung der Knotenpunkte beschrieben. Für die
Formulierung der Bewegungsdifferentialgleichungen wird dann das kontinuier-
liche Verschiebungsfeld durch Interpolation der Knotenpunktsverschiebungen
© Nico Diercks 17
2. Grundlagen
gewonnen. Die Steifigkeit des Kontinuums wird als Steifigkeitsmatrix (bezogen

auf die Knoten) formuliert und alle äußeren Lasten als Einzellasten an den
Knotenpunkten modelliert.
Somit werden die Differentialgleichungen durch Gleichungssysteme (mit der

Steifigkeitsmatrix, dem Vektor der Knotenpunktsverschiebungen und dem Vek-
tor der äußeren Lasten) ersetzt, die dann relativ einfach numerisch behandelt
werden.
Grundsätzlich unterscheidet man bei der FEM zwischen zwei Methoden
• Kraftgrößenmethode: Hier werden die Spannungen im Element appro-

ximiert und somit die Gleichgewichtsbedingungen exakt erfüllt. Da die
Verschiebungen dann aus den Spannungen berechnet werden, sind die
Kompatibilitätsbedingungen nur im Mittel erfüllt
• Verschiebungsmethode: Bei dieser Methode werden die Verschiebungen

im Element approximiert, womit die Kompatibilitätsbedingungen exakt
und die Gleichgewichtsbedingungen nur im Mittel erfüllt werden können.
Da jedoch die Randbedingungen eines Problems in der Regel in Form von

Verschiebungsgrößen vorliegen, ist die Verwendung der Kraftgrößenmethode
äußerst schwierig und daher auch ungewöhnlich. Üblicherweise wählt man aus
diesem Grund die Verschiebungsgrößenmethode, auf die im Folgenden auch
weiter eingegangen wird.
Im Rahmen dieser Arbeit wird das Prinzip der FEM kurz erläutert und einige
Methoden vorgestellt; dabei geht es nur um das prinzipielle Verständnis der
Vorgehensweise. Für nähere Studien wird im Bereich der linearen Theorie die
Lektüre von [4] und für nichtlineare Probleme von [6] und [1] empfohlen.
2.2.2. Die schwache Form der

Bewegungsdifferentialgleichung
2.2.2.1. Das Prinzip von d’Alembert in der Lagrange’schen Fassung
Als Ausgangspunkt für die Formulierung der Bewegungsdifferentialgleichungen

wählt man gern das Prinzip von d’Alembert in der Lagrange’schen Fassung.
© Nico Diercks 18
2. Grundlagen
Hierbei wird auf die betrachtete Struktur ein virtuelles Verschiebungsfeld δu

aufgeprägt und anschließend die dadurch verrichteten Arbeiten bilanziert. Da
diese Arbeitsbilanz nach der Diskretisierung nicht vollständig, sondern nur
im Mittel erfüllt sein muss, spricht man auch von der schwachen Form der
Bewegungsdifferentialgleichung.
δWσ = δWa + δWT (2.16)
wobei δWσ die Arbeit der inneren Spannungen, δWa die Arbeit der äußeren
Lasten, und δWT die Arbeit der Trägheitskräfte bezeichnet. Anschließend gilt
es, Gleichung 2.16 als Funktion der Größen der FEM zu formulieren. Im Fol-
genden soll nun am Beispiel eines längselastischen Stabes exemplarisch die
FE-Formulierung des kontinuumsmechanischen Problems gezeigt werden. Bei
höherwertigen Elementen ist diese Formulierung natürlich etwas umfangrei-
cher; das Prinzip bleibt jedoch immer das Gleiche: am Ende der Formulierung
steht immer anstelle von Differentialgleichungen ein Gleichungssystem für die
Verschiebungen der Knotenfreiheitsgrade.
2.2.2.2. Beispiel: Längselastischer Stab
u1 x qx u2
F1x L F2x
Abbildung 2.7.: Stabelement
Der Stab ist das einfachste Element (Abb. 2.7), da er sich nur in Längsrich-
tung verformt und somit nur über die zwei Freiheitsgrade u1 und u2 , die die
Verschiebung der Knoten 1 und 2 beschreiben, verfügt. Die Belastung erfolgt
dabei durch die Knotenkräfte F1x , F2x und die Streckenlast qx . Zunächst wird
Glg. 2.16 kontinuumsmechanisch ausgewertert. Dabei ist
ZL ZL
δWσ = N δε dx , δWT = − ρAü δu dx
0 0
© Nico Diercks 19
2. Grundlagen
ZL
δWa = F1x δu1 + F2x δu2 + qx δu dx
0
Anschließend folgt die Diskretisierung der Verschiebungen und der Kräfte zu

! !
~û = u1 (t) F1x (t)
r~c =
u2 (t) F2x (t)
Um das kontinuierliche Verschiebungsfeld zu erhalten, wird eine Verschiebungs-

interpolationsfunktion eingeführt. In diesem Fall wird der einfachste Ansatz,
eine lineare Interpolation, gewählt.
!
x x u1 (t)
u(x, t) = (1 − )u1 (t) + u2 (t) = [H1 (x) H2 (x)]
L L u2 (t)
mit den Formfunktionen H1 (x) = 1 − Lx und H2 = Lx und der Matrix der Form-
funktionen oder Verschiebungsinterpolationsmatrix H(x) = [H1 (x) H2 (x)]. So-
mit lautet das Verschiebungsfeld
u(x, t) = H(x) · ~û
Die Dehnung im Stab lautet:
∂ ∂
ε= (u(x, t)) = (H(x)~û) = H0 (x)~û
∂x ∂x
und für die Normalkraft folgt
N = EAε = EA H0 (x) · ~û
Das Beschleunigungsfeld lautet
∂ 2u ∂2
ü = = (H(x)~û) = H(x)~û..
∂t2 ∂t2
© Nico Diercks 20
2. Grundlagen
Anschließend folgt für die virtuelle Arbeit der inneren Spannungen:
ZL ZL
δWσ = ~
N δε dx = δû T
H0T EA H0 dx ~û
0 0
bzw.
~ T k ~û
δWσ = δû (2.17)
mit der Elementsteifigkeitsmatrix
ZL
k= H0T EA H0 dx .
0
Die virtuelle Arbeit der Trägheitskräfte lautet
ZL ZL
δWT = − ~T
ρAü δu dx = −δû HT ρAHdx ~û..
0 0
bzw.
~ T m ~û..
δWT = −δû (2.18)
mit der Elementmassenmatrix
ZL ZL
m= HT ρAHdx = HT IHdx
0 0
wobei die Trägheitsmatrix I eingeführt wird. In diesem Fall ist dies lediglich
ein Skalar ρA.
Schließlich wird die virtuelle Arbeit der äußeren Lasten bestimmt:
ZL ZL
δWa = F1x δu1 + F2x δu2 + ~ T · r~c +
qx δu dx = δû qx δu dx
0 0
© Nico Diercks 21
2. Grundlagen
Die Streckenlast qx muss nun im Intervall 0 ≤ x ≤ L approximiert werden,

wobei die gleiche Formfunktion wie bei den Knotenpunktverschiebungen An-
wendung findet:
!
x x q1
qx = (1 − )q1 + q2 = [H1 (x) H2 (x)] = H~q
L L q2
ZL
~ T (~
δWa = δû rc + HT H dx ~q)
0
bzw.
~ T (~
δWa = δû ~ T · ~r
rc + r~s ) = δû (2.19)
ZL
mit r~s = HT H dx ~q = L~q
0
Zusammengefasst lautet Glg. 2.16 mit den Glgn. 2.17, 2.18, 2.19 also
δWσ = δWa + δWT mit

~ T k ~û
δWσ = δû
~ T m ~û..
δWT = −δû
~ T · ~r
δWa = δû
Die Auswertung dessen liefert
~ T (k ~û + m ~û.. − ~r) = 0

δû
~ T 6= 0 und beliebig folgt

Wegen δû
k ~û + m ~û.. = ~r . (2.20)
Diese Gleichung gilt für ein Element. Zur Lösung des Gesamtproblems werden
diese Elementgleichungen in ein globales Gleichungssystem zusammengeführt,
dessen Form sehr ähnlich zun Glg. 2.20 ist:
~ ~ ~
K Û + M Û .. = R (2.21)
© Nico Diercks 22
2. Grundlagen
Ganz analog zu Glg. 2.20 sind:

K und M die Gesamtsteifigkeits- bzw. Gesamtmassenmatrix sowie
~
Û und R~ die Spaltenmatrix der globalen Knotenfreiheitsgrade
bzw. der globalen Knotenkräfte
Eine genauere Betrachtung der Glg. 2.21 offenbart, dass die Matrizen K und
M lediglich Konstanten enthalten, die von den gewählten Formfunktionen und
von den Materialparametern abhängen. Bei der Betrachtung statischer Proble-
~
me (der Beschleunigungsterm fällt weg wegen Û .. = 0) wurde somit ein Glei-
chungssystem für die Berechnung der Knotenpunktsverschiebungen gewonnen,
das anschließend mit numerischen Methoden gelöst werden kann. Hierbei ist
die statische Bestimmtheit des Systems auch ohne Bedeutung; die Randbedin-
gungen werden einfach an den entsprechenden Knotenpunktsfreiheitsgraden
eingearbeitet.
2.2.3. Nichtlineare FEM
2.2.3.1. Formen der Nichtlinearität
In den vorherigen Ausführungen gelang es, das kontinuumsmechanische Pro-

blem der Statik auf ein lineares Gleichungssystem der Form
~ ~
KÛ = R
zu transformieren. Dies setzt jedoch eine lineare Theorie mit einer konstanten
Matrix K voraus, was in vielen Fällen leider nicht gegeben ist. Die Gründe da-
für können sehr unterschiedlich sein; man unterscheidet in der Regel zwischen
drei Formen der Nichtlinearität:
1. Physikalische Nichtlinearität
beschreibt Nichtlinearitäten im Materialverhalten. Beispiele dafür sind
die bereits vorgestellte Plastizität mit nichtlinearer Spannungs- Deh-
nungsbeziehung oder auch Viskoelastizität. Oft ist mit physikalischer
Nichtlinearität aufgrund großer Verzerrungen geometrische Nichtlinea-
rität verbunden.
© Nico Diercks 23
2. Grundlagen
2. Geometrische Nichtlinearität
führt aufgrund großer Deformationen zu nichtlinearem Verhalten. Als
Beispiel sei ein Biegebalken mit großen Deformationen genannt, bei dem
bei gleich bleibender Belastung die Deformation selbst durch die Ände-
rung des Hebelarmes zu einer Veränderung des Biegemomentes führt.
Auch alle Formen von Durchschlags- oder Verzweigungsproblemen gehö-
ren zu dieser Form der Nichtlinearität.
3. Nichtlinearität infolge wechselnder Randbedingungen

tritt meistens bei Kontaktproblemen auf, wobei sich also die Randbe-
dingungen infolge der Deformation verändern. Genau genommen ist dies
keine eigentliche Nichtlinearität, sondern lediglich abschnittsweise linea-
res Verhalten. Oft kommt aufgrund der großen Deformationen jedoch
zusätzlich noch eine andere Form von Nichtlinearität hinzu.
Da im Rahmen dieser Arbeit ein Materialmodell der Plastizität behandelt

wird, ist die physikalische Nichtlinearität von Bedeutung. Im Folgenden wird
eine Möglichkeit zur Berechnung dessen vorgestellt.
2.2.3.2. Das Newton-Raphson Verfahren
Ausgangspunkt für dieses Verfahren ist die selbe Gleichung wie in der linearen
Theorie, nur dass die Steifigkeitsmatrix eine Funktion des aktuellen Verschie-
bungszustandes ist:
K(~û)~û = ~r (2.22)
Zur Veranschaulichung wird der eindimensionale Fall betrachtet, wobei das

Verschiebungsfeld auf eine Größe u, der Vektor der äußeren Lasten auf eine
Kraft F und die Steifigkeitsmatrix auf einen Eintrag K(u) reduziert wird.
Da aufgrund der Nichtlinearität das Superpositionsprinzip nicht mehr anwend-
bar ist, muss die Belastung in mehreren Schritten inkrementell aufgebracht
werden. Ein Lastschritt wird dabei ∆F und die zugehörige Verschiebung mit
un bezeichnet. Beim Newton-Raphson Verfahren, wie es in Abb. 2.8 darge-
stellt ist, wird zur Berechnung der Verschiebung un+1 zunächst die Steifigkeit
des Vorschrittes Kt (un ) verwendet. In einem ersten Prädiktorschritt wird die
© Nico Diercks 24
2. Grundlagen
F
Fn+1
g0 g1
dg1
Kt = du
∆F
dg0
Kt = du
Fn
∆F
δu0 δu1
un = u0n+1 u1n+1 u2n+1 ukn+1 u
Abbildung 2.8.: Newton-Raphson Verfahren
Verschiebung u1n+1 ermittelt, die bei abschnittsweise linearem Verhalten ent-

stehen würde. Anschließend wird mithilfe eines Residuums g, das aus Glg.
2.22 gewonnen wird, die Verschiebung überprüft und mit mehreren Korrektor-
schritten korrigiert. Dabei wird stets eine neue Tangentensteifigkeit Kt (ukn+1 )
ermittelt.
Da Glg. 2.22 weiterhin ihre Gültigkeit behält, ergibt sich durch die Verwen-
dung der Verschiebungen ukn+1 und uk+1
n+1 eine „out-of-balance-Force", die als
Residuum fungiert:
g = K(ukn+1 )uk+1
n+1 − F
Durch Taylorreihenentwicklung mit Abbruch nach dem ersten Glied entsteht

die Linearisierung:
dg0
g ≈ g0 + δu = 0
du
−1
dg0
⇒ δu0 = − g0 (u0 )
du
© Nico Diercks 25
2. Grundlagen
woraus der nächste Verschiebungswert berechnet werden kann zu
u1n+1 = u0n+1 + δu0
bzw.
uk+1 k
n+1 = un+1 + δuk
Dabei entspricht der Ableitung du dg

genau die Tangentensteifigkeit Kt , die aus
den Materialgesetzen bestimmt werden muss. Nach jeder Iteration wird das
Residuum g ausgewertet. Unterschreitet dies einen gewissen Schwellenwert, so
wird der nächste Lastschritt ∆F aufgebracht und die Iteration über k startet
auf ein Neues.
Das hier vorgestellt Verfahren wird auch Full Newton-Raphson Verfahren ge-
nannt. In Abgrenzung dazu verwendet man bei bestimmten Problemstellungen
auch gern das Modified Newton-Raphson Verfahren. Hierbei wird die Tangen-
tensteifigkeit nicht bei jeder Iteration neu berechnet, sondern nur einmal pro
Lastschritt. Dadurch sind zwar mehr Iteraionen notwendig, jedoch entfällt da-
für der Rechenaufwand für die Neuberechnung der Tangente.
Der Vorteil des Newton-Raphson Verfahrens liegt in seiner schnellen Kon-
vergenz. Allerdings entstehen sehr schnell Schwierigkeiten, sobald die Last-
Verformungskurve nicht monoton steigt. Hier müssen dann andere Verfahren,
wie Bogenlängenverfahren, angewendet werden.
Die Verallgemeinerung auf mehrere Dimensionen ist zwar graphisch recht un-
anschaulich, aber dennoch relativ einfach. Anstelle der skalaren Verschiebung
~
u tritt der Vektor der Knotenverschiebungen Û , anstelle der skalaren Kraft
~
wird der Vektor der Knotenkräfte R̂ verwendet und die Tangentensteifigkeit
Kt muss erweitert werden zur Steifigkeitsmatrix. Ansonsten ist das Vorgehen
ganz analog. Um die Steifigkeitsmatrix für nichtlineare Materialmodelle zu
berechnen, muss aus den Modellgleichungen der Steifigkeitstensor C 4. Stufe
ermittelt werden.
© Nico Diercks 26
2. Grundlagen
2.3. Das Schlimmer-Mahnken Modell
2.3.1. Allgemeines
An dieser Stelle wird das Materialmodell für Klebstoffe, wie es von Michael
Schlimmer und Rolf Mahnken [9] erarbeitet wurde, kurz vorgestellt. Insbe-
sondere die Formulierung der konstitutiven Gleichungen sind zum großen Teil
Auszüge der Veröffentlichung. Zum besseren Verständnis wurden einige Kom-
mentare und Zwischenrechnungen eingeführt.
2.3.2. Kurzvorstellung und Besonderheiten des Modells
Das vorliegende Modell beschreibt elasto-plastisches Verhalten mit isotroper

Verfestigung.
Zur Erarbeitung des Modells wurden Experimente mit dem Klebstoff “Be-
tamate 1496” durchgeführt, wobei die Belastung aus kombinierter Zug- und
Torsionsspannung gewählt wurde. Hierbei wurde deutlich, dass sich, je nach
Verhältnis von Zug- zur Torsionsbelastung, eine andere Fließgrenze einstellt.
Diese von den Autoren “strenght difference effect” genannte Beobachtung ist
in der Formulierung des Modells besonders berücksichtigt werden.
Eine Besonderheit von Klebstoffen in Abgrenzung zur Metallplastizität ist der

Effekt, dass die plastischen Verformungen nicht nur aus den deviatorischen
Spannungen resultieren (also eine reine Gestaltänderung), sondern dass auch
die erste Spannungsinvariante (also der hydrstatische Druck) in die Berech-
nung der plastischen Dehnungen eingeht. In diesem Punkt unterscheidet sich
das von Schlimmer und Mahnken entwickelte Modell ganz grundlegend von an-
deren Modellen der Plastizität. Darin wird auch die Notwendigkeit der Verwen-
dung speziell für Klebstoffe entwickelter Materialmodelle deutlich, um dieser
Besonderheit Rechnung zu tragen.
© Nico Diercks 27
2. Grundlagen
2.3.2.1. Notationen
Bei der Vorstellung dieses Modells wird die Notation von Schlimmer und Mahn-
ken verwendet, wonach eckige Klammern [· ] verwendet werden, um funktio-
nelle Zusammenhänge deutlich zu machen (f [x] bedeutet: f ist eine Funktion
von x). Dies ist geschieht in Abgrenzung zu runden Klammern (· ), mit denen
mathematische Ausdrücke gruppiert werden.
Die drei Invarianten eines beliebigen Tensors 2. Stufe werden wie folgt defi-
niert:
1 1
I1 [A] = 1 : A, I2 [A] = 1 : A2 , I3 [A] = 1 : A3 (2.23)
2 3
wobei 1 den Einheitstensor 2. Stufe mit der Eigenschaft 1· u = u für einen

beliebigen Vektor u darstellt. Wie in Kapitel 2.1.3.2 gezeigt, ist es sinnvoll,
Tensoren in ihren deviatorischen und einen volumetrischen Anteil zu zerlegen.
Der Tensor A wird dabei auf folgende Weise zerlegt:
A = Adev + Avol mit

1 1
Adev = Idev : A, Avol = I1 [A]1, Idev = I − 1 ⊗ 1 (2.24)
3 3
wobei I der Einheitstensor 4. Stufe ist, für den gilt: I : A = A.
2.3.3. Konstitutive Gleichungen
2.3.3.1. Zerlegung der Verzerrungen
Für eine Formulierung mit kleinen Deformationen (also geometrischer Linea-

rität) gilt wie in den vorherigen Betrachtungen:
ε = ε el + ε pl
ε el = C−1 : σ el (2.25)
C = 2GIdev + K 1⊗1
© Nico Diercks 28
2. Grundlagen
Aus Gleichungen 2.25 wird deutlich, dass sich der Verzerrungstensor additiv in
einen elastischen und einen plastischen Anteil zerlegen lässt. Der elastische An-
teil ist eine Funktion des Spannungstensors und des Steifigkeitstensors (Tensor
4. Stufe). Dieser wiederum bestimmt sich direkt aus Materialparametern; dem
Schubmodul G und dem Kompressionsmodul K (Hooke’sches Gesetz). Der pla-
stische Anteil der Verzerrungen hingegen wird aus den Evolutionsgleichungen
berechnet, die im Folgenden eingeführt werden.
2.3.3.2. Fließfunktion und plastisches Potential
Wie bereits in Abschnitt 2.1 erläutert, wird eine Fließbedingung gern als Funk-
tion des Spannungstensors in der Form φ[σ] = 0 formuliert. Für ein isotropes
Material reduziert sich diese Bedingung zu φ[I1 , I2 , I3 ] = 0, wobei die drei In-
varianten I1 , I2 , I3 des Cauchy Spannungstensors analog zu Gleichung 2.23 be-
rechnet werden. Alternativ kann man auch eine Fließbedingung φ[I1 , I20 , I30 ] = 0
formulieren, wobei I20 und I30 die Invarianten des deviatorischen Anteils des
Spannungstensors bezeichnen.
1 1
I1 = 1 : σ I20 = 1 : σ 2dev I30 = 1 : σ 3dev (2.26)
2 3
Wie auch in Schlimmer [11] gezeigt, kann mit einer derartigen Fließfunktion
sowohl der Einfluss des hydrostatischen Druckes als auch der deviatorischen
Spannungen berücksichtigt werden. Im vorliegenden Modell wird die dritte
Invariante I30 nicht berücksichtigt und die Fließfunktion wie folgt formuliert:
1
φ = I20 − ϕ
3
2
ϕ = Y − a1 Y0 I1 − a2 I12
(2.27)
Y = Y0 + R[ev ]
R[ev ] = q(1 − exp[−bev ]) + Hev
Hierbei wird die spannungsartige Variable R[ev ] eingeführt, die (abhängig von
der dehnungsartigen inneren Variable ev , die im nächsten Abschnitt erläutert
wird) die Verfestigung des Materials beschreibt. Weiterhin fällt auf, dass diese
Form der Fließfunktion einige bekannte Modelle der Plastizität einschließt:
© Nico Diercks 29
2. Grundlagen
• a1 = a2 = 0 führt auf eine zylindrische Fließfläche, wie sie von Mises

formulierte (siehe Kap. 2.1.3.3)
a2
• Y = Y0 , µ = a21 6= 0, a2 = − 41 liefert den Kegel 3I20 = Y0 − µI1 von
p
Drucker und Prager
• a1 = 0, a2 6= 0 entspricht dem Ellipsoid von Green
Um eine nichtassoziierte Fließregel zu erhalten wird zusätzlich zur Fließfunkti-

on in Gleichung 2.27 ein plastisches Potential mit der gleichen mathematischen
Struktur eingeführt.
1
φ∗ = I20 − ϕ∗
3 (2.28)
ϕ = Y − a∗1 Y0 I1 − a∗2 I12
∗ 2
wobei a∗1 und a∗2 als zusätzliche Materialparameter eingeführt werden. Somit
wird das Materialverhalten in diesem Modell mit insgesamt 10 Parametern be-
schrieben. Neben denen der Elastizität G und K sind dies: Y0 , q, b, H, a1 , a2 , a∗1 , a∗2 .
2.3.3.3. Evolutionsgleichungen
Der Zuwachs der plastischen Dehnungen wird wie folgt aus dem plastischen
Potential berechnet:
∂φ∗
ε˙pl = λ̇ = λ̇σ dev + λ̇t1
∂σ (2.29)
1
t = (a∗1 Y0 + 2a∗2 I1 )
3
dabei wird der plastische Multiplikator λ̇ aus dem Be- und Entlastungsbedin-
gungen bestimmt
λ̇ ≥ 0, λ̇φ = 0 φ ≤ 0 (2.30)
Zur Formulierung einer Evolutionsgleichung der inneren Variable ev werden 2

Ansätze, die auf der dissipierten Leistung basieren, formuliert:
1. e˙v Y = σ : ε˙pl 2. e˙v Y0 = σ : ε˙pl
© Nico Diercks 30
2. Grundlagen
Da sich später herausstellen wird, dass eine Implementierung des zweiten Fal-
les wesentlich wirtschaftlichere Berechnungen ermöglicht, wird auch nur dieser
implementiert. Aus diesem Grund wird davon abgesehen, den ersten Fall wei-
ter zu betrachten und alle Ausführungen im weiteren Rahmen dieser Arbeit
befassen sich mit dem zweiten Fall, woraus folgt:
1 λ̇
e˙v = σ : ε˙pl = (tI1 + 2I20 ) (2.31)
Y0 Y0
2.3.3.4. Thermodynamische Konsistenz
Ohne an dieser Stelle auf die Details der Thermodynamik einzugehen, muss
nach der Clausius-Duhem-Ungleichung folgende Forderung an das Modell ge-
stellt werden:
εpl − Re˙v ≥ 0
σ : ε̇
mit Gleichung 2.31:

1
εpl (1 − R
σ : ε̇ )≥0 (2.32)
Y0
Einschub: Nach Multiplikation von Gleichung 2.31 mit Y0
εpl = λ̇(tI1 + 2I20 )

σ : ε̇ (2.33)
mit t gemäß Gleichung 2.29 folgt
1
εpl
σ : ε̇ = λ̇( (a∗1 Y0 + 2a∗2 I1 )I1 + 2I20 )
3
1 2
εpl
σ : ε̇ = λ̇( a∗1 Y0 I1 + a∗2 I12 + 2I20 )
3 3
Da sowohl λ̇, als auch I20 und I12 nicht negativ werden können, kann durch die
Forderung a∗1 = 0 und a∗2 ≥ 0 sichergestellt werden, dass
εpl ≥ 0
σ : ε̇
Somit kann aus Gleichung 2.32 gefolgert werden:
1
(1 − R )≥0
Y0
© Nico Diercks 31
2. Grundlagen
Mit der Formulierung für R im Rahmen der Einführung der Fließfunktion

(Gleichungen 2.27) folgt
Y0 − q(1 − exp[−bev ]) + Hev ≥ 0
Hev ≤ Y0 − q(1 − exp[−bev ])
Da (1 − exp[−bev ] ≤ 1) muss ev ≤ Y0H−q gefordert werden. Zusammengefasst

lauten die Restriktionen an die Materialparameter und die innere Variable
also:
Y0 − q
Y0 ≥ 0 a∗1 = 0 a∗2 ≥ 0 H ≥ 0 q ≥ 0 ev ≤ (2.34)
H
Bemerkung: Ein Vergleich der Fließfunktion φ (Glg. 2.27) mit dem plastischen
Potential φ∗ (Glg. 2.28) zeigt, dass die Fließregel nur für a1 = a∗1 und a2 = a∗2
assoziiert ist. Um die beobachteten Effekte modellieren zu können, muss jedoch
a1 > 0 gefordert werden. Damit ist die Fließregel wegen a∗1 = 0 nichtassozia-
tiv.
© Nico Diercks 32
3. Implementierung
3. Implementierung
3.1. Numerische Vorbereitung
3.1.1. Allgemeines
3.1.1.1. Strukturierung des Lösungsweges
Eine FE-Berechnung mit einem eigenen Materialmodell läuft auf mehreren

Ebenen ab. In der obersten Ebene steht ein kommerzielles FE-Programm,
welches die Lösung des globalen Problems mit den in Kapitel 2.2 kurz vor-
gestellten Methoden berechnet. In diese Ebene greift der User nur durch die
gewohnten Befehle (Erstellen der Geometrie, Vernetzung derselbigen, Defini-
tion von Randbedingungen, Lasten und Materialparametern, Vorgabe von zu
verwendenden Verfahren etc.) ein; eine eigentliche Programmierung findet hier
nicht statt.
Diese Programmierung ist in der darunter liegenden Ebene nötig. Um ein
eigenes Materialmodell für diese Berechnung zu verwenden, wird durch das
Hauptprogramm auf Elementebene eine Subroutine aufgerufen, die dann die
ensprechenden Größen berechnet und an das Hauptprogramm zurück gibt. In
diese Subroutine kann jedes beliebige Materialverhalten programmiert werden,
da die Berechnung darin völlig unabhängig vom Hauptprogramm abläuft. Hier
wird dann im Folgenden das in Kapitel 2.3 vorgestellte Materialmodell nach
Schlimmer und Mahnken implementiert.
3.1.1.2. Aufgabe der Subroutine
Bei Aufruf der Subroutine werden vom Hauptprogramm für jeden globalen
Iterationsschritt lediglich die aktuelle Gesamtdehnung ε (n+1) , die plastischen
© Nico Diercks 33
3. Implementierung
(n)
Dehnungen der Vorschrittes ε pl und der letzte Wert für die Verfestigungsgröße
(n)
ev für das jeweilige Element übergeben. Die Subroutine arbeitet also nur auf
Element-Ebene, die Lösung des globalen Problems ist davon völlig losgelöst.
Nun gilt es, aus den übergebenen Größen die aktuellen plastischen Dehnungen
(n+1) (n+1)
ε pl und die Verfestigungsgröße ev zu berechnen und ans Hauptprogramm
zurück zu geben.
Für die Lösung des globalen Systems benötigt das Hauptprogramm zusätzlich
noch die Tangentensteifigkeit C = dσ ε des Elementes, welche ebenfalls durch
dε
die Subroutine berechnet werden muss. Die Spannungen σ, die im Element
entstehen, werden für die Lösung an sich zwar nicht benötigt, dürften für den
Anwender jedoch von großem Interesse sein, weshalb diese ebenfalls auf dieser
Ebene berechnet werden.
3.1.2. numerische Formulierung des Modells
3.1.2.1. Evolutionsgleichungen
Die in Abschnitt 3.1.1.2 geforderten Größen lauten nach Glgn. 2.29 und 2.31:
(n+1)
ε pl ε(n)
=ε εpl
pl + ∆ε
∆εεpl =λ̇σ (n+1)

dev + λ̇(n+1) t1 (3.1)
1 (n+1)
t(n+1) = (a∗1 Y0 + 2a∗2 I1 )
3
und
e(n+1)
v =e(n)
v + ∆ev
λ̇ (n+1) (n+1) (3.2)

0(n+1)
∆ev = (t I1 + 2I2 )
Y0
Schnell wird deutlich, dass zur Auswertung der Gleichungen die aktuellen
(n+1) 0(n+1)
Spannungsinvarianten I1 und I2 benötigt werden, was die Kenntnis
des aktuelle Spannungszustandes erforderlich macht. Um die Notation über-
sichtlich zu halten, wird im Folgenden auf den Index (n+1) verzichtet, d.h. alle
Größen ohne Index beziehen sich auf den aktuellen Iterationsschritt, während
die Größen des Vorschrittes mit (n) gekennzeichnet werden.
© Nico Diercks 34
3. Implementierung
3.1.2.2. Ermittlung der Trial-Spannungen
Da der Spannungszustand zu Beginn der Iteration nicht bekannt ist (und ohne
Kenntnis der plastischen Dehnungen auch nicht vollständig berechnet werden
kann), behilft man sich wie folgt:
• zunächst wird angenommen, dass die zusätzliche Verformung rein ela-

stisch ist, so dass die zugehörigen Spannungen mithilfe des Hooke’schen
Gesetzes berechnet werden können. Da dies im Allgemeinen nicht die
korrekten Spannungen liefert, sondern nur ein erster Versuch ist, nennt
man diese Spannungen Trial-Spannungen σ tr
• mithilfe der Trial-Spannungen und des Verfestigungswertes des Vorschrit-

(n)
tes ev werden die Invarianten berechnet, die Fließfunktion ausgewertet
und die Fließbedingung überprüft.
• ist die Fließbedingung nicht erfüllt (F bzw. φ < 0), so war die An-
nahme korrekt, alle Dehnungsinkremente sind elastisch und die Trial-
Spannungen entsprechen den tatsächlichen Spannungen
• ist die Fließbedingung jedoch erfüllt (F bzw. φ ≥ 0), so treten plastische

Deformationen auf. In einem iterativen Prozess wird nun der plastische
Multiplikator λ̇, die tatsächlichen Spannungen σ, die plastischen Deh-
nungensinkremente ∆ε εpl sowie die Verfestigungsgröße ev nach und nach
angenähert, bis der aktuelle Spannungszustand genau auf der Fließfläche
liegt (F bzw. φ = 0).
© Nico Diercks 35
3. Implementierung
Nach dem Hooke’schen Gesetz ist
(n)
σ tr = C : ε el = C : (ε
ε − ε pl )
(n) (n)
= (2GIdev + K1 ⊗ 1) : (ε
εdev − ε pl,dev + ε vol − ε pl,vol )
wegen Idev : Avol = 0
Idev : Adev = Adev
(1 ⊗ 1) : Adev = 0
(1 ⊗ 1) : Avol = Sp(Avol )1 = Sp(A)1
(n)
σ tr = 2G(εεdev − ε pl,dev ) + K Sp(ε ε − ε (n)
pl )
(n) 1
= 2G(εεdev − ε pl,dev ) + I1tr 1
3
(n)
mit tr
I1 = 3K Sp(ε ε − ε pl ) (3.3)
σ tr = σ tr tr
dev + σ vol
wobei
(n)
σ tr
dev = 2G(ε εdev − ε pl,dev ) (3.4)
1 tr
σ tr
vol = I 1 (3.5)
3 1
Aus Glgn. 3.4 und 3.5 erkennt man, dass die Trial-Spannungen direkt aus den
(n)
vom Hauptprogramm übergebenen Größen ε pl , ε (n+1) und aus Materialpara-
metern berechnet werden.
3.1.2.3. Berechnung des plastischen Multiplikators
Wenn sich nach Auswertung der Fließbedingung ergibt, dass plastische Deh-
nungen auftreten, so gilt es nun, den plastischen Multiplikator zu bestimmen.
In diesem Fall führt dies auf eine Iteration in der zugleich die aktuellen Span-
nungsinvarianten bestimmt werden können. Wäre dies nicht möglich, so müsste
ein „Radial-Return Algorithmus", wie er in [1] oder [6] beschrieben ist, zur An-
wendung gebracht werden, um die Spannungen zurück auf die Fließfläche zu
projizieren.
© Nico Diercks 36
3. Implementierung
Zur Bestimmung der Spannungen wird wieder das Hooke’sche Gesetz heran-
gezogen:
σ = C : ε el
ε − ε pl )
= C : (ε
(n)
ε − ε pl − ∆ε
= C : (ε εpl )
= σ tr − C : ∆ε
εpl
mit ∆ε
εpl aus Glg. 3.1
(n+1)
= σ tr − (2GIdev + K1 ⊗ 1) : (λ̇σ dev + λ̇(n+1) t1)
= σ tr − 2Gλ̇σ dev − λ̇K3t1 (3.6)
Anschließend wird die erste Invariante I1 gebildet:
I1 = 1 : σ
= 1 : (σ tr − 2Gλ̇σ dev − λ̇K3t1)
wegen 1 : Adev = 0
I1 = 1 : (σ tr − λ̇K3t1)
= I1tr − 3λ̇K3t
mit t aus Glg. 3.1
I1 = I1tr − 3λ̇K(a∗1 Y0 + 2a∗2 I1 )
I1tr − 3K λ̇a∗1 Y0
⇒ I1 = (3.7)
1 + 3K λ̇2a∗2
Eine Betrachtung des deviatorischen Anteils von Glg. 3.6 liefert:
σ dev = σ tr
dev − 2Gλ̇σ dev
σ tr
dev = (1 + 2Gλ̇)σ dev
σ tr
dev = δσ dev (3.8)
mit δ = 1 + 2Gλ̇
© Nico Diercks 37
3. Implementierung
Nach Multiplikation beider Seiten von Glg. 3.8 mit sich selbst folgt
δ 2 σ dev : σ dev = σ tr tr
dev : σ dev
δ 2 I20 = I20tr (3.9)
Mit dieser Kenntnis kann nun die innere Variable ev nach Glg. 3.2 entwickelt
werden zu
λ̇
ev = e(n)
v + (tI1 + 2I20 )
Y0
λ̇ 2
= e(n)
v + (tI1 + 2 I20tr ) . (3.10)
Y0 δ
Da bei Erfüllung der Fließbedingung φ ≥ 0 ist, der einzige physikalisch sinn-

volle Wert allerdings φ = 0 lauten muss, kann damit nun ein Residuum r[λ̇]
formuliert werden, in das alle zuvor hergeleiteten Gleichungen eingehen. Da-
mit kann der plastische Multiplikator λ̇ und damit auch die Veränderung der
Spannungen und der Fließfläche bestimmt werden, so dass schließlich φ = 0
gilt.
Zusammengefasst lauten die benötigten Gleichungen:
r[λ̇] = 3I20tr − ϕδ 2 mit

δ = 1 + 2Gλ̇
ϕ = Y 2 − a1 Y0 I1 − a2 I12
Y = Y0 + R
R = q(1 − e−bev ) + Hev
(3.11)
λ̇ 2
ev = e(n)
v + (tI1 + 2 I20tr )
Y0 δ
1 ∗
t = (a1 Y0 + 2a∗2 I1 )
3
I tr − 3K λ̇a∗1 Y0
I1 = 1
1 + 3K λ̇2a∗2
Daraus ergibt sich folgende Iteration für den plastischen Multiplikator:
−1
∂r
λ̇(k+1)
= λ̇ (k)
− r[λ̇(k) ] (3.12)
∂ λ̇
© Nico Diercks 38
3. Implementierung
wobei die Ableitung ∂r

∂ λ̇
aus den Glgn. 3.11 bestimmt wird:
∂r ∂ϕ
= −2δϕ2G − δ 2 mit (3.13)
∂ λ̇ ∂ λ̇
∂ϕ dY ∂∆ev ∂I1
=2Y − (a1 Y0 + 2a2 I1 )
∂ λ̇ dev ∂ λ̇ ∂ λ̇ (3.14)
∂∆e v ∂I1
=2Y (b e−bev + H) − (a1 Y0 + 2a2 I1 )
∂ λ̇ ∂ λ̇
∂I1 −3Ka∗1 Y0 I1tr − 3K λ̇a∗1 Y0

= ∗
− ∗ 2
3K2a∗2
∂ λ̇ 1 + 3K λ̇2a2 (1 + 3K λ̇2a2 )
a∗1 Y0 (1 + 3K λ̇2a∗2 ) + 2a∗2 I1tr − 3K λ̇a∗1 Y0 2a∗2
= − 3K
(1 + 3K λ̇2a∗2 )2
mit I1tr = I1 + 3K λ̇(a∗1 Y0 + 2a∗2 I1 )

∂I1 a∗ Y0 + 2a∗2 I1
= −3K 1 (3.15)
∂ λ̇ 1 + 3K λ̇2a∗2

∂∆ev 1 2 0tr
= tI1 + 2 I2 +
∂ λ̇ Y0 δ

λ̇ dt ∂I1 ∂ 2 0tr
+ I1 + t + I2
Y0 dI1 ∂ λ̇ ∂ λ̇ δ 2
(3.16)
1 2 0tr
= tI1 + 2 I2 +
Y0 δ

λ̇ 2 ∗ ∂I1 8G 0tr
+ a I1 + t − 3 I2
Y0 3 2 ∂ λ̇ δ
3.1.2.4. Spannungsupdate und Bestimmung der plastischen

Dehnungen
Zur Bestimmung der tatsächlichen Spannungen σ werden die bereits im vor-

herigen Kapitel betrachteten Spannungsanteile σ dev und σ vol verwendet. σ dev
folgt direkt aus Glg. 3.8, während σ vol mithilfe der ersten Spannungsinvariante
© Nico Diercks 39
3. Implementierung
I1 berechnet werden kann, die bereits während der Iteration für den plastischen
Multiplikator in Glg. 3.12 bestimmt wurde.
1 1
σ dev = σ tr σ vol = I1 1
δ dev 3
1 1
⇒ σ = σ tr
dev + I1 1 (3.17)
δ 3
Da nun sowohl der plastische Multiplikator als auch der vollständige Span-
nungszustand bekannt ist, lässt sich der Zuwachs an plastischen Dehnungen
nach Glg. 3.1 leicht bestimmen:
εpl = λ̇σ dev + λ̇t1

∆ε
3.1.2.5. Ermittlung der Tangentensteifigkeit
Zuletzt fehlt noch die Berechnung der Tangentensteifigkeit C, die das Haupt-
programm zur Lösung des globalen Systems benötigt. Zum Verständnis dieser
Herleitung ist die Kenntnis der Rechenregeln zur Tensoranalysis im Anhang
(A.2) notwendig.

dσ d 1 tr 1
C= = σ + I1 1
ε
dε ε δ dev 3
dε
1 ∂σ tr

tr ∂ −1 dev ∂ 1
=σ dev ⊗ + + I1 1
∂εε δ2 δ ∂εε ∂εε 3

tr −1 dδ ∂ λ̇ 2G ∂ 1 (n)
=σ dev ⊗ 2 + ε − Sp(ε ε)1 − ε pl,dev +
δ dλ̇ ∂ε ε δ ∂εε 3

d 1
+ I1 1
ε 3
dε
!
−2G ∂ λ̇ 1 1 1 ∂I1 ∂ λ̇ 1 ∂I1
=σ tr
dev ⊗ 2
+ 2G(I − 1 ⊗ 1) + 1 ⊗ +
δ ∂ε ε δ 3 ε 3 ∂ε
3 ∂ λ̇ ∂ε ε
∂I1 3K
mit = 1
∂εε 1 + 3K λ̇2a∗2

2G 1 ∂I1 2G tr ∂ λ̇ K
C= Idev + 1 − 2 σ dev ⊗ + 1⊗1 (3.18)
δ 3 ∂ λ̇ δ ∂εε 1 + 3K λ̇2a∗2
© Nico Diercks 40
3. Implementierung
Die fehlende Ableitung ∂∂ελ̇ε wird nach genauerer Betrachtung des Residuums aus
Glg. 3.11 gewonnen, wobei dieses als implizite Funktion wie folgt aufgefasst
wird:
dr ∂r ∂r ∂ λ̇
ε, ε [λ̇]] = 0 ⇒
r[ε = + =0
ε
dε ε ∂ λ̇ ∂ε
∂ε ε
−1
∂ λ̇ ∂r ∂r
⇒ =− (3.19)
∂εε ∂ λ̇ ∂εε
Die Ableitung ∂∂rλ̇ wurde bereits bei der Iteration des plastischen Multiplikators
in Glg. 3.13 bestimmt. ∂ε ε hingegen wird im folgenden ermittelt:
∂r
∂r ∂
3I20tr − ϕδ 2

=
ε ∂ε
∂ε ε
∂ 1 tr tr 2 2 2
= 3 σ : σ dev − (Y − a1 Y0 I1 − a2 I1 )δ
∂εε 2 dev
∂σ tr

dev tr dY dev ∂I1 ∂I1 2
=3 : σ dev − 2Y − a1 Y 0 − a2 2I1 δ
∂εε de dε ε ∂εε ε
∂ε
v
∂ev ∂I 0tr

dY ∂ev ∂I1
=3σ tr
dev 2G − 2Y + 0tr 2 −
dev ∂I1 ∂ε ε ∂I2 ∂ε ε

∂I1 ∂I1 2
−a1 Y0 − a2 2I1 δ
∂εε ∂εε

dY ∂ev
= 3 − 2Y 2Gσ trdev +
dev ∂I20tr

dY ∂ev ∂I1
+ a1 Y0 + 2a2 I1 − 2Y δ2
dev ∂I1 ∂εε
! !
∂r 4Y λ̇ dY 2Y λ̇t dY δ 2 3K
= 3− 2Gσ tr
dev + a1 Y0 + 2a2 I1 − 1
ε
∂ε Y0 δ 2 dev Y0 dev 1 + 3K λ̇2a∗2
(3.20)
© Nico Diercks 41
3. Implementierung
3.2. ANSYS-spezifische Implementierung
3.2.1. Transformation auf ANSYS-Größen
Bei der Übergabe von Werten zwischen Hauptprogramm und Subroutine sind
einige Besonderheiten zu beachten, da ANSYS teilweise andere Größen verwen-
det, als die, in denen die kontinuumsmechanische Formulierung des Modells
vorgenommen wurde. Daher muss bei der Implementierung darauf geachtet
werden, diese Größen entsprechend zu transformieren, bevor die Berechnungen
durchgeführt werden und diese anschließend wieder zurück zu transformieren.
Tensoren, wie der Spannungs- oder der Verzerrungstensor werden von ANSYS
als Vektoren übergeben. Weiterhin ist zu bedenken, dass ANSYS gemischte
Dehnungen (z.B. ε12 ) immer in Form von Gleitungen (z.B. γ12 ) interpretiert.
Diese Vektoren haben die Form:
   
σ11 ε11
   
 σ22   ε22 
   
 σ   ε 
 33   33 
~σAN SY S =   ~εAN SY S =  
 σ12   γ12 
   
 σ   γ 
 23   23 
σ13 γ13
Weiterhin muss der Steifigkeitstensor 4. Stufe C auf eine 6x6 Matrix transfor-
miert werden. ANSYS fordert diese in der Form:
 
C1111 C1122 C1133 C1112 C1123 C1113
 
 C2211 C2222 C2233 C2212 C2223 C2213 
 
 C
 3311 C3322 C3333 C3312 C3323 C3313 

dsdePl =  
 C1211 C1222 C1233 C1212 C1223 C1213 
 
 C C C C C C 
 2311 2322 2333 2312 2323 2313 
C1311 C1322 C1333 C1312 C1323 C1313
Hierbei wird schnell ersichtlich, dass lediglich 36 Einträge des Tensors C, der
insgesamt über 81 Einträge verfügt, in die ANSYS-Steifigkeitsmatrix eingehen.
© Nico Diercks 42
3. Implementierung
Zu dieser Transformation werden Symmetrieforderungen an C aufgrund der

Symmetrie des Spannungs- und des Verzerrungstensors gestellt.
Ci,j,k,l = Cj,i,k,l
Ci,j,k,l = Ci,j,l,k
Da dies jedoch bei den numerisch berechneten Tensoren nicht immer erfüllt ist,
müssen beide Größen in die Berechnung der 6x6 Steifigkeitsmatrix eingehen.
Weiterhin ist dem Umstand, dass ANSYS Gleitungen stets in der Form des
Schubwinkels γ interpretiert, Rechnung zu tragen. Dazu wird die Berechnung
der Spannungskomponente σi,j betrachtet:
Cij11 ε11 + Cij12 ε12 + Cij13 ε13 +

σi,j = + Cij21 ε21 + Cij22 ε22 + Cij23 ε23 +
+ Cij31 ε31 + Cij32 ε32 + Cij33 ε33
Mit der Symmetrie des Verzerrungstensors εi,j = εj,i folgt:
Cij11 ε11 + (Cij12 + Cij21 )ε12 + (Cij13 + Cij31 )ε13 +

σi,j = + Cij22 ε22 + (Cij23 + Cij32 )ε23 +
+ Cij33 ε33
Der Übergang auf die Schubwinkel 12 γi,j = εi,j liefert
Cij12 +Cij21 Cij13 +Cij31

Cij11 ε11 + 2
γ12 + 2
γ13 +
Cij23 +Cij32
σi,j = + Cij22 ε22 + 2
γ23 +
+ Cij33 ε33
Daraus wird deutlich, dass die Einträge der Schubkomponenten Ci,j,k,l mit k 6=
C +C
l ersetzt werden durch i,j,k,l 2 i,j,l,k . Mithilfe der Symmetrie des Spannungs-
tensors können ähnliche Zusammenhänge gezeigt werden, so dass die Steifig-
keitsmatrix schließlich wie im Anhang unter B.1.3 zusammengesetz wird.
© Nico Diercks 43
3. Implementierung
3.2.2. Die Routine usermat.3d
Zur Verwendung des Materialmodells werden die in Kapitel 3.1 erarbeiteten

Gleichungen in das FEM-Programm ANSYS implementiert. Programmierspra-
che ist dabei Fortran 95; als Entwicklungsumgebung wurde Microsoft Visual
Studio 2003 mit einem Fortran Compiler verwendet. Zum Linken des Codes
sind eine ganze Reihe von ANSYS-Bibiliotheken nötig, die hier nicht weiter
erläutert werden. An dieser Stelle gilt der Dank Herrn Prof. Dr. Armin Fritsch
1
für die Bereitstellung eines Muster-Projektes, mit dem der gesamte Link-
Vorgang automatisch durch Microsoft Visual Studio erfolgt.
Bei der Implementierung selbst ist die Routine so zu programmieren, dass
die aktuellen Gesamtverzerrungen ε , das Lastinkrement ∆ε ε, die plastischen
(n)
Dehnungen des Vorschrittes ε pl und der letzte Wert der Verfestigungsgröße
(n) (n)
ev vom Hauptprogramm in Form der entsprechenden ANSYS-Größen (ev
wird dabei als “user-state-variable” geführt) übergeben werden, die nötigen
Berechnungen durchgeführt und schließlich als Ergebnis die aktuellen plasti-
schen Dehnungen ε pl , den neuen Wert für ev , den Spannungstensor σ und die
Tangentensteifigkeit C in Form der ANSYS-Größen zurück gegeben werden.
Der programmierte Quellcode, der das vorgestellte Materialmodell implemen-
tiert, ist im Anhang B.2 zu finden.
Bei den Tests dieser Implementierung hat sich herausgestellt, dass die Iteration
des plastischen Multiplikators λ̇ in vielen Fällen keine befriedigende Konver-
genz zeigt. Es war zwar zu beobachten, dass die Iteration im Mittel gegen
einen Endwert strebt, der jedoch nicht eingenommen wird; das Ergebnis des
nächsten Iterationsschrittes lag stets über oder unter dem Endwert, ohne dass
sich die Differenz im Laufe der Iteration verringerte.
Um diesem Problem entgegen zu wirken, wurde Glg. 3.12 dahingehend modi-
fiziert, dass die Korrektur des letzten Wertes λ̇(k) mit einem Reduktionsfaktor
multipliziert wird, der im Laufe der Iteration immer kleiner wird.
−1
(k+1) (k) ∂r
λ̇ = λ̇ − r[λ̇(k) ]
∂ λ̇
−1
k ∂r
⇒ λ̇ (k+1)
= λ̇ (k)
− 1 − 0, 95 r[λ̇(k) ] (3.21)
kmax ∂ λ̇
1
Fachhochschule München
© Nico Diercks 44
3. Implementierung
wobei kmax die maximale Zahl der Iterationen bis zum Abbruch bezeichnet.
Durch diese Modifikation wurde erreicht, dass der Korrekturterm anfangs voll
wirksam bleibt, zum Ende der Iteration jedoch auf 20
1
seines Wertes abfällt und
so eine Konvergenz ermöglicht. Die Größe kmax ist standardmäßig auf kmax =
1000 eingestellt, kann durch den User jedoch wie in Kap. 3.2.3.1 dargestellt
verändert werden.
3.2.3. Verwendung der Subroutine durch den Anwender
Um Analysen mit einem auf diese Art programmierten Modell durchzufüh-

ren, sind ein paar Kommandos durch den User notwendig, die im folgenden
aufgeführt sind. Abgesehen davon unterscheidet sich die Arbeit nicht von an-
deren Analysen, die mit bereits vorhandenen Materialmodellen durchgeführt
werden.
3.2.3.1. Aufruf der Routine
Um die usermat.3d mit dem eigenen Inhalt zu verwenden, muss ANSYS über
den ANSYS Product Launcher gestartet werden. Hierbei ist unter der Register-
karte Anpassung/Voreinstellungen die Benutzerdefinierte ANSYS-Programmdatei
zu wählen. Dabei wird als Startdatei die “ANSYS.exe” gewählt, die beim Com-
pilieren und Linken entstanden ist. Anschließend kann die Erstellung des Mo-
dells (Geometrie, Vernetzung etc.) wie gewohnt durchgeführt werden. Der Ver-
weis auf die usermat.3d erfolgt erst bei der Zuweisung des Materials. Wichtig
zu bedenken ist lediglich, dass die Subroutine eine unsymmetrische Tangen-
tensteifigkeit liefert, was durch den User mit dem Kommando
1 NROPT, unsym
einzustellen ist.
© Nico Diercks 45
3. Implementierung
3.2.3.2. Eingabe der Materialparameter
Im folgenden sind die benötigten Kommandos zur Definition des Materials

aufgeführt. Diese Zeilen sind Teil einer Validierungssimulation, die in Kap. 4.2
beschrieben ist. Das vollständige Macro befindet sich im Anhang unter B.3.
1 tb,user ,1,1,10
2 tbtemp,1
3
4 tbdata ,1,1588.7,0.34,36.11,8.08,101.1 !Materialparameter 1−5

5 tbdata ,6,33.96,0.643,0.122,0,0.279 !Materialparameter 6−10
6
7 tb, state ,1,,1
• Mit dem ersten Kommando wird eine Tabelle für ein user defined material
für das Material Nr. 1 mit 1 Temperatur und 10 Datenpunkten deklariert.
• Anschließend wird die erste Temperatur angewählt.
• In den folgenden zwei Zeilen werden die zu dieser Temperatur gehöri-

gen Materialparameter definiert. Dabei ist zu beachten, dass mit einem
tbdata-Kommando maximal 6 Einträge getätigt werden können.
• Zuletzt wird deklariert, dass das Modell über eine user-state-variable

verfügt. Dies ist wichtig, damit die Verfestigungsgröße ev mit übergeben
und abgespeichert wird.
Die Eingabe der Materialparameter erfolgt dabei in folgender Reihenfolge:
Parameter 1-5 : E, ν, Y0 , q, b
(3.22)
Parameter 6-10 : H, a1 , a2 , a∗1 , a∗2
Zusätzlich können mit dem folgenden Kommando

1 tbdata,11,500,20,90,5 !Parameter 11−14
optional noch 4 weitere Parameter eingegeben werden, die zur Optimierung

der Konvergenzgeschwindigkeit und bei Problemen zur Fehlersuche verwendet
werden können. Dabei ist besonders zu beachten, dass bei der Deklarierung
der Tabelle die Anzahl der Datenpunkte auf 14 erweitert werden muss.
Diese zusätzlichen Parameter haben die folgende Bedeutung:
© Nico Diercks 46
3. Implementierung
1. Anzahl der maximalen Iterationen kmax
2. Schrittweite der Iterations-Ausgabe
3. Nr. des Elements, dessen Daten ausgegeben werden
4. Nr. der auszugebenden Integrationspunktes
Die Vorgabe der maximalen Anzahl an Iterationen kann insofern sinnvoll sein,
dass bei einer kleineren Anzahl eine schnellere
Konvergenz
aufgrund des in Glg.
3.21 eingeführten Reduktionsfaktors 1 − 0, 95 kmax in der Iteration auftreten
k
kann.
Sollten Probleme während einer Simulation auftreten, so kann man über die
Parameter 12-14 eine detaillierte Ausgabe aller Zwischenschritte erreichen. In
diesem Beispiel wird vom Element Nr. 90 am Integrationspunkt 5 jeder 20.
Iterationsschritt ausgegeben. Die Ausgabe erfolgt dabei in eine Datei namens
ausgabelog.txt. Ebenfalls wird jede Bisektion mit Angabe des verursachenden
Elements und Integrationspunkt in dieser Datei vermerkt, um die Fehlersuche
zu vereinfachen.
Diese zusätzlichen Modifikationsmöglichkeiten stammen aus dem Debugging
und eine Verwendung sollte eigentlich nicht mehr nötig sein; alle Validierungs-
simulationen wurden ohne diese zusätzlichen Modifikationen erfolgreich durch-
geführt. Dennoch sind diese Optionen nicht deaktiviert worden, damit diese
u.U. bei einer Weiterentwicklung genutzt werden können.
© Nico Diercks 47
4. Modellvalidierung
Nach der ANSYS-Implementierung des Modells ist der nächste Schritt die
Validierung durch Experimente. Dazu wurden im ersten Schritt einige Experi-
mente mit einfacher Belastung ausgewertet, deren Daten aus [12] entnommen
wurde. Anschließend wurden diese Experimente in ANSYS simuliert, um zu
überprüfen, ob die Ergebnisse der FE-Simulation das Experiment nachbilden
können.
Im zweiten Schritt wurden kombinierte Belastungen simuliert. Dazu wurde das
von Schlimmer und Mahnken durchgeführte und in [9] vorgestellte Experiment
mit stumpf geklebten Rohren unter verschiedener Beanspruchung (Zug und
Torsion) nachgebildet, wodurch überprüft werden soll, ob die FE-Rechnung in
der Lage ist, die von Schlimmer und Mahnken beobachteten Effekte (Abhän-
gigkeit der Festigkeit vom Lastfall) zu modellieren.
In allen Experimenten wurde der Strukturklebstoff Betamate 1496˝ verwen-

det, dessen Materialparameter von Schlimmer und Mahnken, wie in Tabelle
4 aufgeführt, bestimmt wurden. Für die Bestimmung der Parameter gibt es
mehrere Verfahren, ein Beispiel dafür kann in [8] gefunden werden.
E [MPa] = 1588,7 H [MPa] = 33,96

ν = 0,34 a1 = 0,643
Y0 [MPa] = 36,11 a2 = 0,122
q [MPa] = 8,08 a∗1 = 0
b = 101,1 a∗2 = 0,279
Tabelle 4.1.: Materialparameter
© Nico Diercks 48
4.1. Experimente mit einfacher

Beanspruchung
4.1.1. Zugversuch
Einer der einfachsten Versuche zur Bestimmung von Materialkennwerten ist

der einachsige Zugversuch. Aus diesem Grund wird er gern angewendet, um
einen ersten Eindruck vom Materialverhalten zu gewinnen.
Dazu wird eine Probe, deren Geometrie in der EN ISO 527-2 vorgegeben ist
≥ 150
60
10 20
L0 =50
106-120
Abbildung 4.1.: Zugprobe nach EN ISO 527-2
(Abmessungen gem. Abb. 4.1, Dicke t = 4mm) , in Längsrichtung verformt und

deren Last-Verformungskurve aufgenommen. Dabei wird die Dehnung des in-
neren Probenteils L0 gemessen und die Spannung auf den inneren Querschnitt
bezogen wird.
Für die FE-Simulation des Versuches kann die mehrfache Symmetrie der Probe
ausgenutzt werden, so dass insgesamt nur 1/8 der Probe mit den entsprechen-
den Randbedingungen modelliert werden muss.
In der Abb. 4.2 ist die Last-Verformungskurven der Simulation, sowie die
Datenpunkte des Experimentes dargestellt. Bei der Betrachtung wird schnell
deutlich, dass die FE-Rechnung in der Lage ist, diesen Versuch gut zu simulie-
ren; lediglich einige Prozent Abweichung sind zu erkennen. Weiterhin bietet die
FE-Simulation auch die Möglichkeit, lokale Effekte aufzuzeigen. Als Beispiel
ist in Abb. 4.3 die lokale Verteilung der plastischen Dehnung dargestellt.
© Nico Diercks 49
35
Zugspannung σz [MPa] 30
25
20
15
10
5 Simulation
Experiment
0
0 0.025 0.05 0.075 0.1
εz
Abbildung 4.2.: Last-Verformungskurve des Zugversuches
Abbildung 4.3.: plastische Dehnung der Zugprobe
4.1.2. Kompressionsversuch
In Ergänzung zum Zugversuch wurde als nächstes ein Kompressionsversuch

betrachtet. Bei diesem Versuch wurde eine zylindrische Probe, wie in Abb. 4.4
dargestellt auf Druck belastet, um die Last-Verformungskurve zu bestimmen.
Die gemessenen Daten sind in Abb. 4.5 mit Kästchen dargestellt. Zum Ver-
gleich sind in derselben Graphik mit der durchgehenden Linie die Ergebnisse
der Simulation dargestellt, wobei hier aufgrund von Symmetriebetrachtungen
wieder nur ein Teil der Probe modelliert werden musste. Die Auswertung dieser
Graphik zeigt wie beim Zugversuch, dass die Simulation die experimentellen
Daten gut nachbilden kann. Die Abweichungen sind in der gleichen geringen
© Nico Diercks 50
20
25
Abbildung 4.4.: Druckprobe
Größenordnung, wodurch das Modell auch in diesem Versuch überzeugen kann.
60
Druckspannung σz [MPa]
50
40
30
20
10 Simulation
Experiment
0
0 0.025 0.05 0.075 0.1
εz
Abbildung 4.5.: Last-Verformungskurve des Kompressionsversuches
4.1.3. Zugscherversuch
Als letztes Experiment mit einfacher Beanspruchung wurde ein Zugscherver-

such modelliert. Dies ist ein Standard-Versuch zur Untersuchung von Klebstof-
fen, der in der DIN EN 1465 genau beschrieben ist. Weitere Informationen zu
diesem Prüfverfahren sind in [2] zu finden.
Die Zugscherprobe, die hierfür verwendet wird, ist in Abb. 4.6 dargestellt, wo-
© Nico Diercks 51
db
dk
mm
25
10 mm
Abbildung 4.6.: Zugscherprobe
bei in diesem Fall eine etwas von der DIN abweichende Geometrie verwendet
wurde. Nach der Norm verfügt eine dünne Zugscherprobe über eine Über-
lappungslänge von 12,5mm, während diese bei einer dicken Probe nur 5mm
beträgt. Im betrachteten Versuch wurde eine dicke Probe untersucht, wobei
die Überlappungslänge abweichend von der Norm mit 10mm ausgeführt war.
Die Dicke der Klebschicht beträgt hierbei dK = 0, 9mm, während die Bleche
mit dB = 6mm ausgeführt sind. Die Verwendung der dicken Probe bringt den
Vorteil, dass die Belastung der Klebschicht weitesgehend in Schubbeanspru-
chung ausfällt. Bei dünnen Blechen kommt es zu Biegedeformation der Bleche,
wodurch die Klebefläche zusätzlich zum Schub noch mit Normalspannung be-
lastet wird
Bei der Simulation des Zugscherversuches konnte lediglich eine Symmetrie-
10
6
F [kN]
2
Simulation
Experiment
0
0 0.2 0.4 0.6 0.8
u [mm]
Abbildung 4.7.: Last-Verformungskurve des Zugscherversuches
© Nico Diercks 52
ebene ausgenutzt werden, wodurch das FE-Modell in diesem Fall über mehr
Freiheitsgrad als bei den letzten Versuchen verfügte. Die Betrachtung der Gra-
phik 4.7 fällt auf, dass die Ergebnisse der Simulation bei diesem Versuch nicht
die gleiche Güte wie im Zug- oder Druckversuch zeigen. Im Bereich der einset-
zenden plastischen Deformationen liegen die Spannungswerte der Simulation
bis zu 15% unter den im Experiment gemessenen. Im Bereich etwas größe-
rer Deformationen stimmen die Werte jedoch wieder überein. Somit sind die
Daten der Simulation dennoch gut verwertbar, da die Spannungen gerade im
kritischen Bereich (die Proben versagten bei ca. 0,9mm Auslenkung) recht gut
zu den MEssdaten passen. Weiterhin ist zu bemerken, dass die Simulation bei
einer vorgegebenen Belastung größere Verschiebungen liefert, so dass die Be-
rechnungen eher konservative Näherungen darstellen.
Als Abschluss dieser Diskussion ist in Abb. 4.8 noch die Verteilung der Span-
Abbildung 4.8.: v. Mises-Spannung der Zugscherprobe
nungen innerhalb der Klebschicht der Zugscherprobe dargestellt. Somit bietet

die Simulation dem Konstrukteur die Möglichkeit, kritische Punkte der Kon-
struktion festzustellen.
© Nico Diercks 53
4.2. Kombinierte Zug-Schubbelastung
Bei diesem Versuch handelt es sich um ein Experiment mit komplexer zweiach-
siger Beanspruchung, dass von Schlimmer und Mahnken durchgeführt wurde.
Bei diesem Experiment wurde die Beobachtung gemacht, dass die Festigkeit ei-
nes Klebstoffs vom Beanspruchungszustand abhängt. Die experimentellen Da-
ten wurden aus [9] entnommen, um sie anschließend mit den aus der ANSYS-
Simulation gewonnenen zu vergleichen.
4.2.1. Versuchsaufbau
Um eine Klebeverbindung sowohl unter Zug als auch unter Schub zu bean-
spruchen, wurden zwei Rohre (Außendurchmesser 60mm, Wandstärke 5mm)
mit einer Klebefuge von 2mm stumpf verklebt und anschließend unter Zug und
Torsion belastet, wie es in Abb. 4.9 skiziert ist. Die Messung der Deformation
50 mm
vo uo
2 mm
vu uu
60 mm
Abbildung 4.9.: Experimenteller Aufbau zur kombinierten Belastung
erfolgt über mehrere Wegaufnehmer. Die axiale Deformation der Klebschicht

v folgt aus der Differenz der Wegaufnehmer vo und vu zu v = vo − vu . Die
© Nico Diercks 54
1/2 γ
yz
Nr. Lastfall α=
εz
1 reiner Zug 0
2 reine Torsion ∞
3 überw. Zug 0,5
4 überw. Torsion 2
Tabelle 4.2.: Lastfälle
tangentiale Verformung u wird ganz analog zu v aus der Differenz von uo und
uu zu u = uo − uu bestimmt. Aus diesen Verschiebungen können anschließend
die Dehnungen unter der Annahme kleiner Deformationen bzw. geometrischer
Linearität berechnet werden
v u
εaxial = εz = γscher = γyz = ,
d d+v
wobei d die Dicke der Klebschicht mit d = 2mm bezeichnet.
4.2.2. Versuchsdurchführung
Zur Untersuchung verschiedener Lastfälle wurden 4 verschiedene Belastungen

(siehe Tabelle 4.2) aufgebracht. Zum einen wurden reine Zug- oder Torsions-
belastungen aufgebracht, zum anderen kombinierte Belastungen mit überwie-
gender Zug- oder überwiegender Torsionsbelastung. Das Verhältnis der Ver-
formungen wird dabei mit α bezeichnet.
4.2.3. Auswerung
Nach Simulation des Versuches in ANSYS sind die Last-Verformungskurven in

Abbildung 4.10 dargestellt. Die durchgezogenen Linien sind dabei die Ergeb-
nisse der Simulation, während die experimentellen Daten von Schlimmer und
Mahnken als Kästchen dargestellt sind.
Bei Betrachtung der Graphen erkennt man schnell die von Schlimmer und
Mahnken beobachtete Abhängigkeit der Festigkeit vom Lastfall; die Fließspan-
nung des Materials zeigt für alle Lastfälle deutlich unterschiedliche Werte. Ein
© Nico Diercks 55
40
Zugspannung σz [MPa]
35
30
25
20
15
10 α=0
5 α = 0, 5
α=2
0
0 0.05 0.1 0.15
εz
30
Schubspannung σyz [MPa]
25
20
15
10
5 α=∞
α=2
α = 0, 5
0
0 0.1 0.2 0.3
γyz
Abbildung 4.10.: Spannungsverläufe
Vergleich der experimentellen Daten mit denen der Simulation, ermöglicht ei-
ne zuverlässige Beurteilung des Modells: Die Spannungen der FE-Rechnung
stimmen tendenziell sehr gut mit denen des Experiments überein. Besonders
eindrucksvoll wird die Abhängigkeit der Festigkeit vom Lastfall nachgebildet.
Einige leichte Schwächen zeigen sich lediglich bei der Berechnung der Zugspan-
nung im torsionsdominierten Lastfall (α = 2) bei größeren Dehnungen. Dies
ist aber durchaus hinnehmbar, da zum einen die hier dominierenden Schub-
spannungen sehr gut nachgebildet werden. Zum anderen ist zu bedenken, dass
im Bereich der Ungenauigkeit die Axialdehnung bereits εz ≥ 6% und somit
die Schubverzerrungen sogar schon im Bereich εyz ≥ 24% liegen. Damit dürfte
deutlich sein, dass die Annahme kleiner Deformationen bereits mehr als frag-
würdig ist.
© Nico Diercks 56
Neben den globalen Last-Verformungskurven bietet die FE-Simulation natür-

lich auch die lokale Verteilung von Spannungen und Dehnungen und ermöglicht
so, kritische Punkte der Konstruktion zu ermitteln. In Abbildung 4.11 wird als
Beispiel die Verteilung der axialen plastischen Dehnung innerhalb der Klebfuge
beim Lastfall α = 0, 5 dargestellt.
Abbildung 4.11.: axiale, plastische Dehnung
© Nico Diercks 57
5. Fazit und kritische Bewertung
5. Fazit und kritische Bewertung
In der Gesamtbetrachtung der Validierungssimulationen wird deutlich, dass

das vorliegende Materialmodell nach Schlimmer und Mahnken die experimen-
tell gewonnenen Daten in großen Bereichen gut simulieren kann. Sowohl die
plastischen Deformationen als auch der Effekt der wechselnden Festigkeit un-
ter kombinierter Beanspruchung werden durch das Modell dargestellt. Ledig-
lich bei der Simulation des Zugscherversuches treten im Anfangsbereich der
Plastizität relevante Abweichungen auf, die im kritischen Bereich jedoch wie-
der kleiner werden. In den weiteren Versuchen weicht die Simulation erst bei
größeren Deformationen von den experimentellen Daten ab.
Daher sollte es das nächste Ziel sein, dieses Modell für die Berücksichtigung
der besonderen Effekte beim Auftreten großer Deformationen zu erweitern. Ein
mögliches Konzept dazu ist mit der Formulierung einer Freie-Energie-Funktion
in [7] zu finden. Es sollte möglich sein, diese Vorgehensweise auf das Modell von
Schlimmer und Mahnken zu übertragen. Weiterhin ließe sich die Genauigkeit
des Modells unter Umständen auch dadurch verbessern, indem man die dritte
Invariante des Spannungstensors in die Fließfunktion mit eingehen lässt.
Ein weiteres Ziel, neben der Verbesserung der Berechnung von Spannungs-
und Verzerrungszuständen, sollte die Formulierung von Versagenshypothesen
sein. Somit ließe sich nicht nur die Festigkeitsberechnung der intakten Struktur
durchführen, sondern zusätzlich das lokale oder globale Versagen der Struktu-
ren berechnen.
Trotz der genannten Verbesserungsmöglichkeiten kann das Ziel dieser Arbeit,
dem Konstruktionsingenieur mit der Implementierung des Schlimmer-Mahnken-
Materialmodells ein zuverlässiges Werkzeug zur Festigkeitsberechnung verkleb-
ter Strukturen zur Verfügung zu stellen, als erreicht gewertet werden. Durch
die Verwendung dieser Routine wird es in Zukunft sicher möglich sein, Ex-
perimente einzusparen und durch Simulationen zu ersetzen, um so die rasant
steigenden Entwicklungskosten komplizierter Strukturen zu senken.
© Nico Diercks 58
Literaturverzeichnis
1 Crisfield, M.A.: Non-linear Finite Element Analysis of Solids and Struc-

tures, Band 1. John Wiley & Sons, 1991. 2.2.1, 3.1.2.3
2 Czarnecki, Dr. Jürgen von: Einführung in die Klebtechnik. Vorle-

sungsscriptum, 2008. 4.1.3
3 Ehlers, W., Prof. Dr.-Ing.: Ergänzungsscriptum: Tensorrechnung -

Eine Einführung. Universität Stuttgart Institut für Mechanik, 2007. A,
A.1.2
4 Emmerling, F.A., Univ.-Prof. Dr.-Ing.: Vorlesungsscriptum Finite

Elemente I. Universität der Bundeswehr München, Fakultät für Luft- und
Raumfahrttechnik, Institut für Mechanik, Oktober 2003. 2.2.1
5 Emmerling, F.A., Univ.-Prof. Dr.-Ing.: Differentialgleichungen der

Technischen Mechanik IV. S. 114-146, Universität der Bundeswehr Mün-
chen, Fakultät für Luft- und Raumfahrttechnik, Institut für Mechanik,
Oktober 2005. 2.1.1
6 J.C. Simo, T.J.R. Hughes: Computational Inelasticity. Springer-Verlag,

1998. 2.2.1, 3.1.2.3
7 Lion, A., Univ.-Prof. Dr.-Ing.: Constitutive modelling in finite thermo-

viscoplasticity: a physical approach based on nonlinear rheological models.
International Journal of Plasticity, 16:469–494, 2000. 5
8 Mahnken, R.: Identification of material parameters for constitutive equa-

tions. In: Stein E, de Borst, Hughes (Herausgeber): Encyclopedia of
Computational Mechanics, Band 2. Wiley, 2004. 4
9 Mahnken R., Schlimmer M.: Simulation of strenght difference in elasto-

plasticity for adhesive materials. International journal for numerical me-
thods in engeneering, 63:1461–1477, März 2005. 2.3.1, 4, 4.2
© Nico Diercks 59
10 Middendorf, Peter, Dr.-Ing.: Ausgewählte Probleme der nichtlinearen

Finite Elemente Methode. S. 71-92, EADS Corporate Research Center
Germany, 2007. 2.1.3.1
11 Schlimmer, M.: Fließverhalten plastisch kompressibler Werkstoffe. Dok-

torarbeit, RWTH Aachen, 1974. 2.3.3.2
12 Wocke, C. P.: Betamate Material Data for Engeneering Analysis. DOW

Automotive, Februar 2006. 4
© Nico Diercks 60
Eidesstattliche Erklärung
Eidesstattliche Erklärung
Ich versichere hiermit, dass ich meine Diplomarbeit mit dem Thema
Validierung und ANSYS-Implementierung des Schlimmer-Mahnken

Materialmodells zur Festigkeitsberechnung von Klebeverbindungen
selbständig verfasst und keine anderen als die angegebenen Quellen und Hilfs-
mittel benutzt habe. Die Arbeit wurde bisher keiner anderen Prüfungsbehörde
vorgelegt und auch nicht veröffentlicht.
Neubiberg, den 15. Dezember 2008
Nico Diercks
© Nico Diercks 61
A. Mathematische Grundlagen
An dieser Stelle wird eine kurze Einführung in die Tensorrechnung gegeben,

wobei die wichtigsten Grundlagen, die zum Verständnis dieser Arbeit benötigt
werden, aufgeführt sind. Dabei wird natürlich kein Anspruch auf eine vollstän-
dige Beschreibung der Tensorrechnung erhoben. Zum tieferen Studium wird
daher die Lektüre von [3] empfohlen, woraus auch diese Beispiele stammen.
A.1. Tensoralgebra
A.1.1. Einige Konventionen
In den folgenden Ausführungen und in der gesamten Arbeit werden einige Kon-
ventionen verwendet, die an dieser Stelle erläutert werden. Bei Verwendung der
Indexschreibweise werden die Ausdrücke gern durch die Einstein’sche Summa-
tionskonvention verkürzt dargestellt. Dabei gilt, dass bei Auftreten von zwei
gleichen Indizes über diese von 1 bis 3 summiert wird:
3
X
ai~ei = ai~ei = a1~e1 + a2~e2 + a3~e3 (A.1)
i=1
Weiterhin wird das Kronecker-Delta eingeführt, mit dem das einfache Skalar-
produkt von Basisvektoren gut beschrieben werden kann:
(
1 für i = k
~ei · ~ek = δik = (A.2)
0 für i =
6 k
Der Vorteil bei der Verwendung des Kronecker-Deltas liegt neben der verkürz-
ten Schreibweise darin, dass bestimmte Indizes ausgetauscht werden können,
© Nico Diercks i
wie im folgenden Beispiel anhand des einfachen Skalarproduktes zweier Vek-

toren gezeigt wird:
~a · ~b = (ai~ei ) · (bk~ek ) = ai bk (~ei · ~ek ) = ai bk δik
Da ohnehin alle Terme mit i 6= k weg fallen, kann i = k gesetzt werden und es
gilt
ai bk δik = ai bi .
A.1.2. Tensorbegriff
Nach [3] ist ein Tensor T 2. Stufe eine lineare Abbildung, die einem Vektor ~u
eindeutig einen Vektor w~ zuordnet:
T~u = w
~ (A.3)
Dabei wird der Tensor mithilfe des dyadischen Produktes ⊗ beschrieben durch
T = (~a ⊗ ~b) mit der Eigenschaft
(~a ⊗ ~b)~c =: (~b · ~c)~a . (A.4)
Die Basisdarstellung des Tensors erfolgt dabei durch
T := (~a ⊗ ~b) = (ai~ei ) ⊗ (bk~ek ) = ai bk (~ei ⊗ ~ek )
Allgemein lässt sich ein Tensor also in folgender Form beschreiben:
T := tik (~ei ⊗ ~ek )
wobei (~ei ⊗ ~ek ) die Tensorbasis und tik die Koeffizientenmatrix mit neun un-
abhängigen Einträgen darstellt mit
 
t11 t12 t13
tik =  t21 t22 t23 
 
t31 t32 t33
© Nico Diercks ii
A.1.3. Produkte von Tensoren
Multiplikation mit einem Vektor
Es gilt die Definition der linearen Abbildung gemäß Glg. A.3 mit
T~u = w
~
In Basisdarstellung lautet dies
T~u = (tik~ei ⊗ ~ek )(uj ~ej ) = tik uj (~ei ⊗ ~ek )~ej
mit der Eigenschaft des dyadischen Produktes nach A.4 und dem Skalarpro-
dukt zweier Vektoren nach A.2
T~u = tik uj δkj ~ei = tik uk~ei
Skalarprodukt von Tensoren
Das Skalarprodukt zweier Tensoren wird im Rahmen dieser Arbeit mit einem
Doppelpunkt [:] gekennzeichnet. Dies geschieht in Abgrenzung zum Tensorpro-
dukt, das hier ohne zusätzliches Zeichen verwendet wird. Diese Nomenklatur
liegt darin begründet, dass das Ergebnis des Skalarproduktes stets ein Skalar
und das des Tensorproduktes stets ein Tensor ist.
In anderer Literatur wird dieses Produkt mitunter auch als doppeltes Skalar-
produkt bezeichnet, während das Tensorprodukt als einfaches Skalarprodukt
eingeführt wird.
A : B = aik (~ei ⊗ ~ek ) : brs (~er ⊗ ~es )
= aik brs (~ei ⊗ ~ek ) : (~er ⊗ ~es )
= aik brs δir δks
= aik bik
© Nico Diercks iii

Das Tensorprodukt
Schließlich soll an dieser Stelle noch das Tensorprodukt (oder einfaches Ska-
larprodukt) vorgestellt werden, dessen Ergebnis stets ein Tensor ist. In Basis-
darstellung lautet dies
AB = aik (~ei ⊗ ~ek )brs (~er ⊗ ~es )

= aik brs (~ei ⊗ ~ek )(~er ⊗ ~es )
= aik brs δkr (~ei ⊗ ~es )
= aik bks (~ei ⊗ ~es )
A.1.4. Tensoren höherer Stufe
Darstellung
Bei der Behandlung von Tensoren höherer Stufe ist es hilfreich, sich am Beispiel
der linearen Abbildung (Glg. A.3) zu orientieren und dieses auf höherwertige
Tensoren zu verallgemeinern. Dabei hat diese Abbildung die gleiche verjün-
gende Eigenschaft, d.h. das Ergebnis ist ein Tensor niedrigerer Stufe. Es gilt
n m n−m
A B= C
wobei die Zahl über dem Symbol die Tensorstufe angibt. Ist die Stufe nicht
angegeben, so handelt es sich um einen Tensor 2. Stufe. Im folgenden Beispiel
wird ein Tensor 4. Stufe mit einem Tensor 2. Stufe multipliziert, wobei das
Ergebnis ein Tensor 2. Stufe ist.
4
A B = aijkl (~ei ⊗ ~ej ⊗ ~ek ⊗ ~el )brs (~er ⊗ ~es )
= aijkl brs δkr δls (~ei ⊗ ~ej )
= aijkl bkl (~ei ⊗ ~ej )
= cij (~ei ⊗ ~ej ) = C
© Nico Diercks iv
Fundamentaltensoren 4. Stufe
Fundamentaltensoren 4. Stufe werden aus dem dyadischen Produkt von 2 Ein-

heitstensoren 2. Stufe gebildet. In dieser Arbeit werden zwei dieser Fundamen-
taltensoren des öfteren verwendet:
1. die Spurbildende Abbildung 1 ⊗ 1

4 23
2. der Einheitstensor 4. Stufe I := 1 = (1 ⊗ 1) T
ij
Wobei (·)T eine spezielle Transposition, d.h. eine Vertauschung des i-ten mit
dem j-ten Basissystems, bezeichnet. Diese Fundamentaltensoren werden wir
folgt gebildet:
1 ⊗ 1 = (~ei ⊗ ~ei ) ⊗ (~ej ⊗ ~ej )

I = ~ei ⊗ ~ej ⊗ ~ei ⊗ ~ej
Die genannten Eigenschaften zeigen sich durch:
(1 ⊗ 1) T = (~ei ⊗ ~ei ⊗ ~ej ⊗ ~ej )trs (~er ⊗ ~es )

= trs δjr δjs (~ei ⊗ ~ei ) = tjj (~ei ⊗ ~ei )
= Sp(T)1
23
(1 ⊗ 1) T T = (~ei ⊗ ~ej ⊗ ~ei ⊗ ~ej )trs (~er ⊗ ~es )
= trs δir δjs (~ei ⊗ ~ej )
= tij (~ei ⊗ ~ej ) = T
© Nico Diercks v
A.2. Tensoranalysis
An dieser Stelle wird eine kurze Einführung in die Tensoranalysis gegeben,

wobei die wichtigsten Vorgehensweisen, wie sie insbesondere bei der Herleitung
der Tangentensteifigkeit in Kapitel 3.1.2.5 benötigt werden, kurz dargestellt
werden.
A.2.1. Einführung des Funktionsbegriffes
Bei der Bezeichnung verschiedener Funktionen der Tensoranalysis gilt:
φ(...) : skalarwertige Funktion  skalaren Variablen

 

 
~v (...) : vektorwertige Funktion von vektoriellen Variablen
T(...) : Tensorwertige Funktion tensoriellen Variablen

 

Beispiel: φ(C) : skalarwertige Tensorfunktion
A.2.2. Tensorielle Ableitungen von Funktionen
Ableitung eines Skalares nach einem Tensor
Das Ergebnis einer Ableitung nach einem Tensor 2. Stufe ist ebenfalls ein
Tensor 2. Stufe. Für die Ableitung einer skalarwerigen Tensorfunktion α(A)
gilt:
∂φ ∂φ
= ~ei ⊗ ~ej
∂A ∂aij
Die Ableitung eines Produktes zweier skalarwertiger Funktionen berechnet sich
in gewohnter Art nach
∂(φ ψ) ∂φ ∂ψ
=ψ +φ
∂A ∂A ∂A
Ableitungen von Tensorfunktionen
Um eine Ableitung der (n × n) Komponenten einer Tensorfunktion T(A) nach

allen (n × n) Komponenten eines anderen Tensors A darzustellen, ist das Er-
© Nico Diercks vi
gebnis dieser Ableitung ein Tensor 4. Stufe. Analog zur Ableitung der skalar-
wertigen Funktion gilt:
∂T ∂tij
= ~ei ⊗ ~ej ⊗ ~ek ⊗ ~el
∂A ∂akl
Bei der Ableitung eines Produktes einer skalarwertigen Funktion mit einer
tensorwertigen Funktion ist folgende Besonderheit zu berücksichtigen:
∂(φT) ∂φ ∂T
=T⊗ +φ
∂A ∂A ∂A
Außerdem gilt:
∂T
= I
∂T
∂I1 (T)
= 1
∂T
© Nico Diercks vii

B. Quellcode
B. Quellcode
B.1. Verzeichnis der verwendeten Variablen
B.1.1. Lokale Variablen
Name Zeichen Kommentar

Materialparameter
E E Elastizitätsmodul
ny ν Querkontraktionszahl
y0 Y0 Anfangsfließspannung
q q
b b
H H
a1 a1
a2 a2
a1s a∗1
a2s a∗2
G G Schubmodul
Ka K Kompressionsmodul
Tensoren
C C Steifigkeitstensor 4. Stufe
dlde ∂ λ̇
∂ε
drde ∂r
∂ε
eps ε 3x3 Verzerrungstensor
epsdev εdev 3x3 Verzerrungstensor, deviatorisch
Fortsetzung auf der nächsten Seite
© Nico Diercks viii

B. Quellcode
Lokale Variablen: Fortsetzung

epspln εnpl 3x3 plastischer Verzerrungstensor bei Aufruf der
Routine
epsplndev εnpl,dev 3x3 plastischer Verzerrungstensor, deviatorisch bei
Aufruf der Routine
epsPlt εpl 3x3 plastischer Verzerrungstensor
ident2 1 Einheitstensor 2. Stufe
ident4 I Einheitstensor 4. Stufe
idev Idev
s σ 3x3 Spannungstensor
sdevtr σ tr
dev 3x3 Spannungstensor, trial, deviatorisch
str σ tr 3x3 Spannungstensor, trial
t1d1 1⊗1
skalare Größen
didl ∂I1
∂ λ̇
dedl ∂e
∂ λ̇
del δ
dphidl ∂ϕ
∂ λ̇
ev ev innere Variable bei Ende der Routine
evn env innere Variable bei Aufruf der Routine
I1 I1 (σ) erste Invariante des Spannungstensors
I1tr I1tr (σ tr ) erste Invariante des Trial-Spannungstensors
I2tr I20tr (σ tr
dev ) zweite Invariante des Trial-Spannungstensors (de-
viatorisch)
jj J = ∂∂rλ̇
lalt λ̇k plastischer Multiplikator
lneu λ̇k+1 plastischer Multiplikator
phi ϕ Teil der Fließfunktion
grphi Φ Fließfunktion
t t(n)
trckeps Sp(εn ) Spur des Verzerrungstensors
trckepspl Sp(εnpl ) Spur des plastischen Verzerrungstensors
© Nico Diercks ix
B. Quellcode
Lokale Variablen: Fortsetzung

tnp1 t
R R(ev ) spannungsartige innere Variable
rdl r[λ̇]
Y Y Teil der Fließfunktion
ys 0
y = de dY
v
programminterne Größen
elemaus Nr. des Elements, dessen Daten ausgegeben wer-
den
intpointout Nr. des auszugebenden Integrationspunktes
i,j,k,l,n Zählindizes
itausgabe Gibt an, wie viele Iterationsschritte ausgegeben
werden
maxit maximale Anzahl an Iterationen
sqtiny kleinste Abweichung zw. 2 Werten zum Beenden
der Iteration
yield binäre Größe, =„true"bei plastischer Deformation
B.1.2. Durch ANSYS definierte Variablen
Name Kommentar
verwendete Input-Größen
matId Material-Nr.
elemId Element-Nr.
kdomIntPt Integrationspunkt-Nr.
ldstep Loadstep-Nr.
isubst Substep-Nr.
Strain Verzerrungen als ANSYS-Vektor
dStrain Verzerrungsinkremente als ANSYS-Vektor
prop Materialparameter als Vektor
verwendete Input-Output-Größen
stress Spannungen als ANSYS-Vektor
© Nico Diercks x
B. Quellcode
Durch ANSYS definierte Variablen: Fortsetzung

Name Kommentar
ustatev User state variables als Vektor
ustatev(1) Verfestigungsgröße ev
verwendete Output-Größen
keycut =1, falls die Iteration keine Konvergenz erreicht
dsdePl 6x6 Steifigkeitsmatrix
© Nico Diercks xi
© Nico Diercks
B. Quellcode
B.1.3. ANSYS-Steifigkeitsmatrix
Der Steifigkeitstensor 4. Stufe C wird wie folgt zur ANSYS-spezifischen 6x6 Steifigkeitsmatrix dsdePl transformiert:
 1 1 1

C1111 C1122 C1133 2
(C1112+ C1121 ) 2
(C1123+ C1132 ) 2
(C1113
+ C1131 )
1 1 1
C2211 C2222 C2233 (C2212+ C2221 ) (C2223+ C2232 ) (C2213
+ C2231 )
 
2 2 2
 
 
1 1 1

 C3311 ! C3322 ! C3333 ! 2
(C3312+ C3321 ) ! 2
(C3323+ C3332 ) ! 2
(C3313
+ C3331 ) ! 

C1211 C1222 C1233 C1212 + C1221 C1223 + C1232 C1213 + C1231
 

 1 1 1 1 1 1 

2 2 2 4 4 4

 +C2111 ! +C2122 ! +C2133 ! +C2112 + C2121 ! +C2123 + C2132 ! +C2113 + C2131 ! 
C2311 C2322 C2333 C2312 + C2321 C2323 + C2332 C2313 + C2331
 

 1 1 1 1 1 1 

2 2 2 4 4 4

 +C3211 ! +C3222 ! +C3233 ! +C3212 + C3221 ! +C3223 + C3232 ! +C3213 + C3231 ! 
C1311 C1322 C1333 C1312 + C1321 C1323 + C1332 C1313 + C1331
 
 1 1 1 1 1 1 
2 2 2 4 4 4
+C3111 +C3122 +C3133 +C3112 + C3121 +C3123 + C3132 +C3113 + C3131
Somit gehen alle 81 Einträge von C in die Matrix dsdePl ein.

xii
B. Quellcode
B.2. Quellcode der usermat.3d

1 EXTERNAL wrinqr, vzero
2
3 INTEGER i , j , k, l , n, maxit, iott , wrinqr, itausgabe, elemaus, intpointout

4 DOUBLE PRECISION E, ny, y0, q, b, H, a1, a2, a1s, a2s, lneu, lalt , ys, dedl, rdl ,
didl , ka, I1 , I1tr , I2tr G, del, phi, grphi, t , tnp1, ev, evn, R, dphidl, Y, jj ,
trckeps , trckepspl , red, epspln(3,3) , epsplndev(3,3), epsPlt(3,3) , eps(3,3) , epsdev
(3,3), str (3,3) , sdevtr(3,3) , s (3,3) , ident2 (3,3) , drde(3,3) , dlde(3,3) , C(3,3,3,3) ,
ident4 (3,3,3,3) , t1d1 (3,3,3,3) , idev (3,3,3,3)
5 Logical yield
6
7 !∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗
8
9 iott = wrinqr(WR_OUTPUT)
10 keycut= 0
11
12 !get material parameters

13 E= prop(1)
14 ny= prop(2)
15 y0= prop(3)
16 q= prop(4)
17 b= prop(5)
18 H= prop(6)
19 a1= prop(7)
20 a2= prop(8)
21 a1s= prop(9)
22 a2s= prop(10)
23
24 ka= E/(3.d0∗(1.d0−2.d0∗ny))
25 G= E/(2.d0∗(ny+1))
26
27 !printout−definition for testing purpose

28 if (nprop<11) then
29 maxit=1000
30 itausgabe=10
31 elemaus=0
32 intpointout=0
33 else
34 maxit=prop(11)
35 itausgabe=prop(12)
36 elemaus=prop(13)
© Nico Diercks xiii

B. Quellcode
37 intpointout=prop(14)
38 end if
39
40 ! Definition of a log−file
41 if (isubst .eq .1. and.kdomIntPt.eq.itpointout.and.elemID.eq.elemaus) then
42 open (unit =20 , file =’ausgabelog.txt’, status =’replace’ , action =’write’)
43 else if (elemID.eq.elemaus.and.kdomIntPt.eq.itpointout) then
44 open (unit =20 , file =’ausgabelog.txt’, status =’old’ , action =’write’ , position=’
append’)
45 end if
46
47 !Define tsstif (1) for hourglass stiffness

48 tsstif (1) = G
49
50 !get total strain tensor at end of increment

51 eps(1,1)= Strain(1)+dStrain(1)
52 eps(2,1)= half∗(Strain(4)+dStrain(4))
54 eps(1,2)= eps(2,1)
57 eps(1,3)= eps(3,1)
58 eps(2,3)= eps(3,2)
60
61 !get plastic strain tensor at beginning of increment

62 epspln(1,1)= epspl(1)
63 epspln(2,1)= half∗epspl(4)
65 epspln(1,2)= epspln(2,1)
71
72 ! Definition of unit−tensors (2. and 4. order)

73 call vzero(ident2,9)
74 do i=1,3
75 ident2( i , i )= 1.d0
76 end do
77
© Nico Diercks xiv

B. Quellcode
78 call vzero(ident4,81)
79 do i=1,3
80 do j=1,3
81 ident4( i , j , i , j )= 1.d0
82 end do
83 end do
84
85 ! Definition of fourth−order tensor 1 dyad 1

86 call vzero(t1d1,81)
87 do i=1,3
88 do j=1,3
89 t1d1(i , i , j , j )= 1.d0
90 end do
91 end do
92
93 ! Definition of fourth−order tensor I, deviatoric

94 do i=1,3
95 do j=1,3
96 do k=1,3
97 do l=1,3
98 idev( i , j ,k, l )= ident4(i , j ,k, l )−third∗t1d1(i,j ,k, l )
99 end do
100 end do
101 end do
102 end do
103
104 !Deviatoric split of strains

105 trckeps= eps(1,1)+eps(2,2)+eps(3,3)
106 do j=1,3
107 do i=1,3
108 epsdev(i, j )= eps(i , j )−third∗trckeps∗ident2(i, j )
109 end do
110 end do
111
112 trckepspl= epspln(1,1)+epspln(2,2)+epspln(3,3)

113 do j=1,3
114 do i=1,3
115 epsplndev(i, j )= epspln(i , j )−third∗trckepspl∗ident2(i, j )
116 end do
117 end do
118
119 !Calculation of deviatoric trial stresses
© Nico Diercks xv
B. Quellcode
120 do j=1,3
121 do i=1,3
122 sdevtr( i , j )= 2.d0∗G∗(epsdev(i,j)−epsplndev(i,j))
123 end do
124 end do
125
126 !Calculation of the invariants of the trial −stress−tensor

127 i1tr = 0.d0 ! First Invariant−Trial
128 do i=1,3
129 i1tr = i1tr +eps(i, i )−epspln(i,i )
130 end do
131 i1tr = 3.d0∗ka∗i1tr
132
133 i2tr = 0.d0 !Second Invariant−Trial,deviatoric

134 do j=1,3
135 do i=1,3
136 i2tr = i2tr +sdevtr(i, j )∗∗2.d0
137 end do
138 end do
139 i2tr = half∗i2tr
140
141 !Calculation of the trial −stress−tensor

142 do j=1,3
143 do i=1,3
144 str ( i , j )= sdevtr(i, j )+third∗i1tr∗ident2(i , j )
145 end do
146 end do
147
148 if (kdomIntPt.eq.intpointout.AND.elemId.eq.elemaus) then

149 ! display of the values for testing purpose
150 write(20,∗)
151 write(20,∗)
152 write(20,∗) "Element: ",elemID
153 write(20,∗) "Loadstep: ",ldstep
154 write(20,∗) "Substep: ", isubst
155 write(20,∗) "Time: ",Time
156 write(20,∗)
157 write(20, ’ (A16,3X $)’) "Lastinkremente:"
158 write(20, ’ (6(EN12.3,2X))’) dStrain(1),dStrain(2),dStrain(3),0.5∗dStrain(4),0.5∗
dStrain(5),0.5∗dStrain(6)
159 write(20,∗)
160 write(20, ’("Ev,n: ",EN14.3)’) evn
© Nico Diercks xvi

B. Quellcode
161 write(20,∗)
162
163 write(20,∗) "Epsilon ges:"

164 do i=1,3
165 write(20, ’ (3(EN14.3))’) eps(i ,1:3)
166 end do
167 write(20,∗)
168
169 write(20,∗) "Epsilon deviatorisch :"

170 do i=1,3
171 write(20, ’ (3(EN14.3))’) epsdev(i ,1:3)
172 end do
173 write(20,∗)
174
175 write(20,∗) "Epsilon plastisch ,n:"

176 do i=1,3
177 write(20, ’ (3(EN14.3))’) epspln(i ,1:3)
178 end do
179 write(20,∗)
180
181 write(20,∗) "Epsilon plastisch , deviatorisch ,n:"

182 do i=1,3
183 write(20, ’ (3(EN14.3))’) epsplndev(i ,1:3)
184 end do
185 write(20,∗)
186
187 write(20,∗) "Trial Spannungen:"

188 do i=1,3
189 write(20, ’ (3(EN14.3))’) str ( i ,1:3)
190 end do
191 write(20,∗)
192
193 write(20,∗) "1. Invariante , trial :", i1tr

194 write(20,∗) "2. Invariante , trial :", i2tr
195 write(20,∗)
196 end if
197
198 !check for yielding

199 R= q∗(1.d0−exp(−b∗evn))+h∗evn
200 y= y0+R
201 phi= y∗∗2.d0−a1∗y0∗i1tr−a2∗i1tr∗∗2.d0
202 grphi= i2tr −third∗phi
© Nico Diercks xvii

B. Quellcode
203 if (grphi>0) then

204 yield= .true.
205 else
206 yield= . false .
207 end if
208
209 if (kdomIntPt.eq.itpointout.AND.elemId.eq.elemaus) then

210 write(20,∗) "Gross−Phi:", grphi
211 write(20,∗)
212 end if
213
214 if ( yield ) then ! plastic deformation

215 ! Iteration of the plastic multiplier
216 lalt = 0.d0 !no meaning of the value, just to start the loop
217 lneu= 1e−7 !start−value
218 n= 0
219

221 write(20,∗) " "
222 write(20, ’ (A5,2(4X,A10),7(3X,A11))’) "N", "Lambda", "Phi", "J", "I1", "kl.
Phi", "dphi/dl", "di/dl", "de/dl", "r[lambda]"
223 end if
224
225 do while(abs( lalt −lneu)>sqtiny)

226 n= n+1
227 lalt = lneu
228 i1= ( i1tr −3.d0∗ka∗lalt∗a1s∗y0)/(1.d0+3.d0∗ka∗lalt∗2.d0∗a2s)
229 t= third∗(a1s∗y0+2.d0∗a2s∗i1)
230 del= 1.d0+2.d0∗G∗lalt
231 ev= evn+(lalt/y0)∗(t∗i1+(2.d0/(del∗∗2.d0))∗i2tr)
232 R= q∗(1.d0−exp(−b∗ev))+h∗ev
233 y= y0+R
234 phi= y∗∗2.d0−a1∗y0∗i1−a2∗i1∗∗2.d0
235
236 ! r [ delta Lambda]

237 rdl= 3.d0∗i2tr−phi∗del∗∗2.d0
238
239 !Ableitung d−y/d−ev

240 ys= q∗b∗exp(−b∗ev)+H
241
242 !Ableitung d−I1/d−Lambda

243 didl= −3.d0∗ka∗((i1∗2.d0∗a2s+a1s∗y0)/(1.d0+3.d0∗ka∗lalt∗2.d0∗a2s))
© Nico Diercks xviii

B. Quellcode
244
245 !Ableitung d−e/d−Lambda

246 dedl= (1.d0/y0)∗(t∗i1+(2.d0/del∗∗2.d0)∗i2tr)+(lalt/y0)∗(2.d0∗(a2s/3.d0)∗
didl∗i1+t∗didl−(8.d0∗G/(del∗∗3.d0))∗i2tr)
247
248 !Ableitung d−phi/d−Lambda

249 dphidl= 2.d0∗y∗ys∗dedl−(a1∗y0+2.d0∗a2∗i1)∗didl
250
251 !J
252 jj = −(del∗∗2∗dphidl+2.d0∗phi∗del∗2.d0∗G)
253
254 !Konvergenzverbesserungen
255 red=1−.95d0∗(real(n)/real(maxit))
256
257 lneu= lalt −red∗(1.d0/jj)∗rdl

258 grphi= i2tr ∗(1.d0/(del∗∗2.d0))−third∗phi
259
260 !Ausgabe der Zwischenschritte

261 if (kdomIntPt.eq.itpointout.and.mod(n,itausgabe).eq.0.and.elemID.eq.elemaus)
then
262 write(20, ’ (I5 ,X,2(EN15.5,X),7(EN13.4,X))’)
263 n, lneu, grphi, jj , i1 , phi, dphidl, didl , dedl, rdl
264 end if
265
266 !abort−criterion, if no convergence is reached

267 if (n>maxit) then
268 keycut= 1
270 write(20,∗) "Keine Konvergenz nach ",n,"Iterationen"
271 write(20,∗)
272 end if
273
274 ev= evn

275 open (unit =20 , file =’ausgabelog.txt’, status =’old’ , action =’write’ ,
position=’append’)
276 write(20,∗)
277 write(20,∗) "Bisection at Element ",ElemID,", Substep ",isubst
278 write( iott ,∗) "Bisection at Element ",ElemID,", Substep ",isubst
279 write(20,∗)
280 write(20,∗) "∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗∗"
281 close (20)
282 goto 100
© Nico Diercks xix

B. Quellcode
283 end if
284 end do
285
286 !Update of the stresses

287 do j=1,3
288 do i=1,3
289 s( i , j )= (1.d0/del)∗sdevtr(i , j )+third∗(i1∗ident2(i , j ))
290 end do
291
292 end do
293
294 !Calculation of the Tangent−operator

295 do i=1,3
296 do j=1,3
297 drde(i , j )= (a1∗y0+2.d0∗a2∗i1−2.d0∗y∗ys∗lneu∗t/y0)∗(del∗∗2.d0∗3.d0
∗ka/(1.d0+3.d0∗ka∗lneu∗2.d0∗a2s))∗ident2(i,j) +(3.d0−4.d0∗y∗ys∗
lneu/(y0∗del∗∗2.d0))∗2.d0∗g∗sdevtr(i,j)
298 dlde( i , j )= (−1.d0/jj)∗drde(i , j )
299 end do
300 end do
301
302 call vzero(C,81)

303 do i=1,3
304 do j=1,3
305 do k=1,3
306 do l=1,3
307 c( i , j ,k, l )= (2.d0∗g/del)∗idev(i, j ,k, l )+(third∗didl∗ident2(i ,
j )−(2.d0∗g/(del∗∗2.d0))∗sdevtr(i,j ))∗dlde(k, l )+(ka/(1.d0
+3.d0∗ka∗lneu∗2.d0∗a2s))∗t1d1(i,j,k,l)
308 end do
309 end do
310 end do
311 end do
312
313 ! calculation of the new plastic strains

314 do j=1,3
315 do i=1,3
316 epsPlt(i , j )= epspln(i , j )+lneu∗(1.d0/del)∗sdevtr(i,j )+lneu∗t∗ident2(i,j )
317 end do
318 end do
319
320 else ! elastic deformation
© Nico Diercks xx
B. Quellcode
321 ev= evn

322 i1= i1tr
323 n= 0.d0
324 lneu= 0.d0
325
326 ! calculation of the tangent−operator

327 call vzero(C,81)
328 do i=1,3
329 do j=1,3
330 do k=1,3
331 do l=1,3
332 c( i , j ,k, l )= 2.d0∗G∗Idev(i,j,k, l )+ka∗t1d1(i,j,k, l )
333 end do
334 end do
335 end do
336 end do
337 end if
338
339 !Display of the results for testing purpose

341 write(20,∗)
342 write(20,∗)"Nach ",n," Iterationen : delta lambda= ", lneu
343 write(20,∗)
344
345 write(20,∗) "Spannungstensor:"

346 do i=1,3
347 if ( yield ) then
348 write(20, ’ (3(EN14.3))’) s( i ,1:3)
349 else
350 write(20, ’ (3(EN14.3))’) str ( i ,1:3)
351 end if
352 end do
353 write(20,∗)
354
355 write(20,∗) "1. Invariante :", i1

356 write(20,∗)
357
358 write(20,∗) "Epsilon plastisch ,n+1:"

359 do i=1,3
361 write(20, ’ (3(EN14.3))’) epsPlt(i ,1:3)
362 else
© Nico Diercks xxi

B. Quellcode
363 write(20, ’ (3(EN14.3))’) epspln(i ,1:3)

364 end if
365 end do
366 write(20,∗)
367
368 write(20,∗) "ev, n+1:"

369 write(20, ’ ((EN14.3))’) ev
370 write(20,∗)
371 end if
372
373 !Update of the hardening user statate variable

374 ustatev(1)=ev
375

377 !Transformation of the stress −tensor to the ansys−vector
378 stress (1)= s (1,1)
379 stress (2)= s (2,2)
380 stress (3)= s (3,3)
381 stress (4)= s (2,1)
382 stress (5)= s (3,2)
383 stress (6)= s (3,1)
384
385 !Transformation of the plastic strain−tensor to the ansys−vector

386 epsPl(1)= epsPlt(1,1)
389 epsPl(4)= 2.d0∗epsPlt(1,2)
392
393 else
394 stress (1)= str (1,1)
395 stress (2)= str (2,2)
396 stress (3)= str (3,3)
397 stress (4)= str (2,1)
398 stress (5)= str (3,2)
399 stress (6)= str (3,1)
400 end if
401
402 !Transformation of the fourth−order tangent−operator

403 dsdePl(1,1)= C(1,1,1,1)
404 dsdePl(2,1)= C(2,2,1,1)
© Nico Diercks xxii

B. Quellcode
405 dsdePl(3,1)= C(3,3,1,1)

406 dsdePl(4,1)= HALF ∗ (C(1,2,1,1) + C(2,1,1,1))
407 dsdePl(5,1)= HALF ∗ (C(2,3,1,1) + C(3,2,1,1))
408 dsdePl(6,1)= HALF ∗ (C(1,3,1,1) + C(3,1,1,1))
409
410 dsdePl(1,2)= C(1,1,2,2)

411 dsdePl(2,2)= C(2,2,2,2)
412 dsdePl(3,2)= C(3,3,2,2)
413 dsdePl(4,2)= HALF ∗ (C(1,2,2,2) + C(2,1,2,2))
414 dsdePl(5,2)= HALF ∗ (C(2,3,2,2) + C(3,2,2,2))
415 dsdePl(6,2)= HALF ∗ (C(1,3,2,2) + C(3,1,2,2))
416
417 dsdePl(1,3)= C(1,1,3,3)

418 dsdePl(2,3)= C(2,2,3,3)
419 dsdePl(3,3)= C(3,3,3,3)
420 dsdePl(4,3)= HALF ∗ (C(1,2,3,3) + C(2,1,3,3))
421 dsdePl(5,3)= HALF ∗ (C(2,3,3,3) + C(3,2,3,3))
422 dsdePl(6,3)= HALF ∗ (C(1,3,3,3) + C(3,1,3,3))
423
424 dsdePl(1,4)= HALF ∗ (C(1,1,1,2) + C(1,1,2,1))

425 dsdePl(2,4)= HALF ∗ (C(2,2,1,2) + C(2,2,2,1))
426 dsdePl(3,4)= HALF ∗ (C(3,3,1,2) + C(3,3,2,1))
427 dsdePl(4,4)= (C(1,2,1,2)+C(1,2,2,1)+C(2,1,1,2)+C(2,1,2,1))/4.d0
428 dsdePl(5,4)= (C(2,3,1,2)+C(2,3,2,1)+C(3,2,1,2)+C(3,2,2,1))/4.d0
429 dsdePl(6,4)= (C(1,3,1,2)+C(1,3,2,1)+C(3,1,1,2)+C(3,1,2,1))/4.d0
430
431 dsdePl(1,5)= HALF ∗ (C(1,1,2,3) + C(1,1,3,2))

432 dsdePl(2,5)= HALF ∗ (C(2,2,2,3) + C(2,2,3,2))
433 dsdePl(3,5)= HALF ∗ (C(3,3,2,3) + C(3,3,3,2))
434 dsdePl(4,5)= (C(1,2,2,3)+C(1,2,3,2)+C(2,1,2,3)+C(2,1,3,2))/4.d0
435 dsdePl(5,5)= (C(2,3,2,3)+C(2,3,3,2)+C(3,2,2,3)+C(3,2,3,2))/4.d0
436 dsdePl(6,5)= (C(1,3,2,3)+C(1,3,3,2)+C(3,1,2,3)+C(3,1,3,2))/4.d0
437
438 dsdePl(1,6)= HALF ∗ (C(1,1,1,3) + C(1,1,3,1))

439 dsdePl(2,6)= HALF ∗ (C(2,2,1,3) + C(2,2,3,1))
440 dsdePl(3,6)= HALF ∗ (C(3,3,1,3) + C(3,3,3,1))
441 dsdePl(4,6)= (C(1,2,1,3)+C(1,2,3,1)+C(2,1,1,3)+C(2,1,3,1))/4.d0
442 dsdePl(5,6)= (C(2,3,1,3)+C(2,3,3,1)+C(3,2,1,3)+C(3,2,3,1))/4.d0
443 dsdePl(6,6)= (C(1,3,1,3)+C(1,3,3,1)+C(3,1,1,3)+C(3,1,3,1))/4.d0
444
445 !Display of the results

© Nico Diercks xxiii

B. Quellcode
447 write(20,∗) "C−Tensor:"

448 do i=1,3
449 do j=1,3
450 write(20, ’ (3(3(EN12.3,2X),4X))’) C(i,1,j ,1:3) , C(i ,2, j ,1:3) ,C(i ,3, j ,1:3)
451 end do
452 write(20,∗)
453 end do
454 write(20,∗)
455
456 write(20,∗) "Tangentenoperator:"

457 do i=1,6
458 write(20, ’ (6(EN12.3,2X))’) dsdePl(i,1:6)
459 end do
460 write(20,∗)
461 write(20,∗) "−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−"
462 write(20,∗)
463 end if
464
465 100 if (elemID.eq.elemaus.and.kdomIntPt.eq.itpointout) close(20)

466
467 return
468 end
© Nico Diercks xxiv

B. Quellcode
B.3. ANSYS-Macro zur

Validierungssimulation
1 fini
2 /clea
3 /VIEW,1,1,2,3
4
5 alpha=0.5
6
7 dk=2 !Dicke d. Klebschicht

8 lr =8 !Länge der Rohre
9 ri =25 !Innenradius
10 ra=30 !Außenrdius
11
12 epsz=0.25 !Längsdehnung bezogen auf Klebschichtdicke

13 epsyz=2∗alpha∗epsz∗(1+epsz) !Gleitung bezogen auf Klebschichtdicke
14
15 rsize =1.25 !Elementgröße in Radialrichtung

16 tsize =2 !Elementgröße in Tangentialrichtung
17 zsize =0.5 !Elementgröße in Klebung in Dickenrichtung
18 lsize =4 !Elementgröße in Rohr in Dickenrichtung
19
20 /prep7
21
22 !_________________________________________
23 ! Erstellung der Geometrie
24
25 cyl4 ,0,0, ri ,0, ra ,180,dk !Klebung

26 cyl4 ,0,0, ri ,180, ra ,360,dk
27
28 cyl4 ,0,0, ri ,0, ra,180,−lr !Rohr

29 cyl4 ,0,0, ri ,180, ra,360,−lr
30
31 wpoffs ,0,0, dk
32 cyl4 ,0,0, ri ,0, ra ,180, lr
33 cyl4 ,0,0, ri ,180, ra ,360, lr
34 wpoffs,0,0,−dk
35
36 csys,1 ! Definition von globalen Polarkoordinaten

37
38 nummrg,kp !Verbindung der Keypoints
© Nico Diercks xxv

B. Quellcode
39
40 !_________________________________________
41 ! Definition der Elementgrößen
42
43 lsel , s , loc ,x,( ri +ra)/2

44 lesize , all , rsize
45
46 lsel , s , loc ,y,45

47 lesize , all ,,,2∗3.1417∗(( ri +ra)/2)/tsize
48
49 lsel , s , loc ,z,dk/2

50 lesize , all , zsize
51
52 lsel , s , loc ,z,−lr/2

53 lsel ,a, loc ,z,dk+lr/2
54 lesize , all , lsize
55
56 !_________________________________________
57 !Meshing
58
59 vsel , s ,volu ,,1,2

60 mat,1
61 et ,1,185
62 keyopt,1,2,0 ! Full Integration
63 vmesh,1,2
64
65 vsel , s ,volu ,,3,6

66 mat,2
67 et ,2,185
68 keyopt,2,2,1 !Uniform Integration
69 vmesh,3,6
70
71 nsel , all
72 nrotat, all !Drehung der lokalen Koordinatensysteme
73
74 !_________________________________________
75 !Material−Definition
76
77 tb,user ,1,1,14
78 tbtemp,1
79
80 tbdata ,1,1588.7,0.34,36.11,8.08,101.1 !Materialparameter 1−5
© Nico Diercks xxvi

B. Quellcode
81 tbdata ,6,33.96,0.643,0.122,0,0.279 !Materialparameter 5−10

82
83 tbdata,11,500,20,90,5 ! optional : Maximale Iterationen, Schrittweite der

Ausgabe, Auszugebendes Element und Integrationspunkt
84
85 tb, state ,1,,1

86
87 mp,ex,2,69e3
88 mp,nuxy,2,0.3
89
90 !_________________________________________
91 !Randbedingungen
92
93 nsel , s , loc ,z,−lr

94 d, all ,uy,0 ! Unterseite , tangentiale Sperrung
95 d, all ,uz,0 ! Unterseite , z−Sperrung
96
97 !_________________________________________
98 !Solution
99
100 /solu
101
102 alls
103
104 DELTIM,0.01,,0.01
105 TIME,1
106
107 AUTOTS,on
108
109 nsel , s , loc ,z,dk+lr+lsize,dk+lr

110 d, all ,uz,epsz∗dk !Belastungsvorgabe (Zug)
111
112 nsel , s , loc ,z,dk+lr

113 nsel , r , loc ,x, ri
114 d, all ,uy,epsyz∗dk∗ri/((ri+ra)/2) !Belastungsvorgabe (Torsion)
115

117 nsel , r , loc ,x ,(( ri +ra)/2)
118 d, all ,uy,epsyz∗dk !Belastungsvorgabe (Torsion)
119

121 nsel , r , loc ,x,ra
© Nico Diercks xxvii

B. Quellcode
122 d, all ,uy,epsyz∗dk∗ra/((ri+ra)/2) !Belastungsvorgabe (Torsion)

123

125 nsel , r , loc ,x, ri +.25∗(ra−ri)
126 d, all ,uy,epsyz∗dk∗(ri+.25∗(ra−ri))/((ri+ra)/2) !Belastungsvorgabe (Torsion)
127

129 nsel , r , loc ,x, ri +.75∗(ra−ri)
130 d, all ,uy,epsyz∗dk∗(ri+.75∗(ra−ri))/((ri+ra)/2) !Belastungsvorgabe (Torsion)
131
132 alls
133
134 outres, all , all

135 outres, svar , all
136
137 NROPT, unsym

138
139 eqslv ,,4

140
141 solve
142
143 !_________________________________________
144 !Postprocessing
145
146 /post1
147 /VIEW,1,1,2,3
148 rsys ,1
149 PLDISP,2
150
151 !_________________________________________
152 !Time−History−Postprocessing
153
154 /POST26
155 !Ermittlung der Verschiebung
156 NSOL,2,147,u,z, Uz_Unten
157 NSOL,3,62,u,z, Uz_Oben
158
159 !Ermittlung der Reaktionskraft

161 ∗get,anzahl,node,0,count
162 ∗get, erster ,node,0,num,min
163
© Nico Diercks xxviii

B. Quellcode
164 rforce ,4, erster , f ,z

165 knoten=erster
166
167 ∗do,i ,1, anzahl−1

168 ∗get,knoten,node,knoten,nxth
169 rforce ,5, knoten,f ,z
170 add ,4,4,5,, zKraft
171 ∗enddo
172
173 !Darstellung der Zugkraft über der Axialverschiebung

174 XVAR,2
175 plvar,4
© Nico Diercks xxix

Beiträge zur Materialtheorie
ISSN 1862-5703
___________________________________________________________________
In dieser Reihe sind bisher folgende Hefte erschienen:
1/07 RETKA, JOHANNES / HÖFER, PHILIPP

Numerische Simulation aushärtender Klebstoffe
2/07 RANZ, THOMAS

Ein feuchte- und temperaturabhängiger anisotroper Werkstoff: Holz
3/07 THURNER, FLORIAN / RANZ, THOMAS

Rheologisches Verhalten von Holz unter Biegebelastung
4/07 LION, ALEXANDER

Einführung in die Lineare Viskoelastizität
5/07 RANZ, THOMAS

Elementare Materialmodelle der Linearen Viskoelastizität im Zeitbereich

BZM Heft 1 09

Hochgeladen von

Copyright:

Verfügbare Formate

BZM Heft 1 09

Hochgeladen von

Dokumentinformationen

Copyright

Verfügbare Formate

Dieses Dokument teilen

Dokument teilen oder einbetten

Freigabeoptionen

Stufen Sie dieses Dokument als nützlich ein?

Sind diese Inhalte unangemessen?

Copyright:

Verfügbare Formate

BZM Heft 1 09

Hochgeladen von

Copyright:

Verfügbare Formate

BEITRÄGE

Herausgeber: Univ.-Prof. Dr.-Ing. habil. Alexander Lion

Dipl.-Ing. Nico Diercks

Validierung und ANSYS-

Heft 01/09 ISSN 1862-5703

Dipl.-Ing. Nico Diercks

Universität der Bundeswehr München

Heft 01/09 der Reihe

Beiträge zur Materialtheorie

Herausgeber: Univ.-Prof. Dr.-Ing. habil. Alexander Lion

Univ.-Prof. Dr.-Ing. habil. Alexander Lion

Dipl.-Ing. Nico Diercks:

Beiträge zur Materialtheorie 01/09

Postanschrift des Verfassers:

2.2. Grundlagen der Finite-Elemente-Methode . . . . . . . . . . . . 17

3.2.3.1. Aufruf der Routine . . . . . . . . . . . . . . . . 45

5. Fazit und kritische Bewertung 58

© Nico Diercks III

2.1. σ(ε) eines Zugversuches . . . . . . . . . . . . . . . . . . . . . . 4

4.1. Zugprobe nach EN ISO 527-2 . . . . . . . . . . . . . . . . . . . 49

Bei modernen Fertigungsprozessen, sowohl in der Automobilindustrie als auch

Um diese Zuverlässigkeit zu erzielen, sind bei den Klebeverbindungen, zum

1.2. Ziel der Arbeit

Im Rahmen dieser Arbeit wird ein Materialmodell der Elastoplastizität in die

1.3. Aufbau der Arbeit

Der Hauptteil dieser Arbeit gliedert sich in 3 Abschnitte:

3. Das Schlimmer-Mahnken Materialmodell

2.1.2. eindimensionale Plastizitätstheorie

Um die Effekte und die Herangehensweise an selbige zu verdeutlichen, wird

Abbildung 2.1.: σ(ε) eines Zugversuches

trachtet wird. In Abbildung 2.1 ist eine solche Spannungs-Dehnungskurve skiz-

1. Entlastung: dσ < 0, hier gilt das elastische Materialgesetz (Hooke’sches

2. Belastungssteigerung: dσ > 0, wobei weitere Deformationen auftreten

Zur mathematischen Formulierung der Fallunterscheidung wird eine Fließfunk-

2.1.2.2. Beispiel: „Prandtl-Reuss-Plastizität", eindimensional

Abbildung 2.2.: σ(ε) nach Prandtl-Reuss

zugehörige Fließfunktion lautet

Wird σ weiter gesteigert bis |σ| = σF , so folgt F (σ) = 0. Auch wenn in

2.1.2.3. Elasto-Plastizität mit linearer Verfestigung

Nach der relativ einfachen „Prandtl-Reuss-Plastizitätfolgt nun die Betrachtung

Abbildung 2.3.: σ(ε) mit Verfestigung

F (σ, εpl ) = 0 ⇒ (σ − Bεpl )2 = σF2

Nach Bildung des Differentials und Umstellen nach dεpl

mithilfe des Materialgesetzes der Elastizität folgt

Nach Integration erhält man

Die Bestimmung der Integrationskonstante C erfolgt durch

2.1.2.4. Zusammenhang zwischen den Moduln

In den vorherigen Kapiteln wurde das Materialverhalten mithilfe folgender

3. dem plastischen Tangentenmodul

Abbildung 2.4.: infinitesimales Stück der σ(ε) -Kurve

Zur Bestimmung der Konstanten B wird Glg. 2.4 erneut betrachtet:

Weiterhin lautet Glg. 2.5

Aus den Glgn. 2.4 und 2.8 folgt

wobei man durch Koeffizientenvergleich sofort erkennt

2.1.3. Verallgemeinerung auf 3 Dimensionen am Beispiel

σ31 σ32 σ33 ε31 ε32 ε33

Wegen der Symmetrie der Tensoren

6. bei Verfestigung: Evolutionsgleichungen für die Spannungen

2.1.3.2. Zerlegung in deviatorische und volumetrische Anteile

Zunächst werden Spannungs- und Verzerrungstensor in die deviatorischen und