1. No category

Indholdsfortegnelse
INDHOLDSFORTEGNELSE
DEL I
FORSØG ................................................................................... 3
A
Elastiske konstanter...........................................................................................5
A.1
A.2
A.3
A.4
B
Dataopsamling ..............................................................................................5
Brudstyrkemåling på massivt aluminiumsemne ...........................................5
Elasticitetsmodul og Poissons forhold for massivt materiale .......................8
Elasticitetsmodul for porøst materiale ........................................................15
Forsøg med massiv og porøs konsol ...............................................................21
B.1 Dataopsamling ............................................................................................21
B.2 Databehandling ...........................................................................................23
B.3 Resultater ....................................................................................................26
DEL II
ANALYTISKE MODELLER ..........................................................31
C
Simpel bjælkemodel........................................................................................33
C.1
C.2
C.3
C.4
C.5
C.6
C.7
C.8
D
Airys spændingsfunktion ................................................................................41
DEL III
E
Tværsnitsdata ..............................................................................................33
Virkelige snitkræfter ...................................................................................34
Virkelige tøjninger ......................................................................................35
Virtuelle snitkræfter ....................................................................................35
Virtuelle Kræfters Princip ...........................................................................36
Udbøjning....................................................................................................37
Kontrol af udbøjning ...................................................................................37
Spændinger..................................................................................................38
NUMERISKE MODELLER ...........................................................43
Bestemmelse af materialeparametre ved beregninger.....................................45
E.1 Elasticitetsmodul og Poissons forhold ........................................................45
E.2 Forskydningsmodul.....................................................................................48
F
Elementmetodeteori .........................................................................................51
F.1
F.2
F.3
F.4
G
CST- og LST-elementopbygning....................................................................57
G.1
G.2
G.3
G.4
G.5
H
Potentiel energi ...........................................................................................51
Variation af potentiel energi .......................................................................53
Konvergens af potentiel energi ...................................................................55
Konvergens af lastens arbejede...................................................................56
Opbygning af CST-elementer .....................................................................57
Opstilling af materialestivhedsmatricen......................................................64
Opbygning af LST-model ...........................................................................66
Lastfordeling på knuder ..............................................................................70
Spændinger i konsol beregnet ud fra CST- og LST-model ........................73
Eliminering af indre frihedsgrader ..................................................................75
1
Indholdsfortegnelse
H.1
H.2
H.3
H.4
H.5
I
Topologioptimering ......................................................................................... 81
I.1
I.2
I.3
I.4
I.5
I.6
I.7
J
Oprindelig målfunktion............................................................................... 81
Sidebetingelser ............................................................................................ 82
Endelig målfunktion.................................................................................... 82
Optimeringsbetingelse ................................................................................ 83
Oprindelig følsomhed ................................................................................. 83
Filtreret følsomhed...................................................................................... 83
Opdateringsrutine........................................................................................ 85
Singularitet af spændinger ............................................................................... 87
J.1
J.2
J.3
J.4
2
Potentiel energi ........................................................................................... 75
Adskillelse af indre og ydre frihedsgrader.................................................. 76
Variation med hensyn til indre frihedsgrader ............................................. 77
Potentiel energi som funktion af ydre frihedsgrader................................... 78
Bestemmelse af indre frihedsgraders flytninger ......................................... 80
Problematik ................................................................................................. 87
Afhængigheder............................................................................................ 88
Spændingsvariation..................................................................................... 88
Singularitet for konsol................................................................................. 89
Del I
Forsøg
A Elastiske konstanter
A
Del I Forsøg
ELASTISKE KONSTANTER
Her er gennemgået databehandlingen til forsøgene, udført til bestemmelse af elastiske konstanter for
aluminiumsmaterialet. Princippet for dataopsamlingen i forsøgene er beskrevet generelt for alle forsøgene, mens databehandlingen på de enkelte forsøg er beskrevet for hver af de forskellige forsøg.
A.1
DATAOPSAMLING
Under forsøgene er der løbende opsamlet data om kraft, tøjning og flytning. Dataene bliver lagret i
datafiler til videre bearbejdelse, og der vil i de enkelte forsøgsafsnit blive henvist til mapper på den
vedlagte cd-rom, hvor de respektive data kan findes.
Den påførte kraft i forsøgene registreres som et analogt spændingsfald, konverteres til digitalt signal
og opsamles i Volt, hvor 10 V = 25 kN.
Tøjninger er registreret vha. strain gauges. Disse registrerer en ændring i modstand, som konverteres
til et digitalt signal. Outputtet fra strain gauge målingerne er [μm/m].
Flytninger er målt med flytningsmålere, og outputtet på disse data er [mm]
Alle data er opsamlet med en frekvens på 5 Hz.
A.2
BRUDSTYRKEMÅLING PÅ MASSIVT ALUMINIUMSEMNE
I det følgende er databehandlingen til brudforsøget på aluminium beskrevet. Forsøgsbeskrivelsen af
trækbrudforsøget kan findes i hovedrapport afsnit 1.1. I databehandlingen vil der først blive gennemgået principperne for databehandlingen og de anvendte formler, dernæst følger optegning af
grafer og beregning af dertilhørende resultater.
Dokumentation på data og databehandling findes på den vedlagte cd-rom i følgende mapper:
•
•
Elastiske konstanter forsøg\Brudforsøg 1-2
Elastiske konstanter forsøg \Brudforsøg 3
A.2.1
Databehandling
De målte trækkræfter og flytninger på forsøgsemnerne afbildes med flytningen som ordinat, hvorved
materialets arbejdskurve fremkommer. Dog kan arbejdskurven ikke bruges til at fastlægge elasticitetsmodulen, da dette kræver, at ordinaten angiver den relative flytning i aksialretningen, det vil sige
5
Del I Forsøg
A Elastiske konstanter
tøjningen. Til det aktuelle formål er den faktiske flytning som ordinat dog tilstrækkelig til at konstatere, hvornår materialet overgår til flydning.
Ud fra arbejdskurven aflæses kraften ved begyndende flydning, Fy , hvor arbejdskurven overgår fra
lineær til krum. Kraften ved brud, Fu , aflæses som maksimalværdien af F på arbejdskurven. De
tilhørende flyde- og brudspændinger beregnes herefter ved (A.1), idet emnerne kun udsættes for
normalkraft.
fy =
Fy
A
, fu =
Fu
A
(A.1)
hvor
fy
Fy
fu
Fu
A
er flydespænding [MPa]
er normalkraft ved flydning [N]
er brudspænding [MPa]
er normalkraft ved brud [N]
er tværsnitsarealet [mm2]
I figur 1, figur 2 og figur 3 er vist de fremkomne arbejdskurver for hver af de tre belastningsramper.
Desuden er angivet de aflæste værdier for flyde- og brudkræfter med tilhørende beregnede spændinger. Præcisionen af disse aflæsninger er ikke afgørende, idet formålet med trækbrudforsøget som
nævnt kun tjener som kontrol for gyldigheden af andre undersøgelser i projektet. Bemærk, at figur 2
og figur 3 kun er udtryk for materialets arbejdskurve til og med brudpunktet. Herefter fortsætter
kurven, da måleudstyret fortsat er aktivt efter brud, men ud fra et materialeanalytisk synspunkt er
kurven uinteressant efter brudstadiet.
Forsøg 1: Forsøgsemne 1 belastet til flydning.
Under forsøget blev der gennemført en belastningsrampe med en flytning på 0 − 14 mm over 200
sekunder.
12
10
F [kN]
8
6
4
2
0
0
5
u [mm]
Figur 1: Arbejdskurve for delforsøg 1.
Følgende værdier er aflæst på figur 1:
6
10
15
A Elastiske konstanter
Del I Forsøg
Fy = 8400 N
fy =
8400 N
49,19 mm 2
= 171MPa
Forsøg 2: Forsøgsemne 1 belastet til brud.
Under forsøget blev der gennemført en belastningsrampe med en flytning på 0 − 20 mm over 400
sekunder.
12
10
F [kN]
8
6
4
2
0
0
5
u [mm]
10
15
Figur 2: Arbejdskurve for delforsøg 2.
Følgende værdier er aflæst på figur 2:
Fu = 10650 N
fu =
10650 N
49,19 mm 2
= 216 MPa
Forsøg 3: Forsøgsemne 2 belastet til brud.
Under forsøget blev der gennemført en belastningsrampe med en flytning på 0 − 25 mm over 400
sekunder.
7
Del I Forsøg
A Elastiske konstanter
12
10
F [kN]
8
6
4
2
0
0
5
u [mm]
10
15
Figur 3: Arbejdskurve for delforsøg 3.
Følgende værdier er aflæst på figur 3:
Fy = 8900 N
fy =
8900 N
49,19 mm 2
= 181MPa
Fu = 10700 N
fu =
10700 N
49,19 mm 2
= 217 MPa
A.3
ELASTICITETSMODUL OG POISSONS FORHOLD FOR MASSIVT
MATERIALE
Dette er databehandlingen tilhørende forsøget til bestemmelse af elasticitetsmodul og Poissons forhold for det massive aluminiumsmateriale. Forsøgsbeskrivelsen til forsøget kan findes i hovedrapport afsnit 1.2. I databehandlingen vil der først blive gennemgået principperne for databehandlingen
og de anvendte formler, dernæst følger optegning af grafer og beregning af resultater.
Dokumentation på data og databehandling findes på den vedlagte cd-rom i følgende mapper:
•
•
•
8
Elastiske konstanter forsøg \Elastiske konstanter massiv forsøg 1
Elastiske konstanter forsøg \Elastiske konstanter massiv forsøg 2
Elastiske konstanter forsøg \Elastiske konstanter massiv forsøg 3
A Elastiske konstanter
A.3.1
Del I Forsøg
Databehandling
Ud fra data på kraften findes spændingen i tværsnittet af forsøgsemnet ved (A.1). Da tøjningen i
længderetningen måles, kan arbejdskurven for materialet optegnes. E er afbildet som hældningen på
arbejdskurven, og udregnes for de enkelte målinger som
E=
σ
ε længde
(A.2)
Fra (A.2) optegnes E som funktion af ε længde . Ud fra denne graf sorteres ”dårlige data” fra for at få
en mere præcis værdi for E. De dårlige data omfatter de målte værdier for de små spændinger og
tøjninger i starten og slutningen af belastningsrampen, det vil sige. ved den initielle belastning og
den afsluttende aflastning af emnet. For disse målinger er der fundet store udsving for E.
Poissons forhold ν findes som forholdet mellem tværkontraktionen ε tvær og længdetøjningen ε længde .
Denne udregnes også for de enkelte målinger som
ν=
ε tvær
ε længde
(A.3)
Ud fra (A.3) optegnes en graf med ν som funktion af ε længde .
Dernæst udregnes middelværdien for både E og ν som
x=
x1 + x2 + ...xn
n
(A.4)
hvor
x
er middelværdien
x
er de uafhængige observationer
n
er antal observationer
[Teknisk Ståbi 2004, p 42]
For at kunne vurdere fejlen på de beregnede middelværdier beregnes spredningen s ved
s=
( x1 − x )2 + ( x2 + x ) 2 + ...( xn − x )2
n −1
(A.5)
[Teknisk Ståbi 2004, p 42]
For at få den procentvise fejl beregnes variationskoefficienten δ ved
δ=
s
x
(A.6)
[Teknisk Ståbi 2004, p 42]
I det følgende er ovenstående databehandling gennemført for de enkelte forsøg. Det skal bemærkes
at de ”dårlige data” er frasorteret i figurerne.
9
Del I Forsøg
A Elastiske konstanter
Forsøg 1: Forsøgsemne 3 - belastning på langs af valseretningen
Forsøget er gennemført med en belastningsrampe med kræfter på 0-18 000 N og 18 000-0 N på 200
s.
80
70
60
σ [MPa]
50
40
30
20
10
0
0
0.2
0.4
0.6
εlængde [-]
0.8
1
1.2
-3
x 10
0.8
1
1.2
-3
x 10
Figur 4: Arbejdskurve for forsøg 1.
4
7.02
x 10
E [MPa]
7
6.98
6.96
6.94
6.92
0
0.2
0.4
0.6
εlængde [-]
Figur 5: E som funktion af ε længde for forsøg 1.
10
A Elastiske konstanter
Del I Forsøg
-0.33
-0.332
ν [-]
-0.334
-0.336
-0.338
-0.34
0
0.2
0.4
0.6
0.8
εlængde [-]
1
1.2
-3
x 10
Figur 6: v som funktion af ε længde for forsøg 1.
Forsøg 2: Forsøgsemne 4 - belastning på langs af valseretningen
Forsøget er gennemført med en belastningsrampe med kræfter på 0-6000 N på 100 s.
25
σ [MPa]
20
15
10
5
0
0
2
εlængde [-]
4
x 10
-4
Figur 7: Arbejdskurve for forsøg 2.
11
Del I Forsøg
A Elastiske konstanter
4
6.85
x 10
E [MPa]
6.84
6.83
6.82
6.81
6.8
0
2
εlængde [-]
4
-4
x 10
Figur 8: E som funktion af ε længde for forsøg 2.
-0.336
-0.338
ν [-]
-0.34
-0.342
-0.344
-0.346
-0.348
0
2
εlængde [-]
4
-4
x 10
Figur 9: v som funktion af ε længde for forsøg 2.
Forsøg 3: Forsøgsemne 5 - belastning på tværs af valseretningen
Forsøget er gennemført med en belastningsrampe med kræfter på 0-18 000 N og 18 000-0 N på 200
s.
12
A Elastiske konstanter
Del I Forsøg
80
70
σ [MPa]
60
50
40
30
20
10
0
0
0.2
0.4
0.6
εlængde [-]
0.8
1
1.2
-3
x 10
0.8
1
1.2
-3
x 10
Figur 10: Arbejdskurve for forsøg 3.
4
7.26
x 10
7.24
E [MPa]
7.22
7.2
7.18
7.16
7.14
0
0.2
0.4
0.6
εlængde [-]
Figur 11: E som funktion af ε længde for forsøg 3.
13
Del I Forsøg
A Elastiske konstanter
-0.346
-0.348
ν [-]
-0.35
-0.352
-0.354
-0.356
-0.358
0
0.2
0.4
0.6
0.8
εlængde [-]
1
1.2
-3
x 10
Figur 12: v som funktion af ε længde for forsøg 3.
For alle de ovenstående figurer er de ”dårlige data” sorteret fra.
Arbejdskurverne på figur 4, figur 7, figur 10 viser at materialet opfører sig lineært elastisk. Det at
arbejdskurverne forbliver rette linier under hele forløbet indikerer, at forsøgsemnerne ikke er blevet
belastet til flydning.
Graferne for E på figur 5, figur 8, figur 11 viser at E varierer en lille smule for de enkelte målinger.
Variationen er dog ikke stor efter de ”dårlige data” er sorteret fra, og det vurderes derfor rimeligt at
regne E for værende en konstant, beregnet som middelværdien af de enkelte målinger. Det kan ses at
målingerne på E har lidt lavere værdier under afbelastning, hvilket skyldes fænomenet hysterese.
Dette fænomen er behandlet under diskussion af forsøgsresultater i hovedrapport afsnit 1.2.4.
Graferne for ν på figur 6, figur 9 og figur 12 viser at også denne værdi varierer en lille smule for de
enkelte målinger. Variationen på denne er minimal efter at de ”dårlige data” er frasorteret, så også
her vurderes det rimelig at regne den for værende en konstant værdi, udregnet som middelværdien af
de enkelte målinger.
Tabel 1 viser middelværdien for E for de enkelte forsøg, samt en spredning og en variationskoefficient. Spredningen og variationskoefficienten for E viser at usikkerheden ved at regne E som konstant
er minimal.
Tabel 1: Middelværdier for E, samt
spredning og variationskoefficient.
Forsøg
1
2
3
E [MPa]
s [MPa]
δ [-]
69 589
68 181
71 872
226
63
276
0.003
0.001
0.004
Tabel 2 viser middelværdien for υ for de enkelte forsøg, samt en spredning og en variationskoefficient. Også her viser spredningen og variationskoefficienten, at usikkerheden ved at regne υ som konstant er minimal.
14
A Elastiske konstanter
Del I Forsøg
Tabel 2: Middelværdier for υ,
samt spredning og variationskoefficient.
A.4
Forsøg
υ [-]
s [-]
δ [-]
1
2
3
0.336
0.344
0.351
0.001
0.002
0.002
0.004
0.005
0.005
ELASTICITETSMODUL FOR PORØST MATERIALE
Dette er databehandlingen tilhørende forsøget til bestemmelse af et ækvivalent elasticitetsmodul for
det porøse aluminiumsmateriale. Forsøgsbeskrivelsen til forsøget kan findes i hovedrapport afsnit
1.3. I databehandlingen vil der først blive gennemgået principperne for databehandlingen og de anvendte formler, dernæst følger optegning af grafer og beregning af resultater.
Dokumentation på data og databehandling findes på den vedlagte cd-rom i følgende mapper:
•
•
•
Elastiske konstanter forsøg \Elasticitetsmodul porøs forsøg 1
Elastiske konstanter forsøg \Elasticitetsmodul porøs forsøg 2
Elastiske konstanter forsøg \Elasticitetsmodul porøs forsøg 3
A.4.1
Databehandling
For den påførte belastningsrampe udregnes en ækvivalent spænding ved (A.1), hvor tværsnitsarealet
er indsat, som var materialet massivt. Tøjningen kan bestemmes ud fra flytningsmåleren, som måler
et gennemsnit for længdetøjningerne mellem to målepunkter. Tøjningen kan bestemmes som
ε længde =
ΔL L − L0
=
L0
L0
(A.7)
hvor
L0
er afstanden mellem de to målepunkter i udeformeret tilstand
L
er afstanden mellem de to målepunkter i deformeret tilstand
Da der er målt en sammenhæng mellem belastning og tøjninger, kan en ækvivalent arbejdskurve
optegnes for det porøse materiale. Elasticitetsmodulen for de enkelte målinger udregnes ved (A.2),
og E optegnes som funktion af ε længde . Ud fra denne graf sorteres ”dårlige data” fra for at få en mere
præcis værdi for E. Da de målte værdier ved aflastningen gav højere værdier for E end ved belastning, er der her valgt at se bort fra disse data. Dette er gjort ud fra betragtningen om, at målingerne
ved aflastning burde have de samme eller lavere værdier for stivheden på grund af hysterese. Der må
altså være fejl på målingerne ved aflastning af forsøgsemnet. Data for små tøjninger er ligeledes
sorteret fra, da målingerne her viser store udsving.
Det ækvivalente elasticitetsmodul findes som middelværdien af de beregnede hældninger på den
ækvivalente arbejdskurve. Middelværdien udregnes som (A.4), spredningen udregnes som (A.5), og
variationskoefficienten udregnes som (A.6).
I det følgende er ovenstående databehandling gennemført for de enkelte forsøg. Det skal bemærkes
at de ”dårlige data” er frasorteret i figurerne.
15
Del I Forsøg
A Elastiske konstanter
Forsøg 1: Forsøgsemne 6
Forsøget er gennemført med en belastningsrampe med kræfter på 0-3000 N og 3000-0 N på 200 s.
12
10
σ [MPa]
8
6
4
2
0
0
2
εlængde [-]
4
-4
x 10
Figur 13: Arbejdskurve for forsøg 1.
4
3
x 10
2.95
E [MPa]
2.9
2.85
2.8
2.75
2.7
0
2
εlængde [-]
4
-4
x 10
Figur 14: Eækv som funktion af ε længde for forsøg 1.
Forsøg 2: Forsøgsemne 6
Forsøget er gennemført med en belastningsrampe med kræfter på 0-6000 N og 6000-0 N på 200 s.
16
A Elastiske konstanter
Del I Forsøg
25
σ [MPa]
20
15
10
5
0
0
2
4
εlængde [-]
6
8
x 10
-4
Figur 15: Arbejdskurve for forsøg 2.
3
x 10
4
E [MPa]
2.95
2.9
2.85
2.8
0
2
4
εlængde [-]
6
8
-4
x 10
Figur 16: Eækv som funktion af ε længde for forsøg 2.
Forsøg 3: Forsøgsemne 6
Forsøget er gennemført med en belastningsrampe med kræfter på 0-9000 N og 9000-0 N på 200 s.
17
Del I Forsøg
A Elastiske konstanter
35
30
σ [MPa]
25
20
15
10
5
0
0
0.2
0.4
0.6
εlængde [-]
0.8
1
1.2
-3
x 10
0.8
1
1.2
-3
x 10
Figur 17: Arbejdskurve for forsøg 3.
4
2.94
x 10
E [MPa]
2.92
2.9
2.88
2.86
2.84
0
0.2
0.4
0.6
εlængde [-]
Figur 18: Eækv som funktion af ε længde for forsøg 3.
Arbejdskurverne på figur 13, figur 15 og figur 17 viser, at det porøse materiale opfører sig lineært
elastisk. Det at arbejdskurverne forbliver rette linier under hele forløbet indikerer, at forsøgsemnerne
ikke er blevet belastet til flydning.
Graferne for E på figur 14, figur 16 og figur 18 viser at Eækv varierer en lille smule for de enkelte
målinger. Variationen er dog ikke stor efter de ”dårlige data” er sorteret fra, og det vurderes derfor
rimeligt at regne Eækv for værende en konstant, beregnet som middelværdien af de enkelte målinger.
Tabel 3 viser middelværdien for Eækv for de enkelte forsøg, samt en spredning og en variationskoefficient. Spredningen og variationskoefficienten for Eækv viser at usikkerheden ved at regne Eækv som
konstant er minimal.
18
A Elastiske konstanter
Del I Forsøg
Tabel 3: Middelværdier for Eækv , samt
spredningen s og variationskoefficienten δ.
Forsøg
1
2
3
E [MPa]
s [MPa]
δ [-]
28 978
29 014
28 844
210
134
105
0.007
0.005
0.004
19
B Forsøg med massiv og porøs konsol
B
B.1
Del I Forsøg
FORSØG MED MASSIV OG PORØS
KONSOL
DATAOPSAMLING
Under forsøgene er der løbende opsamlet data om kraft, tøjning og flytning. Dataene bliver lagret i
datafiler til videre bearbejdelse. Der videre behandling af dataene er udført i MatLab og er vedlagt
på cd-rom i følgende mapper:
•
•
Forsøg\Massiv
Forsøg\Porøs
Den påførte kraft registreres som et analogt spændingsfald, og konverteres til digital signal og opsamles i volt. For at montere konsollen var en forbelastning nødvendig. Denne får dog ingen indflydelse på resultaterne af kraft, tøjninger og nedbøjninger, da disse alle blev nulstillet efter forbelastningen var påført.
Tøjninger er registreret vha. strain gauges. Der er i forsøget anvendt tre enkeltgauges med en gaugefaktor kg=2.15 og fire rosettegauges med en gaugefaktor på kg=2.12. Gauge 1 - 12 anvendes til at
registrere tøjninger til brug til sammenligning med de numeriske og analytiske undersøgelser, og
gauge 18 – 20 anvendes til at kontrol. På figur 19 ses placeringen af flytnings- og tøjningsmålere for
det massive forsøg. For den massive og porøse konsol er der syv understøtningspositioner. For hver
understøtningsposition gennemføres et delforsøg, hvor forsøgsemnet henholdsvis belastes og aflastes, mens data opsamles. Til sidst kan de endelige resultater findes ved at superponere resultaterne
fra de 7 delforsøg. Dette er tilladt, så længe spændingerne holdes inden for det lineære elastiske
område, og så længe deformationerne er små.
21
Del I Forsøg
B Forsøg med massiv og porøs konsol
Forside
1
2
2
3
1
10
5
11
12
55.0
20.0
20.0
13.3°
8
20.0
4
9
2 0.0
4
7
28.1
6
2 8 .1
28.1
18
13.3°
3
135.0
210.0
15.0
Bagside
1
9
20.0
13.3°
28.1
13.3°
135.0
5
20
55.0
7
6
45.0
45.0
8
45.0
60.0
435.0
Figur 19: Placering af flytningsmålere og strain gauges. Tal i cirkler refererer til strain gauges, og tal i firkanter refererer til flytningsmålere. Mål i mm.
Alle gauges er monteret på halvbroer, som vist på figur 20, men da der kun er brug for en kvartbro,
er der placeret en konstant modstand på den ene af de aktive modstande i Wheatstone broen. Dermed
er ε13 = 0 , og den bliver til en kvartbro.
Figur 20: Wheatstone bro med to aktive modstande (halvbro).
Outputtet V0 konverteres fra et analog signal til digital signal, og omregnes til en tøjning inden lagring i data-filen. Sammenhængen mellem tøjning og ændring i spændingsfald for halvbroen er givet
ved
ΔV0 =
hvor
22
Kg
4
(ε12 − ε13 )Ve
(B.1)
B Forsøg med massiv og porøs konsol
Del I Forsøg
ΔV0
er ændringen i spændingsfaldet over 1 – 4 [V]
Kg
er gaugefaktoren for strain gaugen [-]
Ve
er et konstant spændingsfald fra 2 – 3 [V]
ε12 , ε13 er målte tøjninger [-]
[Hansen 1998]
Idet ε13 = 0, isoleres ε12 i (B.1) til:
ε12 =
K g Ve
4 ΔV0
(B.2)
Flytninger er registreret ved hjælp af induktive flytningsmålere af fabrikat HF Jensen, og data opsamles som ændring i spændingsfald og konverteres til en flytning, idet der er en lineær sammenhæng mellem ændring i spændingsfald og flytning. Flytningsmålerne kan måle i intervallet ± 1 mm
[Hansen 1998, notat 6].
B.2
DATABEHANDLING
I det følgende er håndteringen af dataene beskrevet. I datafilerne er registreret sammenhørende værdier af tid, kraft, flytninger og tøjninger. Tiden vil ikke blive beskrevet i det følgende, da denne ikke
har nogen relevans. Der er ikke skelnet mellem det massive og det porøse forsøgsemne.
B.2.1
Kraft
Som beskrevet tidligere bliver den påførte kraft lagret i datafilerne som en ændring i spænding i
Volt. Sammenhængen mellem ændring i spændingsfald ΔU og kraft F er angivet i (B.3).
F = ΔU ⋅
B.2.2
6 kN
10 V
(B.3)
Flytninger
Flytningerne er registreret i de 8 flytningsmålere for hvert af de 7 delforsøg. Ved at betragte delforsøg 1, hvor understøtningerne står i position 1 for den massive konsol, ses det jf. figur 21, at belastnings- og aflastningskurven ikke er sammenfaldende som forventet.
23
Del I Forsøg
B Forsøg med massiv og porøs konsol
5
Kraft [kN]
4
3
2
1
0
0
0.1
0.2
0.3
Flytning [mm]
0.4
0.5
Figur 21: Belastnings- og aflastningskurve for flytningsmåler 1 for placering 1.
Dette skyldes fænomenet hysterese, hvor der sker dissipation af tøjningsenergien i form af afgivelse
af varme. Da dette energitab ikke ønskes medtaget, vælges der at se bort fra aflastningskurven for
både det massive og det porøse forsøgsemne, og dermed vil kun belastningskurven indgå i den endelige beregning af flytningen.
Flytningsmålerne er upræcise, og det er derfor nødvendigt at kalibrere disse for at opnå en så præcis
flytning som muligt. På figur 22 ses kalibreringsgrafen for flytningsmåler 1. Kalibreringen er udført
for hver af de anvendte flytningsmålere ved at benytte et måleapparat med større præcision end flytningsmålerne. Sammenhørende værdier af de målte værdier og de eksakte værdier er fundet, og ved
hjælp af regression er der bestemt en lineær sammenhæng. Denne er efterfølgende brugt til at korrigere de målte værdier, og dermed er en bedre nøjagtighed opnået.
Kalibrerede flytninger [mm]
2
1.5
1
0.5
0
0
0.5
1
1.5
Målte flytninger [mm]
2
Figur 22: Lineær kalibreringskurve for flytningsmåler 1 samt
kalibreringspunkter.
Som det ses på figur 19, sidder de 8 flytningsmålere parvis to og to. Ved databehandlingen kunne det
konstateres, at de parvise flytningsmålere ikke har vist den samme flytning efter korrigeringen, hvilket ses på figur 23, hvor flytningsmåler 1 og 5 er plottet for delforsøg 1 med understøtninger i position 1.
24
B Forsøg med massiv og porøs konsol
Del I Forsøg
5
Kraft [kN]
4
3
2
1
0
0
0.1
0.2
0.3
Flytning [mm]
0.4
0.5
Figur 23: Belastningskurve for flytningsmåler 1 og 5 for
delforsøg 1.
Dette kan skyldes, at forsøgsemnet ikke har stået helt lodret, og derfor er det valgt at beregne den
resulterende flytning som et gennemsnit af de parvise flytningsmålere.
På figur 23 ses det, at hverken kurverne for flytningsmåler 1 og 5 ligger på en ret linie som forventet.
Grunden til dette er målingsunøjagtigheder, hvilket medfører, at der ikke er linearitet. For at kunne
superponere alle delforsøg og i øvrigt kompensere for de nævnte måleunøjagtigheder, er der udført
lineær regression på måleresultaterne for flytningerne.
Hermed haves flytningsgrafer for de 4 målepunkter for hvert af de 7 delforsøg.
B.2.3
Tøjninger
For at måle tøjninger blev der monteret 4 rosette gauges hvis placering ses på figur 19, for den massive konsol. På figur 24 ses det, at der som for tøjningskurverne også optræder hysterese, hvorfor der
igen kun betragtes belastningskurven.
Belastnings- og aflastningskurve mht. tøjning
5
Kraft [kN]
4
3
2
1
0
0
20
40
60
Tøjning [μm/m]
80
100
Figur 24: Tøjningskurve for gauge 1.
Denne ligger heller ikke på en ret linie som forventet, hvorfor der ligeledes foretages en lineær regression, og dermed haves tøjningsgrafer for alle gauges for hvert af de 7 delforsøg.
25
Del I Forsøg
B Forsøg med massiv og porøs konsol
B.3
RESULTATER
For at opnå de samlede flytninger og tøjninger for konsollen, benyttes superpositionsprincippet, hvor
resultaterne fra de syv belastningspositioner adderes. I det følgende vil resultaterne for den massive
og dernæst den porøse konsol blive angivet. Der vil for den massive konsol blive angivet flytninger
og tøjningerne for en fladelast på 4.5 MPa. For den porøse konsol vil der blive angivet de fire flytninger for en fladelast på 1.5 MPa. Ved at vælge fladelasten til ovenstående for henholdsvis den
massive og porøse konsol, er der sikret at denne ikke flyder.
B.3.1
Massiv konsol
For en fladelast på 4.5 MPa bliver de samlede flytninger og tøjninger angivet i det følgende, hvor
flytningerne ses i tabel 4.
Tabel 4:Flytninger for den massive konsol. Placering af flytningsmålere er angivet på figur 19.
Flytningsmåler
1 og 5
2 og 6
3 og 7
4 og 8
Samlede flytning [mm]
0.2900
0.1633
0.0896
0.0062
Da flytningerne ønskes beskrevet for konsollen, korrigeres de således, at flytningen ved indspændingen er u2 = 0. Det antages her, at ændringen i flytning fra fastholdelsen til flytningsmåler 1-5 er
Δu2 ≈ 0, og dermed kan flytningerne angivet i tabel 5 med hensyn til x1 findes.
Tabel 5: Flytninger mht. x1. Koordinatsystem som
figur 25.
x1 [mm]
0
15
60
105
150
u2 [mm]
0
0.0000
0.1267
0.2004
0.2838
På figur 25 ses placeringen af de fire rosettegauges samt det globale koordinatsystem. I det følgende
vil tøjningerne for rosette gauges blive transformeret til det globale koordinatsystem, hvorefter
spændingerne vil blive beregnet.
x2
x1
A(13.5, −17.5)
B (70.4, −12.4)
D (79.6, −51.4)
C (26.5, −72.3)
Figur 25: Placering af rosette gauges (A, B, C og D).Koordinater i mm.
Tøjninger målt i rosettegauges for konsollen er angivet i tabel 6, hvor gauge nummer refererer til
figur 19.
Tabel 6: Tøjninger for den massive konsol.
26
B Forsøg med massiv og porøs konsol
Nr.
Tøjning ⎣⎡ ⎦⎤
μm
m
Del I Forsøg
1
2
3
4
5
6
7
8
9
10
129.16
-55.78
-541.77
-141.76
298.05
365.75
136.33
-45.52
-493.7
-135.63
Nr.
Tøjning ⎡⎣ ⎤⎦
μm
m
11
12
18
19
20
266.59
342.53
-527.02
-468.55
410.84
Tøjningerne målt i rosettegauges ønskes angivet for ε11 , ε 22 og ε12 i det globale koordinatsystem
der ses på figur 25. Til dette benyttes transformationsformlerne angivet ved (B.4).
ε i = ε11 ⋅ cos 2 βi + ε 22 ⋅ sin 2 βi + ε12 sin βi ⋅ cos βi
(B.4)
hvor
i
β
er hhv. gauge 1, 2 og 3 for rosettegaugen
er vinklen mellem den enkelte strain gauge og det globale koordinatsystem
[Hansen 1998]
Ud fra figur 26 ses geometrien af de fire rosettegauges.
x2
1
βC
βB
βA
3
x1
13.3°
2
Figur 26: Geometri af rosettegauges.
Ved indsættelse af værdier i formel (B.4) fås følgende tøjninger:
Tabel 7: Tøjninger i det globale koordinatsystem.
A
B
C
D
B.3.2
ε11 [ mm ]
ε 22 [ mm ]
ε12 [ mm ]
3.9027e-4
4.222e-4
-4.0073e-4
-4.3886e-4
-1.8337e-4
-1.9821e-4
4.3359e-5
2.6247e-5
3.8822e-5
5.2741e-5
2.6012e-4
2.8481e-4
Spændinger
For at beregne spændingerne i konsollen benyttes Hookes lov, der er givet ved følgende:
{σ } = [ D ] ⋅ {ε }
hvor
[ D]
{σ }
er spændingerne for henholdsvis σ 11 , σ 22 og σ 12
{ε }
er tøjningerne for henholdsvis ε11 , ε 22 og 2ε12
(B.5)
er materialets stivhedsmatrice
Materialet stivhed regnes ved
27
Del I Forsøg
B Forsøg med massiv og porøs konsol
[ D] =
⎡1 ν
E ⎢
ν 1
1 −ν 2 ⎢
⎢⎣ 0 0
0⎤
0 ⎥⎥
1−ν
⎥
2 ⎦
(B.6)
hvor
E
er elasticitetsmodulen beregnet i afsnit A.3
er Poissons forhold beregnet i afsnit A.3
ν
Ved indsættelse af tøjningerne i (B.5) fås følgende spændinger der ses i tabel 8.
Tabel 8: Spændinger i den massive konsol for en fladelast på 4.5 MPa.
Rosettegauge
A
B
C
D
B.3.3
ε 11 [MPa]
ε 22 [MPa]
ε 12 [MPa]
25.54
27.64
-30.06
-33.49
-3.95
-4.26
-7.24
-9.58
1.00
1.36
6.69
7.32
Effektive spændinger
Der er i dette projekt foruden σ 11 , σ 22 og σ 12 spændinger anvendt effektive Von Mises spændinger.
De effektive spændinger er defineret til følgende
σe =
3
2
sij sij
sij = σ ij − 13 σ kk δ ij
(B.7)
hvor
σe
er de effektive spændinger
er spændingsdeviatoren
sij
[Byskov 2002, p107-108]
Udtrykket i (B.7) kan omskrives til følgende udtryk ved udnyttelse af regneregler for indeksnotation
σ e = σ 112 + σ 222 + 3σ 122 − 2σ11σ 22
(B.8)
I udtrykket i (B.8) er alle indgangsværdier af σ angivet i tabel 8. De effektive Von Mises spændinger
anvendes til at kontrollere at der ikke sker flydning idet idet Von Mises flydekriterium er givet ved
σe = σ y
hvor
σy
er flydespændingen
Indsættes værdierne fra tabel 8 i formel (B.8) fås følgende værdier for den effektive spænding.
28
(B.9)
B Forsøg med massiv og porøs konsol
Del I Forsøg
Tabel 9: Effektive spændinger for massiv konsol.
Rosettegauge
29.5
32.0
25.6
27.0
A
B
C
D
B.3.4
σ e [ MPa ]
Porøs konsol
For en fladelast på 2.5 MPa bliver de samlede flytninger som angivet i tabel 10.
Tabel 10: Flytninger for porøs konsol.
Flytningsmåler
1 og 5
2 og 6
3 og 7
4 og 8
Samlede flytning [mm]
0.5033
0.2906
0.1438
0.006
Det er igen interessant at kende flytningerne på den porøse konsol mht. x1-aksen, og de er fundet på
samme måde som for den massive konsol og kan ses i tabel 11.
Tabel 11: Flytninger i x1-retning med koordinatsystem som i figur 25.
x1
0
15
60
105
150
u2
0
0.0000
0.2127
0.3595
0.4973
29
Del II
Analytiske modeller
C Simpel bjælkemodel
C
Del II Analytiske modeller
SIMPEL BJÆLKEMODEL
I dette bilag er der foretaget en bestemmelse af udbøjningen og spændingen for bjælkemodellen,
som er præsenteret i hovedrapport afsnit 5.1 og vist på figur 27. Bjælken er påvirket af en jævnt
fordelt linielast p.
p
x2 , w
x1
Figur 27: Bjælkemodel med varierende tværsnit påvirket af en jævnt
fordelt linielast p.
Udbøjningen er beregnet for både en Bernouli-Euler og en Timoshenko-bjælke efter Virtuelle Kræfters Princip på følgende fremgangsmåde:
•
•
•
•
•
•
Tværsnitsdata defineres
Virkelige snitkræfter bestemmes
Virkelige tøjninger bestemmes vha. konstitutive betingelser
Virtuel krafttilstand indføres og virtuelle snitkræfter bestemmes
Virtuelle Kræfters Princip opstilles
Udbøjningen fås ved at løse arbejdsligningen
Herefter er spændingerne beregnet for en Bernouli-Euler-bjælke ved hjælp af Naviers og Grashofs
formel. Alle beregningerne er vedlagt på cd-rom som bjælke.mws.
C.1
TVÆRSNITSDATA
Højden varierer lineært med bjælkeaksen på følgende måde
h ( x1 ) =
L x1
− , 0 ≤ x1 ≤ L
2 2
(C.1)
hvor
h
L
er bjælkens tværsnitshøjde [m]
er bjælkens længde [m]
Inertimomentet om bjælkeaksen for det rektangulære tværsnit med konstant tykkelse og varierende
højde bliver
33
Del II Analytiske modeller
C Simpel bjælkemodel
I ( x1 ) =
3
1
t ( h ( x1 ) ) ,
12
1 ⎛L x ⎞
= t⎜ − 1⎟
12 ⎝ 2 2 ⎠
0 ≤ x1 ≤ L
(C.2)
3
hvor
I
er inertimomentet om bjælkeaksen for et rektangulært tværsnit [ m 4 ]
t
er bjælkens tykkelse [m]
[Teknisk Ståbi 2004, p34]
Det effektive areal for et rektangulært tværsnit er
5
t h ( x1 )
6
5 ⎛L x ⎞
= t⎜ − 1⎟
6 ⎝2 2⎠
Ae ( x1 ) =
(C.3)
hvor
er det effektive areal for et rektangulært tværsnit [ m 2 ]
Ae
[Byskov 2002, p137]
C.2
VIRKELIGE SNITKRÆFTER
For at bestemme snitkræfterne betragtes et bjælkeudsnit i afstanden ( L − x1 ) fra bjælkespidsen, som
vist på figur 28. Normalkraften er ikke vist, fordi den er lig nul når bjælken kun er påvirket af lodrette kræfter.
p
M
V
x1
L − x1
x1
L
Figur 28: Bjælkeudsnit med forskydningskraften V, momentet M,
og den ydre linielast p.
Snitkræfterne bestemmes ud fra henholdsvis lodret ligevægt og momentligevægt
V ( x1 ) = p ( L − x1 )
1
2
M ( x1 ) = p ( L − x1 )
2
hvor
34
V
M
er den indre forskydningskraft [N]
er det indre moment [Nm]
p
er den ydre linielast ⎡⎣ mN ⎤⎦
, 0 ≤ x1 ≤ L
(C.4)
C Simpel bjælkemodel
C.3
Del II Analytiske modeller
VIRKELIGE TØJNINGER
Materialet antages at være lineært elastisk, hvorfor de konstitutive feltbetingelser kan indføres ved
Hookes lov
V ( x1 ) = G Ae ( x1 ) γ ( x1 )
, 0 ≤ x1 ≤ L
M ( x1 ) = E I ( x1 ) κ ( x1 )
(C.5)
hvor
G
er forskydningsmodulen [Pa]
γ
er forskydningstøjningen [-]
E
er elastitetsmodulen [Pa]
κ
er krumningstøjningen [-]
[Byskov 2002, p136]
Der er set bort fra aksialtøjningen, idet bjælken ikke har normalkræfter. Ved omskrivning af (C.5)
fås følgende forskrift for krumningstøjningen for både Timoshenko-bjælken og Bernoulli-Eulerbjælken
κ ( x1 ) =
=
M ( x1 )
E I ( x1 )
6 p ( L − x1 )
2
⎛L x ⎞
Et⎜ − 1 ⎟
⎝2 2⎠
48 p
=
( L − x1 ) E t
3
(C.6)
For Timoshenko-bjælken fås endvidere den konstante forskydningstøjning
γ ( x1 ) =
V ( x1 )
G Ae ( x1 )
6 p ( L − x1 )
⎛L x ⎞
5G t ⎜ − 1 ⎟
⎝2 2⎠
12 p
=
5G t
=
C.4
(C.7)
VIRTUELLE SNITKRÆFTER
Der indføres en virtuel krafttilstand, hvor bjælken i stedet for sin virkelige belastning påføres en
virtuel opadrettet enkeltkraft, som vist på figur 29. Den virtuelle kraft er placeret i en vilkårlig afstand a fra understøtningen i intervallet 0 ≤ a < L . Herved kan udbøjningen i angrebspunktet bestemmes.
35
Del II Analytiske modeller
C Simpel bjælkemodel
δp
a
L
Figur 29: Virtuel krafttilstand
med en opadrettet enkeltkraft
δP placeret i en vilkårlig afstand a fra understøtningen.
De virtuelle snitkræfter af dette system kan bestemmes ved at betragte et bjælkeudsnit i afstanden
( L − x1 ) fra bjælkespidsen, som vist på figur 30.
δP
δM
δV
a − x1
L − x1
x1
Figur 30: Bjælkeudsnit med den
virtuelle forskydningskraft δV, det
virtuelle moment δM, og den ydre
virtuelle enkeltkraft δP.
De virtuelle snitkræfter fås ved henholdsvis lodret ligevægt og momentligevægt
⎧δ P
δ V ( x1 ) = ⎨
⎩0
0 ≤ x1 < a
a < x1 ≤ L
⎧δ P ( a − x1 )
δ M ( x1 )= ⎨
0
⎩
0 ≤ x1 < a
a ≤ x1 ≤ L
(C.8)
hvor
δV
δM
δP
C.5
er den indre virtuelle forskydningskraft [N]
er det indre virtuelle moment [Nm]
er den ydre virtuelle enkeltkraft [N]
VIRTUELLE KRÆFTERS PRINCIP
De virtuelle kræfters princip med virkelige tøjninger og de virtuelle snitkræfter formuleres således
L
L
L
0
0
0
w ( a ) ⋅ δ P = ∫ κ ( x1 ) δ M ( x1 ) dx1 + ∫ γ ( x1 ) δ V ( x1 ) dx1 + ∫ ε ( x1 ) δ N ( x1 ) dx1
hvor
w
δN
ε
36
er den opadrettede udbøjning [m]
er den indre virtuelle normalkraft [Nm]
er den virkelige aksialtøjning [-]
(C.9)
C Simpel bjælkemodel
Del II Analytiske modeller
[Byskov 2002, p130]
Da aksialtøjningen er lig nul, bliver det tilhørende virtuelle indre arbejde lig nul. Da de virtuelle
snitkræfter er lig nul til højre for den virtuelle enkeltkraft jf. (C.8), bliver de øvre grænser for integralerne lig a. Arbejdsligningen (C.9) reduceres derfor til
a
a
0
0
w ( a ) ⋅ δ P = ∫ κ ( x1 ) δ M ( x1 ) dx1 + ∫ γ ( x1 ) δ V ( x1 ) dx1
C.6
(C.10)
UDBØJNING
Ved indsættelse af (C.6) og (C.8) i (C.10), integration og variabel substitution a = x1 fås udbøjningen w(x1) for Bernoulli-Euler-bjælken, idet sidste led i (C.10) udelades
w ( x1 )
BE
=
48 p ⎛
⎛ L − x1 ⎞ ⎞
⎜ x1 + ( L − x1 ) ln ⎜
⎟ ⎟ , 0 ≤ x1 < L
Et ⎝
⎝ L ⎠⎠
(C.11)
hvor
wBE
er den opadrettede udbøjning for Bernoulli-Euler-bjælken [m]
Udbøjningen i bjælkespidsen fås som grænseværdien af (C.11) når x1 → L
w( L)
BE
= lim w ( x1 )
BE
x1 → L
=
(C.12)
48 p L
Et
På tilsvarende vis indsættes (C.6)-(C.8) i (C.10) for Timoshenko-bjælken, hvis løsning kan skrives
som summen af udbøjningen hidrørende fra Bernoulli-Euler-bjælken (C.11) og et forskydningsbidrag
w ( x1 )
TS
= w ( x1 )
BE
12 p
x1
5G t
+
, 0 ≤ x1 < L
(C.13)
forskydningsbidrag
hvor
wTS
er den opadrettede udbøjning for Timoshenko-bjælken [m]
Spidsudbøjningen for Timoshenko-bjælken fås som grænseværdien af (C.13) når x1 → L
w( L)
TS
= lim w ( x1 )
TS
x1 → L
=
=
C.7
48 p L 12 p L
+
Et
5G t
(C.14)
12 p L ( 20 G + E )
5EGt
KONTROL AF UDBØJNING
Ved indsættelse af de virtuelle snitkræfter i virtuelle kræfters pincip fås et flytningsfelt, der opfylder
de kinematiske felt- og randbetingelser. I det følgende afsnit eftervises dette for løsningen af Ber37
Del II Analytiske modeller
C Simpel bjælkemodel
noulli-Euler-bjælken. De kinematiske randbetingelser kan indses af det statiske system på figur 27,
hvor indspændingen forhindrer udbøjning og vinkeldrejning
w ( 0) = 0
dw ( x1 )
∧
dx1
=0
(C.15)
x1 = 0
Ved indsættelse af løsningen (C.11) i randbetingelserne (C.15), ses at løsningen giver hverken udbøjning eller vinkeldrejning i understøtningen
dw BE ( x1 )
48 p
⎛L⎞
w BE ( 0 ) =
L ln ⎜ ⎟ = 0 ∧
Et
⎝L⎠
dx1
x1 = 0
ok
⎛L⎞
48 p ln ⎜ ⎟
⎝L⎠ =0
=
Et
(C.16)
ok
Bjælkens kinematiske feltbetingelse er
d 2 w ( x1 )
dx1
= κ ( x1 )
(C.17)
[Byskov 2002, p130]
Ved differentiation af (C.11) to gange fås
⎛
⎛ L ⎞⎞
⎜ 48 p ln ⎜
⎟⎟
d w ( x1 ) d ⎜
48 p
⎝ L − x1 ⎠ ⎟ =
=
⎜
⎟
dx1
dx1 ⎝
Et
⎠ E t ( L − x1 )
2
BE
(C.18)
Ved sammenligning med (C.6), ses at løsningen ligeledes opfylder den kinematiske feltbetingelse.
C.8
SPÆNDINGER
Spændingerne er beregnet for Bernoulli-Euler-bjælken. Normalspændingen findes ved hjælp af Naviers formel for ren bøjning:
σ 11 = −
M ( x1 )
I ( x1 )
⋅ x2
(C.19)
hvor
σ 11
er normalspændingen i x1-retningen [Pa]
[Byskov 2005, p148]
Ved indsættelse af forskrifterne for inertimomentet (C.2) og bøjningsmomentet (C.4) i (C.19) fås
normalspændingen til
σ 11 = −
48 p
⋅x
( L − x1 ) t 2
Forskydningsspændingen findes ved hjælp af Grashofs formel:
38
(C.20)
C Simpel bjælkemodel
Del II Analytiske modeller
σ 12 =
V ( x1 ) Su ( x2 )
(C.21)
I ( x1 ) t
hvor
σ 11
er forskydningsspændingen i x2-retningen [Pa]
Su
er det statiske moment omkring bjælkeaksen af et delareal af tværsnittet, hvor forskyd-
ningsspændingen ønskes bestemt [ m3 ]
[Byskov 2005, p163]
Det statiske moment af et delareal for et rektangulært tværsnit er givet ved
Su =
2
⎞
1 ⎛ ⎛ h ( x1 ) ⎞
2
t⋅⎜⎜
⎟ − x2 ⎟
⎟
2 ⎜⎝ 2 ⎠
⎝
⎠
(C.22)
[Byskov 2005, p163]
Ved indsættelse af (C.2), (C.4) og (C.22) i (C.21) fås forskydningsspændingen til
σ 12 =
3 p ( −4 x2 + L − x1 ) ( 4 x2 + L − x1 )
( L − x1 )
2
(C.23)
t
For en plan, ret bjælke kan den effektive spænding findes ved hjælp af normal- og forskydningsspændingen.
σ eff = σ 112 + 3σ 12 2
(C.24)
hvor
σ eff
er den effektive spænding [Pa]
Forskriften for de effektive spændinger er udledt i bilag G, blot er σ 22 = 0.
Ved indsættelse af spændingernes forskrifter (C.20) og (C.23) i (C.24) fås
σ eff
⎛ ( −4 x2 + L − x1 ) ⋅ ( 4 x2 + L − x1 ) ⎞
3p
=
⋅ 256 x2 2 + 3 ⎜
⎟
L − x1
( L − x1 ) t
⎝
⎠
2
(C.25)
39
D Airys spændingsfunktion
D
Del II Analytiske modeller
AIRYS SPÆNDINGSFUNKTION
Fra hovedrapport afsnit 5.2 haves Airys spændingsfunktion for konsollen på følgende form.
1 r 2 ⋅ σ ⋅ (1 − cos ( 2α ) + cos ( 2θ ) − cos ( 2θ ) ⋅ cos ( 2α ) )
4
α sin ( 2α ) + cos ( 2α ) − 1
2
1 r ⋅ σ ⋅ ( −2α sin ( 2α ) + 2θ sin ( 2α ) − sin ( 2α ) ⋅ sin ( 2θ ) )
− ⋅
4
α sin ( 2α ) + cos ( 2α ) − 1
Φ ( r ,θ ) = −
(D.1)
Spændingsfunktionen kan reduceres væsentligt ved at udnytte følgende trigonometriske sammenhænge
sin 2 a = 1 − cos 2 a
sin 2a = 2 sin a cos a
(D.2)
cos 2a = 2 cos 2 a − 1
1 − cos 2 2a = 2sin 2 a
Herved reduceres (D.1) til
Φ ( r ,θ ) = −
−
= −
−
(
)
2
2
1 r ⋅ σ ⋅ 2sin (α ) + cos ( 2θ ) ⋅ (1 − cos ( 2α ) )
4
α sin ( 2α ) + cos ( 2α ) − 1
2
1 r ⋅ σ ⋅ ( ( −2α + 2θ ) sin ( 2α ) − 2sin (α ) ⋅ cos (α ) ⋅ 2sin (θ ) ⋅ cos (θ ) )
α sin ( 2α ) + cos ( 2α ) − 1
4
2
2
2
1 r ⋅ σ ⋅ ( 2sin (α ) + cos ( 2θ ) ⋅ 2sin (α ) )
4
α sin ( 2α ) + cos ( 2α ) − 1
2
1 r ⋅ σ ⋅ ( ( −α + θ ) ⋅ 4sin (α ) ⋅ cos (α ) − 4sin (α ) ⋅ cos (α ) ⋅ sin (θ ) ⋅ cos (θ ) )
(D.3)
α sin ( 2α ) + cos ( 2α ) − 1
4
(
)
2
2
1 r ⋅ σ ⋅ 2sin (α ) ⋅ (1 + cos ( 2θ ) )
= −
α sin ( 2α ) + cos ( 2α ) − 1
4
2
1 r ⋅ σ ⋅ ( ( −α + θ ) ⋅ 4sin (α ) ⋅ cos (α ) − 4sin (α ) ⋅ cos (α ) ⋅ sin (θ ) ⋅ cos (θ ) )
−
α sin ( 2α ) + cos ( 2α ) − 1
4
41
Del II Analytiske modeller
= −
−
D Airys spændingsfunktion
2
2
2
1 r ⋅ σ ⋅ ( 2sin (α ) ⋅ 2cos (θ ) )
4 α ⋅ 2sin (α ) cos (α ) − 2sin 2 (α )
2
1 r ⋅ σ ⋅ ( ( −α + θ ) ⋅ 4sin (α ) ⋅ cos (α ) − 4sin (α ) ⋅ cos (α ) ⋅ sin (θ ) ⋅ cos (θ ) )
α ⋅ 2sin (α ) cos (α ) − 2sin 2 (α )
4
2
2
1 r ⋅ σ ⋅ ( 4sin (α ) ⋅ cos (α ) ⋅ tan (α ) ⋅ cos (θ ) )
=−
8
α ⋅ sin (α ) cos (α ) − sin 2 (α )
2
1 r ⋅ σ ⋅ ( ( −α + θ ) ⋅ 4sin (α ) ⋅ cos (α ) − 4sin (α ) ⋅ cos (α ) ⋅ sin (θ ) ⋅ cos (θ ) )
−
8
α ⋅ sin (α ) cos (α ) − sin 2 (α )
= −
=
2
2
1 r ⋅ σ ⋅ 4sin (α ) ⋅ cos (α ) ⋅ ( tan (α ) ⋅ cos (θ ) − α + θ − sin (θ ) cos (θ ) )
8
α ⋅ sin (α ) cos (α ) − sin 2 (α )
2
2
1 r ⋅ σ ⋅ cos (α ) ⋅ (α − θ + sin (θ ) cos (θ ) − tan (α ) ⋅ cos (θ ) )
α ⋅ cos (α ) − sin (α )
2
Der indføres en konstant
C=
1
cos (α )
2 α ⋅ cos (α ) − sin (α )
(D.4)
Herved kan spændingsfunktionen endeligt omskrives til en simplere form
Φ ( r ,θ ) = C ⋅ r 2 ⋅ σ ⋅ (α − θ + sin (θ ) ⋅ cos (θ ) − cos 2 (θ ) ⋅ tan (α ) )
42
(D.5)
Del III
Numeriske modeller
E Bestemmelse af materialeparametre ved beregninger
E
Del III Numeriske modeller
BESTEMMELSE AF
MATERIALEPARAMETRE VED
BEREGNINGER
E.1
ELASTICITETSMODUL OG POISSONS FORHOLD
E.1.1
Udskrivning af spændingsrelation
I det følgende betragtes relationen mellem lokale og globale spændinger givet ved (E.1).
σ ij =
1
V
∫σ
S
lokal
ik
nk x j dS
(E.1)
[Jensen 2006, note 6]
For at illustrere anvendelsen af (E.1) tages der udgangspunkt i figur 31, som viser de enkelte spændingskomposanter σ ij på RVE samt enhedsnormalvektorerne n på de to ydre flader givet ved
⎡1⎤
nvandret = ⎢ ⎥
⎣0⎦
⎡0⎤
nlodret = ⎢ ⎥
⎣1⎦
for den vandrette rand
for den lodrette rand
45
Del III Numeriske modeller
E Bestemmelse af materialeparametre ved beregninger
σ 22
nlodret
σ 21
σ 12
a
σ 11
x2
nvandret
x1
a
Figur 31: Spændingskomposanter, tøjningstilstand og normalvektorer på udsnit af
RVE. Stiplet omrids angiver udeformeret tilstand.
I det følgende gennemgås i detaljer, hvordan den globale spændingskomposant σ11 bestemmes.
Da problemet regnes plant, skal indeksnotationen i (E.1) kun evalueres for to dimensioner. Udtrykket skal gennemløbes for både den vandrette og den lodrette rand. Udskrevet generelt for σ11 på en
rand bliver (E.1)
σ 11 =
⎤
1⎡
lok
lok
⎢ ∫ σ 11 n1 x1 dS + ∫ σ 12 n2 x1 dS ⎥
V ⎣S
S
⎦
(E.2)
Her er n1 og n2 hhv. x1 - og x2 -komposanten i den enhedsnormalvektor, som tilhører den betragtede rand.
For den lodrette rand er komposanterne af enhedsnormalen hhv. n1 = 1 og n2 = 0 . Derfor forsvinder
andet led for denne rand. Da x1 -koordinaten er konstant a på den lodrette rand varierer denne således ikke i summationen. Da problemet modelleres vha. elementmetoden er det en diskret summation
af knudekræfter der udføres i modsætning til en egentlig integration. σ11 -bidraget fra den lodrette
rand bliver således
a
σ 11 =
1
σ 11lok a dx2
V ∫0
Tilsvarende er der kun ét bidrag fra den vandrette rand, idet dennes første komposant i normalenhedsvektoren er 0, hvorfor kun andet led i (E.2) er forskellig fra nul. I dette tilfælde er stedkoordinaten x1 ikke konstant, men varierer fra 0 til a, hvorfor dette bidrag ser således ud
46
E Bestemmelse af materialeparametre ved beregninger
Del III Numeriske modeller
a
σ 11 =
1
σ 12lok x1 dx1
∫
V 0
For σ 22 gennemføres samme operation, hvorfor de samlede udtryk for σ 11 og σ 22 kommer til at se
ud som følger
a
a
a
⎤ 1 ⎡ a lok
⎤
1⎡
lok
σ
x
dx
+
σ
x
dx
=
σ
a
dx
+
σ 12 x1 dx1 ⎥
⎢ ∫ 11 1 2 ∫ 12 1 1 ⎥
⎢ ∫ 11
2
∫
V ⎣0
0
0
⎦ V ⎣0
⎦
a
a
a
a
⎤ 1⎡
⎤
1⎡
= ⎢ ∫ σ 21lok x2 dx2 + ∫ σ 22 x2 dx1 ⎥ = ⎢ ∫ σ 21lok x2 dx2 + ∫ σ 22 a dx1 ⎥
V ⎣0
0
0
⎦ V ⎣0
⎦
σ 11 =
σ 22
E.1.2
Konvergensundersøgelse af spændingsrelation
Figur 32 og figur 33 viser resultatet af en konvergensundersøgelse af værdierne for σ 11 og σ 22 som
MPa
funktion af antallet af CST-elementer. Som det ses, opnås tilnærmelsesvis stationære værdier ved
anvendelse af 10.000 elementer.
29700
29600
29500
29400
29300
29200
29100
29000
28900
0
2000
4000
6000
8000
10000
Antal elementer
Figur 32: Konvergensundersøgelse af σ 11 for stigende antal CST-elementer.
6500
MPa
6450
6400
6350
6300
0
2000
4000
6000
8000
10000
Antal elementer
Figur 33: Konvergensundersøgelse af σ 22 for stigende antal CST-elementer.
47
Del III Numeriske modeller
E Bestemmelse af materialeparametre ved beregninger
Det vælges således at regne med følgende værdier for n = 10.000 bliver således
σ11 = 31 895 MPa
σ 22 = 7 471MPa
E.2
FORSKYDNINGSMODUL
E.2.1
Transformation af forskydningsspændinger
I dette afsnit omregnes spændingerne i det drejede koordinatsystem til forskydningsspændinger i det
oprindelige koordinatsystem. Figur 34 viser spændingerne udregnet ved elementmetoden
'
, σ 12' ) og de spændinger der ønskes fundet vha. transformationsformlerne ( σ 11 , σ 22 , σ 12 ).
( σ 11' , σ 22
σ 22'
σ 12'
'
2
x
σ 11'
σ 12'
σ 12'
σ 11'
σ 12'
11
σ
12
σ
σ
σ
{n}
σ 22
σ
σ 22'
12
{m}
θ
12
12
σ
σ 22
11
x1'
Figur 34: Udsnit af RVE til bestemmelse af forskydningsspændinger.
Retningsvektorerne {n} og {m} kan skrives som
⎡ n ⎤ ⎡ − sin(θ ) ⎤
{n} = ⎢ 1 ⎥ = ⎢
⎥
⎣ n2 ⎦ ⎣cos(θ ) ⎦
(E.3)
⎡ m ⎤ ⎡cos(θ ) ⎤
{m} = ⎢ 1 ⎥ = ⎢
⎥
⎣ m2 ⎦ ⎣sin(θ ) ⎦
'
Transformationsformlerne benyttes til at beskrive sammenhængen mellem spændingerne σ αβ
i det
drejede koordinatsystem, og spændingerne σ αβ i det oprindelige koordinatsystem. Transformationsformlerne er givet ved (E.4) og udskrevet for de enkelte spændingskomposanter.
48
E Bestemmelse af materialeparametre ved beregninger
Del III Numeriske modeller
σ 11' = σ αβ nα nβ
σ 11' = σ 11n12 + σ 22 n22 + 2σ 12 n1n2
σ 22' = σ αβ mα mβ
σ 22' = σ 11m12 + σ 22 m22 + 2σ 12 m1m2
(E.4)
σ 12' = σ 21' = σ αβ nα mβ
σ 12' = σ 11n1m1 + σ 22 n2 m2 + σ 12 n1m2 + σ 12 n2 m1
[Jensen 2006, Note 4]
3 ligninger med 3 ubekendte løses, og giver følgende løsninger til spændingerne udtrykt ved de oprindelige koordinater. Beregninger vedlagt på cd-rom som Transformationsformler.mw
σ 11 = −2 cos(θ )sin(θ )σ 12' + cos2 (θ )σ 22' + σ 11' − cos2 (θ )σ 11'
σ 22 = 2 cos(θ )sin(θ )σ 12' + cos 2 (θ )σ 11' + σ 22' − cos 2 (θ )σ 22'
σ 12 = cos(θ )sin(θ )(σ 22' − σ 11' ) + 2 cos 2 (θ )σ 12' − σ 12'
Det er kun forskydningsspændingen σ 12 der anvendes til bestemmelse af forskydningsmodulen.
Derfor er det kun den der behandles i det følgende. σ 12 omskrives vha. følgende trigonometriske
funktioner.
sin(2θ ) = 2sin(θ ) cos(θ )
1
sin(2θ ) = cos(θ )sin(θ )
2
(E.5)
2 cos 2 (θ ) = (1 + cos 2(θ ))
(E.6)
Indføres de trigonometriske funktioner (E.5) og (E.6) i udtrykket for forskydningsspændingen fås
1
2
1
'
= sin(2θ )(σ 22
− σ 11' ) + σ 12' + cos(2θ )σ 12' − σ 12'
2
1
'
= sin(2θ )(σ 22
− σ 11' ) + cos(2θ )σ 12'
2
σ 12 = sin(2θ )(σ 22' − σ 11' ) + (1 + cos(2θ ))σ 12' − σ 12'
(E.7)
49
F Elementmetodeteori
F
Del III Numeriske modeller
ELEMENTMETODETEORI
I dette bilag foretages en gennemgang af teorien bag elementmetoden. Der tages udgangspunkt i
forskriften for potentiel energi, og det vises hvordan der ved brug af variationsprincippet kan opstilles en matrixligning, som kan anvendes i elementmetoden. Gennemgangen bliver foretaget ud fra en
opstilling af den potentielle energi for konsollen i dette projekt.
F.1
POTENTIEL ENERGI
Den potentielle energi beskriver den samlede energi i et system i form af tøjningsenergi, samt potentiel energi fra volumenkræfterne og kræfterne langs randen af legemet. Den potentielle energi for
kinematisk linearitet og lineær elasticitet er opstillet generelt i følgende udtryk
Π P (ui ) =
∫
1
2 V
Dijkl ε ij ε kl dV − ∫ qi ui dV − ∫ τ i ui dS
V
S
(F.1)
hvor
ΠP
er den potentielle energi
u
D
er flytningsfeltet for legemet
er en 4. ordens tensor, der beskriver den konstitutive relation mellem tøjning og spænding
er tøjningstensoren
er volumenet af det udeformerede legeme
er volumenkræfter
ε
V
q
τ
er kræfter på randen af volumenet
S
er arealet af randen
[Byskov 2002, p90]
Det første led på højresiden af (F.1) beskriver tøjningsenergien, det midterste led den potentielle
energi fra volumenkræfterne, mens det sidste led beskriver den potentielle energi fra kræfterne langs
randen af legemet. Π P er en funktion af ui, da den kan opstilles ud fra at vilkårligt flytningsfelt beskrevet ved ui.
I dette projekt ses der bort fra volumenkræfter som for eksempel egenlasten, da denne har lille variation gennem legemet. Dermed kan de i stedet påføres den ydre rand som en fladelast. Dermed bliver
den potentielle energi i dette tilfælde
51
Del III Numeriske modeller
F Elementmetodeteori
Π P (ui ) =
F.1.1
∫
1
2 V
Dijkl ε ij ε kl dV − ∫ τ i ui dS
S
(F.2)
Notation af Dijkl
For at den potentielle energi kan opstilles på matrixform er det nødvendigt, at notationen af den 4.
ordens tensor Dijkl opstilles på matrixform. Dette gøres ved at samle den 2. ordens tensor ε i en vektor. Dette medfører, at relationen mellem Dijkl og εij kan noteres ved følgende
Dijkl ε ij ε kl → Dij ε i ε j
[Byskov 2002, p93]
Dij kaldes materialets stivhedsmatrice.
F.1.2
Opstilling af flytningsfelt
Den potentielle energi kan opstilles for et vilkårligt valg af flytningsfelt, så længe det overholder de
kinematiske betingelser. For dette projekt omhandlende et plant tilfælde er de kinematiske betingelser, at der er kontinuitet i flytningerne i x1 og x2-aksens retninger.
Flytningsfeltet i (F.1) kan beskrives ved at indføre elementer med knudepunkter med tilhørende
frihedsgrader som styrer flytningsfeltet. Hermed kan flytningsfeltet beskrives ved en relation mellem
flytningerne i frihedsgraderne og det totale flytningsfelt som opstillet i det følgende
{u} = [ N ]{V }
(F.3)
hvor
er flytningsfeltet i x1 og x2-aksens retninger
er flytningerne i frihedsgraderne
er en flytningsinterpolationsmatrice som knytter relationen mellem frihedsgrader og
flytningsfelt sammen
[Byskov 2002, p371]
u
V
N
Det skal bemærkes, at de indførte frihedsgrader kan betragtes som koefficienter til et funktionsudtryk, hvis værdier indtil videre kan vælges frit, da (F.1) gælder for et vilkårligt valg af flytningsfelt,
som overholder de kinematiske betingelser.
F.1.3
Opstilling af tøjningsfelt
Tøjningstensoren i (F.1) kan ligeledes beskrives ved de førnævnte frihedsgrader, da tøjningerne i en
skive afhænger lineært af flytningerne. Sammenhængen mellem tøjning og flytning udtrykkes ved
følgende sammenhæng
{ε } = [ B]{V }
(F.4)
hvor
B
er tøjningsinterpolationsmatricen
[Byskov 2002, 372]
F.1.4
Opstilling af last
Kraften på randen af legemet i (F.1) omskrives til ækvivalente punktlaster virkende i frihedsgraderne, hvilket giver følgende relation
52
F Elementmetodeteori
Del III Numeriske modeller
∫ τ u dS = {V } { f }
T
S
(F.5)
i i
hvor
f
er den ækvivalente kraftvektor virkende i frihedsgraderne
[Byskov 2002, p372]
Bemærk at frihedsgradsvektoren er transponeret således at matrixdimensionerne i matrixproduktet
stadig stemmer overens.
F.1.5
Potentiel energi på matrixform
Ved indsættelse af (F.4) og (F.5) i (F.2) kan den potentielle energi nu opstilles på matrixform ved
følgende
Π P ({V }) =
1
2
=
1
2
=
1
2
=
1
2
∫ ([ D][ B]{V }) [ B]{V }dA −{V } { f }
∫ ([ B]{V }) [ D] [ B]{V }dA −{V } { f }
∫ {V } [ B] [ D][ B]{V }dA −{V } { f }
{V } ∫ [ B] [ D][ B]dA{V } −{V } { f }
T
T
A
T
T
T
A
T
T
T
T
T
A
T
(F.6)
A
= 12 {V }T [ K ]{V } − {V }T { f }
hvor
[K]
er konstruktionens stivhedsmatrice defineret ved [ K ] = ∫ [ B ]T [ D ][ B ]dA , idet tykkelsen
A
er indeholdt i [D]
[Byskov 2002, p373]
Ved omskrivningerne i (F.6) er anvendt, at den transponerede af et matrixprodukt er lig de enkelte
matricer transponeret og i omvendt rækkefølge. Ligeledes er det anvendt at materialets stivhedsmatrix er symmetrisk, da Dijkl = Dklij for lineær elastisk materiale [Byskov 2002, p60]. Dette medfører at
[ D]T = [D]. Volumenintegralet er blevet ændret til et fladeintegral, da der i dette projekt behandles
et plant tilfælde. Den potentielle energi er nu en funktion af frihedsgraderne.
F.2
VARIATION AF POTENTIEL ENERGI
Til at bestemme de korrekte koefficienter i flytningsfeltet beskrevet ved frihedsgraderne kan princippet om stationær potentiel energi benyttes som angivet i (F.7). Princippet går ud på at variere den
potentielle energi ved variation af koefficienterne i det estimerede flytningsfelt således, at den mindste potentielle energi opnås. Princippet udnytter, at den potentielle energi er stationær, når variationen er nul. Dette giver det bedste flytningsfelt indenfor den valgte funktionsform. Princippet om
stationær potentiel energi kan defineres ved følgende
δΠ P (α ) =
∂Π P (α + εδα )
=0
∂ε
ε =0
(F.7)
hvor
δΠ P
er variationen af den potentielle energi
α
ε
δα
er et felt, eksempelvis et flytningsfelt
er amplituden af en lille afvigelse fra den potentielle energi
er formen af en lille variation fra den potentielle energi
53
Del III Numeriske modeller
F Elementmetodeteori
[Byskov 2002, p517]
Udtrykket i (F.7) anvendes i det følgende, hvor feltet α substitueres med flytningsfeltet beskrevet
ved {V}. Det skal bemærkes, at der er igennem hele gennemgangen er anvendt stivhedsmatricer og
frihedsgrader på konstruktionsniveau. De anvendte principper gælder ikke på elementniveau, idet
variation af frihedsgrader skal ske globalt for at sikre konstruktionens sammenhæng.
Anvendes (F.6) til at opstille den potentielle energi for flytningen V + εδ V i stedet for V, fås følgende
Π P ({V } + ε
{δ V }) =
1
2
({V } + ε{δ V })T [ K ] ({V } + ε{δ V }) − ({V } + ε{δ V })T { f }
= 12 {V }T [ K ]{V } + 12 {V }T [ K ]ε{δ V } + 12 ε{δ V }T [ K ]{V }
+ 12 ε 2 {δ V }T [ K ]{δ V } − {V }T { f } − ε{δ V }T { f }
= 12 {V }T [ K ]{V } + 12 ε ({V }T [ K ]{δ V }) + 12 ε{δ V }T [ K ]{V }
T
+ 12 ε 2 {δ V }T [ K ]{δ V } − {V }T { f } − ε{δ V }T { f }
= 12 {V }T [ K ]{V } + 12 ε ([ K ]{δ V }) {V } + 12 ε{δ V }T [ K ]{V }
T
+ 12 ε 2 {δ V }T [ K ]{δ V } − {V }T { f } − ε{δ V }T { f }
= 12 {V }T [ K ]{V } + 12 ε{δ V }T [ K ]{V } + 12 ε{δ V }T [ K ]{V }
+ 12 ε 2 {δ V }T [ K ]{δ V } − {V }T { f } − ε{δ V }T { f }
= 12 {V }T [ K ]{V } + ε{δ V }T [ K ]{V }
+ 12 ε 2 {δ V }T [ K ]{δ V } − {V }T { f } − ε{δ V }T { f }
(F.8)
I (F.8) er det udnyttet, at ε er en skalar, som således kan flyttes frit rundt i leddet. Der er ligeledes
anvendt, at hvert led i udtrykket for potentiel energi er en skalar, da den potentielle energi er en skalar. Dermed er et leds transponerede lig ledet inden transponering, hvilket er anvendt til omskrivning
af leddet 12 {V }T [ K ]ε{δ V } = 12 ε ({V }T [ K ]{δ V }) .
T
Den afledte af den potentielle energi i (F.8) med hensyn til ε bliver således
∂Π P ({V } + ε{δ V })
= {δ V }T [ K ]{V } + ε{δ V }T [ K ]{δ V } − {δ V }T { f }
∂ε
(F.9)
Ved nu at anvende definitionen af princippet om stationær potentiel energi i (F.7), hvor ε i udtrykket i (F.9) sættes lig nul fås følgende
δΠ P (V ) = {δ V }T [ K ]{V } − {δ V }T { f } = 0
{δ V }T ([ K ]{V } − { f }) = 0
(F.10)
Princippet om stationær potentiel energi skal gælde for alle variationer af {V}, hvilket medfører
følgende matrixligning
54
F Elementmetodeteori
Del III Numeriske modeller
{δ V }T ([ K ]{V } − { f }) = 0
(F.11)
∀{δ V }
[ K ]{V } = { f }
Hermed er matrixligningen, som anvendes i den anvendte elementmetode opstillet ud fra den potentielle energi.
F.3
KONVERGENS AF POTENTIEL ENERGI
Det er muligt at vise, at den potentielle energi for et estimeret flytningsfelt altid er større end den
potentielle energi for det eksakte flytningsfelt. Dermed kan der foretages en konvergensanalyse af de
estimerede flytningsfelter svarende til forskelligt antal elementinddelinger af legemet, og den potentielle energi vil aftage monotont for forøgelse af antal elementinddelinger.
Det vises i det følgende, hvorledes dette udsagn kan bevises. Det er kun gældende for lineært elastisk
materiale. Den eksakte potentielle energi er udledt i (F.6). Den potentielle energi for et estimeret
flytningsfelt kan udtrykkes som følgende
({ }) = {V} [ K ]{V} − {V} { f }
Π P V
T
1
2
T
(F.12)
hvor
er flytningerne i frihedsgraderne svarende til det estimerede flytningsfelt
V
[Byskov 2002, p519]
Det estimerede flytningsfelt kan også udtrykkes som det eksakte flytningsfelt adderet med en lille
afvigelse som angivet i (F.13)
V = V + εδ V
[Byskov 2002, p519]
(F.13)
Ved indsættelse af (F.13) i (F.12) kan den potentielle energi for et estimeret flytningsfelt omskrives
til
({ } ) =
1
2
=
1
2
Π P V
{V + εδ V }
T
[ K ]{V + εδ V } − {V + εδ V }
T
{f}
{f}
+ε ({V } [ K ]{δ V } − {δ V } { f })
{V } [ K ]{V } − {V }
T
T
T
(F.14)
T
+ 12 ε 2 {δ V } [ K ]{δ V }
T
Det første led i resultatet fra (F.14) er ved sammenligning med (F.6) den potentielle energi for det
eksakte flytningsfelt. Parentesen i andet led svarer til virtuelle flytningers princip. Dette kan ses ved
at regne tilbage i (F.6) og samtidig anvende definitionerne i (F.5) og (F.4). Dermed kan andet led
skrives som
{V } [ K ]{δ V } − {δ V }
T
T
{f}= ∫
A
Dij ε iδε j dA − ∫ τ iδ ui dS = 0
S
(F.15)
[Byskov 2002, p515]
Hvis den potentielle energi for den eksakte løsning er minimum, skal udtrykket givet ved den øverste
linie i (F.16) være sandt. Ved indsættelse af værdier for den potentielle energi for de to flytningsfelter og ved anvendelse af (F.15) fås udsagnet i nederste linie af (F.16). Dette angiver, at tøjningsener55
Del III Numeriske modeller
F Elementmetodeteori
gien fra et vilkårligt flytningsfelt altid vil være positiv, hvilket er et udsagn, der altid er sandt, da det
ellers vil betyde, at der kan genereres energi af et system ved elastisk deformation.
({ }) > Π ({V }) ,
Π P V
1
2
{V } [ K ]{V } − {V }
T
T
P
V ≠ V
{ f } + ε ({V } [ K ]{δ V } − {δ V } { f })
T
T
+ 12 ε 2 {δ V } [ K ]{δ V } >
T
1
2
1
2
{V } [ K ]{V } − {V }
T
T
{ }
(F.16)
f
ε 2 {δ V } [ K ]{δ V } > 0
T
[Byskov 2002, p519]
F.4
KONVERGENS AF LASTENS ARBEJEDE
I det følgende vises, at lastens arbejde for det eksakte flytningsfelt altid er større end ved det estimerede flytningsfelt fundet elementberegninger ud fra potentiel energi. Der tages udgangspunkt i (F.16)
, der ved indsættelse af (F.6) giver
1
2
{V }T [ K ]{V } − {V }T { f } < 12 {V }T [ K ]{V } − {V }T { f }
(F.17)
(F.17) omskrives ved at anvende princippet om stationær potentiel energi som formuleret i (F.10).
Her udnyttes det, at princippet skal gælde for vilkårlige værdier af δ V og δ V , så længe flytningsfeltet overholder de kinematiske betingelser. Ved at sætte δ V = V og δ V = V fås ved indsættelse i
(F.10) følgende relation
{V }T [ K ]{V } − {V }T { f } = 0
{V }T [ K ]{V } = {V }T { f }
(F.18)
{V }T [ K ]{V } − {V }T { f } = 0
{V }T [ K ]{V } = {V }T { f }
Ved indsættelse af relationerne fra (F.18) i (F.17) fås følgende
1
2
{V }T { f } − {V }T { f } < 12 {V }T { f } − {V }T { f }
− 1 {V }T { f } < − 1 {V }T { f }
2
2
(F.19)
{V }T { f } > {V }T { f }
Uligheden i (F.19) viser, at kræfternes arbejde for det eksakte flytningsfelt er større end for det estimerede flytningsfelt. Dette er ensbetydende med, at konstruktionen modelleres for stiv ved anvendelse af elementmetoden for potentiel energi. Det mest eksakte flytningsfelt ved anvendelse af elementmetoden fås dermed ved anvendelse af den elementfordeling, der giver den største værdi af
kræfternes arbejde.
56
G CST- og LST-elementopbygning
G
Del III Numeriske modeller
CST- OG LSTELEMENTOPBYGNING
I dette bilag opstilles en CST-model (Constant Strain Triangle) og en LST-model (Linear Strain
Triangle) for projektkonsollen. Formålet ved opstillingen af CST- og LST-modellen er at analysere
antallet af frihedsgrader, der skal udføres for at opnå tilpas nøjagtighed på nedbøjningen samt beregnings- og opstillingstiden for de to elementtyper. Inddelingen af konsollen i N antal lag svarende til 4
elementer samt det globale koordinatsystem ses på figur 35.
b
N
b
a
a
N
x2
x1
Figur 35: Systemopdeling af konsol.
G.1
OPBYGNING AF CST-ELEMENTER
Dette afsnit omhandler opbygningen af CST-modellen for projektkonsollen. Der vil, hvor ikke andet
er angivet, blive anvendt [Byskov, 2002]. Der opstilles en lokal stivhedsmatrice for et enkelt delelement, hvorefter det beskrives, hvordan de lokale stivhedsmatricer placeres i den globale stivhedsmatrice. Der vil sideløbende være henvisninger til Maple og Matlab, der er anvendt til beregningerne.
Til beregning af elementets stivhedsmatrice benyttes:
[k ] = ∫ [ B]T [ D][ B]dA
(G.1)
A
hvor
[B]
er tøjningsfordelingsmatricen
[D]
er materialets stivhedsmatrice
[k]
er elementets stivhedsmatrice
[Byskov 2002, p400]
57
Del III Numeriske modeller
G CST- og LST-elementopbygning
For at opstille tøjningsfordelingsmatricen betragtes et delelement af konsollen, som ses på figur 36,
hvor antallet af frihedsgrader er vist. For en fuldstændig modellering af konsollen er det nødvendigt
med to forskellige elementer. Et element med vandret underkant og et element med vandret overkant. Nedenfor er gennemgået opstillingen af stivhedsmatricen for elementet med vandret overkant.
x2
v4
v6
v3
v5
1
v2
x1
v1
Figur 36: Lokale frihedsgrader på CST-element.
G.1.1
Flytningsinterpolationsmatrice [N]
CST-elementer har lineære flytningsfelter igennem elementet. Det er en forudsætning, at der ikke
sker overlap mellem delelementet og tilstødende delelementer, hvorfor det entydige flytningsfelt skal
beskrives ved lineære sammenhænge mellem flytningerne i frihedsgraderne. Der vil i det følgende
blive opstillet enhedsflytninger af de enkelte frihedsgrader hver for sig. På figur 37 er enhedsflytningerne i x1 -retningen vist.
v5 = 1
v3 = 1
2
3
2
3
1
1
2
3
1
v1 = 1
Figur 37: Mulige flytninger i x1-retning.
Til at beskrive sammenhængen mellem frihedsgrader og flytninger i elementet anvendes relationen
{u} = [ N ]{v}
(G.2)
hvor
{u}
[N ]
er flytningsvektoren i x1 og x2 retningen
er flytningsinterpolationsmatricen
{v}
er frihedsgradsvektoren
[Byskov, 2002, p378]
Formen af [ N ] ses af (G.3), hvor indgangene kan bestemmes ved at se på flytningstilstandene i figur
37.
58
G CST- og LST-elementopbygning
Del III Numeriske modeller
⎡N
N12
N 22
[ N ] = ⎢ N11
⎣
21
N13
N 23
N14
N 24
N15
N 25
N16 ⎤
N 26 ⎥⎦
(G.3)
Fremgangsmåden til bestemmelse af [ N ] er at give en frihedsgrad en enhedsflytning i en af hovedretningerne, hvor de resterende frihedsgrader fastholdes, hvilket er angivet i figur 37. Det ses, at der
ved de opstillede tilfælde kun er en indgang i N der er forskellig fra nul, da kun en frihedsgrad er
forskellig fra nul og flytningen for hvert tilfælde er nul i x2 retningen. Det er derfor muligt at bestemme indgangene enkeltvis ved at opstille en forskrift for flytningstilstanden for hvert tilfælde som
så bliver lig den respektive indgang i [N]. Flytningsfeltet for CST-elementet ved enhedsflytning af
frihedsgrad v1 kan beskrives på følgende form.
N11 ( x1 , x2 ) = C ⋅ l23
(G.4)
hvor
C
er en konstant
l23
er ligningen for linien mellem knude 2 og 3 sat lig nul
[Byskov, 2002, p394]
Et flytningsfelt af denne form er gældende for alle enhedsflytninger, idet alle kinematiske betingelser
på denne måde sikres. En flytning lig nul langs den udeformerede rand sikres ved l23 og en enhedsflytning i v1 sikres ved bestemmelse af konstanten C. Lineære flytninger sikres idet x indgår ved 1.
orden.
For en enhedsflytning af u1 i x1 -retningen, som ses på figur 37, beregnes indgang N11 . For den
udeformerede rand l23 , findes forskriften for delelement 1 på figur 36
l23 = x2 −
a
=0
N
hvor
a
N
er højde af hele konsollen
er antallet af rækker ved inddelingerne af delelementerne
Til bestemmelse af konstanten C benyttes, at enhedsflytningen i x1 -retningen er lig 1 i punktet
( 0, 0 ) , hvorfor følgende relation benyttes:
N11 ( 0,0 ) = C ⋅ l23 = 1
Ved indsættelse af l23 beregnes:
C=−
N
a
Herved er indgang N11 følgende:
N11 ( x1 , x2 ) = −
a
⋅ x2 + 1
N
Ved at betragte figur 37 findes, at følgende gælder:
59
Del III Numeriske modeller
G CST- og LST-elementopbygning
N11 = N 22 , N13 = N 24 , N15 = N 26
Dette er gældende, da den udeformerede rand er den samme for hvert tilfælde. Analogt beregnes de
resterende indgange i [ N ] , hvorved den er følgende:
N ⋅ ( ax1 − bx2 )
⎡ N
−
0
⎢ − ⋅ x2 + 1
a ⋅b
[N ] = ⎢ a
⎢
N
− ⋅ x2 + 1
0
0
⎢
a
⎣
0
−
hvor
a
b
N
N ⋅ ( ax1 − bx2 )
N ⋅ x1
b
a ⋅b
0
⎤
0 ⎥
⎥
N ⋅ x1 ⎥
⎥
b ⎦
er højden af konsollen, se figur 35
er længden af konsollen, se figur 35
er antal inddelinger af konsollen i højden eller bredden, se figur 35
G.1.2
Tøjningsfordelingsmatrice [B]
Til beregningen af tøjningsfordelingsmatricen [ B ] anvendes flytnings-tøjningsrelationen i (G.5),
hvor der er set bort fra ikke-lineære led, da det antages at flytningerne er små
ε αβ =
1
( uα , β + uβ ,α )
2
(G.5)
hvor
er tøjningstensoren
εαβ
[Byskov 2002, p392]
Tøjningsfordelingsmatricen beregnes ud fra følgende:
⎡ ε1 ⎤ ⎡ ε11 ⎤ ⎡ u1,1 ⎤
{ε ( xα )} ≡ ⎢⎢ε 2 ⎥⎥ ≡ ⎢⎢ ε 22 ⎥⎥ = ⎢⎢ u2,2 ⎥⎥
⎢⎣ε 3 ⎥⎦ ⎢⎣ 2ε12 ⎥⎦ ⎢⎣u1,2 + u2,1 ⎥⎦
(G.6)
[Byskov 2002, p395]
Ud fra ovenstående relationer kan flytningsinterpolationsmatricen indsættes, hvor (G.2) udnyttes
⎡ u1 ⎤ ⎡ N1 j ⎤
⎥ ⋅ vj
⎥=⎢
⎣u2 ⎦ ⎣ N 2 j ⎦
{u} = ⎢
(G.7)
Således bliver tøjningsfordelingsmatricen [ B ] følgende
∂N1 j
⎡
⎤
⎢
⎥
∂x1
⎢
⎥
⎡ B1 j ⎤ ⎢
⎥
∂
N
2j
⎡⎣ B ( xα ) ⎤⎦ = ⎢⎢ B2 j ⎥⎥ = ⎢
⎥
∂x2
⎢
⎥
⎢⎣ B3 j ⎥⎦ ⎢
∂N
∂N ⎥
⎢ 1j + 2 j ⎥
⎢⎣ ∂x2
∂x1 ⎥⎦
(G.8)
[Byskov 2002, p395]
Ved differentiation af [ N ] , vedlagt på cd-rom som "bestemmelse af k1 og k2 homogen bjaelke.mw”
fås:
60
G CST- og LST-elementopbygning
Del III Numeriske modeller
⎡ N
⎢− b
⎢
⎡⎣ B ( xα ) ⎤⎦ = ⎢ 0
⎢
⎢
⎢ 0
⎣⎢
−
0
N
b
0
0
N
b
−
0
N
a
N
b
−
N
a
0
0
N
a
⎤
0⎥
⎥
N⎥
a⎥
⎥
0⎥
⎦⎥
(G.9)
Det bemærkes, at tøjningsfordelingsmatricen ikke afhænger af x1 og x2 , hvorfor denne metode kaldes Constant Strain Triangle.
Materialets stivhedsmatrice [D]
G.1.3
For isotrope lineære elastiske materiale, kan materialets stivhedsmatrice [D] benyttes og er udledt i
afsnit G.2 og givet ved følgende
0 ⎤
⎡1 ν
⎢ν 1
0 ⎥⎥
E ⋅t ⎢
⋅
[ D] =
1 −ν ⎥
1 −ν 2 ⎢
⎢0 0
⎥
2 ⎦
⎣
(G.10)
[Byskov, 2002, p192]
Elementets stivhedsmatrice [k]
G.1.4
Ved indsættelse af de fundne matricer, [B] og [D], i formel (G.1), fås den lokale stivhedsmatrice for
CST-elementet med vandret overside
bEt
⎡
⎢ 4 (ν + 1) a
⎢
⎢
0
⎢
⎢
⎢
bEt
⎢
⎢− 4 ν + 1 a
(
)
⎢
k=⎢
Et
⎢
⎢ 4 (ν + 1)
⎢
⎢
0
⎢
⎢
⎢
⎢ − Et
⎢ 4 (ν + 1)
⎣
−
bEt
2a ( −1 + ν 2 )
−
Etν
2 ( −1 + ν 2 )
bEt
4 (ν + 1) a
Et
4 (ν + 1)
0
Etν
2 ( −1 + ν 2 )
bEt
2a ( −1 + ν 2 )
Etν
2 ( −1 + ν 2 )
Et
4 ( −1 + ν )
Eta
2b ( −1 + ν 2 )
−
0
−
−
Et ( 2 a 2 + b 2 − b 2ν )
4ab ( −1 + ν 2 )
Et ( −2b 2 − a 2 + a 2ν )
bEt
2a ( −1 + ν 2 )
Et
4 ( −1 + ν )
Etν
2 ( −1 + ν 2 )
Eta
2b ( −1 + ν 2 )
−
0
Et
(
4 ν + 1)
−
4ab ( −1 + ν 2 )
Etν
2 ( −1 + ν 2 )
Eta
4 (ν + 1) b
−
−
Etν
2 ( −1 + ν 2 )
Eta
2b ( −1 + ν 2 )
0
Et
⎤
4 (ν + 1) ⎥
⎥
⎥
0
⎥
⎥
⎥
Et
⎥
4 (ν + 1) ⎥⎥
⎥
Eta ⎥
−
4 (ν + 1) b ⎥
⎥
⎥
0
⎥
⎥
⎥
Eta
⎥
4 (ν + 1) b ⎥⎦
−
Ved udregning af fladeintegralet er det udnyttet, at samtlige indgange i matricen [ B]T [ D][ B] er
konstante, jf. (G.9) og (G.10), og dermed er værdierne af fladeintegralet fundet som indgangsværdien multipliceret med arealet af trekanten.
G.1.5
Global stivhedsmatrice [K]
Ved inddelingen af konsollen i N 2 antal delelementer, fås tilsvarende N 2 antal lokale stivhedsmatricer, der skal indsættes i den globale stivhedsmatrice. Indsætningen er foretaget ved hjælp af Calfem, fil vedlagt på cd-rom som CSTmassivkonsol.m.
61
Del III Numeriske modeller
G CST- og LST-elementopbygning
Ved en manuel indsættelse af de lokale stivhedsmatricer i den globale stivhedsmatrice, skal den
globale stivhedsmatrice være en n x n-matrice svarende til antallet af frihedsgrader. Ud fra det delelement der betragtes, placeres den lokale stivhedsmatrice i den globale stivhedsmatrice, så der er
overensstemmelse imellem enhedsflytningerne og indgangene, hvor den lokale stivhedsmatrice bliver indsat.
Ved at opstille den globale stivhedsmatrice betragtes konsollen inddelt i 4 delelementer som vist på
figur 38. Der er på konsollen indtegnet de globale frihedsgrader.
V8
V7
V12
V10
V11
V9
2
4
3
V4
V3
V6
V5
1
V2
V1
Figur 38: Konsollens globale frihedsgrader.
Den globale stivhedsmatrice for konsollen er følgende
⎡ K1;1
⎢ ⋅
[ K ] = ⎢⎢ ⋅
⎢
⎢⎣ K12;1
⋅
⋅
⋅
⋅
K1;12 ⎤
⎥
⎥
⎥
⎥
K12;12 ⎥⎦
Ved at betragte deleelement 2 på figur 38, er dennes globale frihedsgrader som vist på figur 39.
V8
V10
V7
V9
2
V4
V3
Figur 39: Element 2 med globale frihedsgrader.
Ved at beregne den lokale stivhedsmatrice for delelement 2 ud fra ovenstående teori fås følgende:
62
G CST- og LST-elementopbygning
Del III Numeriske modeller
⎡ K 3;3
⎢K
⎢ 4;3
⎢ K 7;3
⎡⎣ kdelelement 2 ⎦⎤ = ⎢
⎢ K8;3
⎢ K9;3
⎢
⎣⎢ K10;3
K 3;4
K 4;4
K 7;4
K8;4
K9;4
K10;4
K 3;7
K 4;7
K 7;7
K8;7
K 9;7
K10;7
K 3;8
K 4;8
K 7;8
K8;8
K9;8
K10;8
K 3;9
K 4;9
K 7;9
K8;9
K9;9
K10;9
K 3;10 ⎤
K 4;10 ⎥⎥
K 7;10 ⎥
⎥
K8;10 ⎥
K 9;10 ⎥
⎥
K10;10 ⎦⎥
Den lokale stivhedsmatrices indgange indføres i den globale stivhedshedsmatrice ved addering, da
der er forudsat lineær elasticitet, hvor superpositionsprincippet gælder, hvorfor indgangene fra den
lokale stivhedsmatrice placeres som følger:
⎡⋅
⎢⋅
⎢
⎢⋅
⎢
⎢⋅
⎢⋅
⎢
⋅
[ K ] = ⎢⎢
⋅
⎢
⎢⋅
⎢
⎢⋅
⎢⋅
⎢
⎢⋅
⎢⎣⋅
G.1.6
⋅
⋅
⋅
⋅
⋅ + K3;3
⋅ + K 4;3
⋅
⋅
⋅
⋅
⋅ + K 7;3
⋅ + K8;3
⋅ + K 9;3
⋅ + K10;3
⋅
⋅
⋅
⋅
⋅
⋅
+ K 3;4
+ K 4;4
⋅
⋅
+ K 7;4
+ K8;4
+ K9;4
+ K10;4
⋅
⋅
⋅ ⋅
⋅
⋅
⋅
⋅
⋅
⋅
⋅
⋅
⋅
⋅
⋅
⋅
⋅
⋅
⋅ + K3;7
⋅ + K 4;7
⋅
⋅
⋅
⋅
⋅ + K 7;7
⋅ + K8;7
⋅ + K 9;7
⋅ + K10;7
⋅
⋅
⋅
⋅
⋅
⋅
⋅
⋅
+ K3;8
+ K 4;8
⋅
⋅
+ K 7;8
+ K8;8
+ K9;8
+ K10;8
⋅
⋅
⋅
+ K 3;9
+ K 4;9
⋅
⋅
+ K 7;9
+ K8;9
+ K 9;9
+ K10;9
⋅
⋅
⋅
+ K3;10
+ K 4;10
⋅
⋅
+ K 7;10
+ K8;10
+ K 9;10
+ K10;10
⋅
⋅
⋅ ⋅⎤
⋅ ⋅⎥⎥
⋅ ⋅⎥
⎥
⋅ ⋅⎥
⋅ ⋅⎥
⎥
⋅ ⋅⎥
⋅ ⋅⎥
⎥
⋅ ⋅⎥
⎥
⋅ ⋅⎥
⋅ ⋅⎥
⎥
⋅ ⋅⎥
⋅ ⋅⎥⎦
Matriceligningen til bestemmelse af flytningerne
Ved opstilling af elementmetoden er det muligt at bestemme flytningerne ud fra en given belastning,
hvilket er givet ved matriceligningen (G.11), som er udledt i bilag F
{ f } = [ K ] ⋅ {V }
(G.11)
hvor
{f}
[K ]
er lastvektoren der virker på konsollen i knuderne
er den global stivhedsmatrice
{V }
er den globale frihedsgradsvektor
[Byskov 2002, p386]
Lastvektoren for et CST-element opbygges således, at linielasten fra konsollen fordeles gennemsnitlig ud på hver knude ved konsollens overflade som vist i bilag G.4. Flytningsvektoren bestemmes ud
fra konsollens randbetingelser således, at frihedsgraderne ved indspændingen sættes lig nul. Ved
løsning af ligningssystemet bestemmes flytningerne af samtlige frihedsgrader i systemet, som så kan
anvendes til at finde de tilhørende spændinger og tøjninger.
G.1.7
Indsættelse af randbetingelser
For at tage højde for de kinematiske randbetingelser, svarende til at frihedsgraderne i understøtningerne har flytninger lig nul, skal [K] {V} og {f} modificeres. Dette kan gøres ved helt at slette de
63
Del III Numeriske modeller
G CST- og LST-elementopbygning
rækker og søjler som svarer til en randbetingelse lig nul samtidig med at de slettes i {V} og {f}. Herved mistes informationen omkring hvilke frihedsgrader der har flytning lig nul. En anden modificering hvor informationen og flytninger lig nul bevares er i det følgende gennemgået.
Modificeringen bygger på, at en matriceligning skal have alle indgange i {V} som variable og alle
andre indgange som konstanter, før denne kan løses ved traditionelle metoder. Fastholdelsen af visse
frihedsgrader i {V} kan derfor sikres ved at ændre indgange i [K] og {f}. Princippet i modificeringen
er, at hvis frihedsgrad Vi skal være fastholdt, sættes alle værdier i søjle i og række i i [K] lig nul,
bortset fra indgang Kii som sættes lig 1. Samtidig sættes fi lig nul. Princippet er vist i (G.12) for i = 3.
⎡ K11
⎢K
⎢ 21
⎢ K31
⎢
⎣K 41
K12
K 22
K32
K13
K23
K33
K 42
K 43
K14 ⎤ ⎡V1 ⎤ ⎡ f1 ⎤ ⎡ K11
⎢ ⎥
K 24 ⎥⎥ ⎢⎢V2 ⎥⎥ ⎢ f 2 ⎥ ⎢⎢K 21
⋅
=
→
K34 ⎥ ⎢ 0 ⎥ ⎢ f3 ⎥ ⎢ 0
⎥ ⎢ ⎥ ⎢ ⎥ ⎢
K 44 ⎦ ⎣V4 ⎦ ⎣⎢ f 4 ⎦⎥ ⎣K 41
K12
K 22
0
K 42
0 K14 ⎤ ⎡V1 ⎤ ⎡ f1 ⎤
⎢ ⎥
0 K 24 ⎥⎥ ⎢⎢V2 ⎥⎥ ⎢ f 2 ⎥
⋅
=
1 0 ⎥ ⎢V3 ⎥ ⎢ 0 ⎥
⎥ ⎢ ⎥ ⎢ ⎥
0 K 44 ⎦ ⎣V4 ⎦ ⎣⎢ f 4 ⎦⎥
(G.12)
Det ses i (G.12), at ved modificeringen giver V3 intet bidrag til ligningerne forskellig fra række i.
Samtidig giver ligningen for række i V3 = 0, hvilket er det ønskede. Ved at modificere [K] og
{f}
som vist i (G.12), sikres det dermed, at de kinematiske randbetingelser er overholdt, samtidig med at
matriceligningen har alle indgange i {V} som variable og indgange i [K] og { f } som konstanter.
G.2
OPSTILLING AF MATERIALESTIVHEDSMATRICEN
Den konstitutive betingelse i form af Hookes lov skal omskrives fra den oprindelige form med indeksnotation til matriceform, så den kan anvendes i elementmetoden. Dette projekt omhandler en
konsol alene udsat for membranpåvirkning svarende til plan spændingstilstand. Hermed kan den
konstitutive betingelse for konsollen skrives på formen
{N } = [ D]{ε }
N i = Dij ε j
(G.13)
hvor
N
er kraftvektoren
D
er materialestivhedsmatricen
ε
er tøjningsvektoren
[Byskov 2002, p190]
For at bestemme indgangene i [D] tages der udgangspunkt i plan spændingstilstand for isotropt materiale, da tykkelsen af konsollen antages lille i forhold til andre dimensioner. I hovedrapport kapitel
11 er der foretaget en nærmere analyse af, hvor god antagelsen om isotrope egenskaber er. Hermed
bliver Hookes lov angivet ved følgende:
E ⎛
ν
⎞
ε γγ δαβ ⎟
⎜ εαβ +
1 +ν ⎝
1 −ν
⎠
≡0
σ αβ =
σ3 j
hvor
σ
ε
E
ν
64
er spændingerne i planen
er tøjningerne i planen
er elasticitetsmodulen
er Poissons forhold
(G.14)
G CST- og LST-elementopbygning
Del III Numeriske modeller
δ
er Kroneckers delta
[Byskov 2002, p99]
Det ses af (G.14), at samtlige spændinger i x3-retningen er sat lig nul, svarende til plan spændingstilstand.
Det er ved (G.14) nu muligt at udregne spændingerne hver for sig. Da spændingstensoren altid er
symmetrisk, er der kun tre forskellige indgange i σαβ, som bestemmes til følgende
E ⎛
ν
(ε11 + ε 22 ) ⎞⎟
⎜ ε11 +
1 +ν ⎝
1 −ν
⎠
E
σ 11 =
( ε11 + νε 22 )
1 −ν 2
σ 11 =
ν
E ⎛
(ε11 + ε 22 ) ⎞⎟
⎜ ε 22 +
1 +ν ⎝
1 −ν
⎠
E
=
(νε11 + ε 22 )
1 −ν 2
σ 22 =
σ 22
(G.15)
E
ε12
1 +ν
E
σ 12 =
(1 − ν )ε12
1 −ν 2
σ 12 =
Den 2. ordens spændingstensor σαβ opstilles på vektorform. Da σ12 = σ21 defineres denne som
σ αβ
⎡ σ 1 ⎤ ⎡ σ 11 ⎤
⎡σ 11 σ 12 ⎤
⎢
⎥ ⎢ ⎥
=⎢
⎥ → σ i = ⎢σ 2 ⎥ ≡ ⎢σ 22 ⎥
σ
σ
⎣ 21
22 ⎦
⎢⎣σ 3 ⎥⎦ ⎢⎣σ 12 ⎥⎦
(G.16)
For at kunne opstille den tilsvarende 2. ordens tøjningstensor εαβ på vektorform skal det sikres, at
tøjning og spænding stadig er arbejdskonjugerede, således at det virtuelle arbejde stadig er gældende. Dette betegnes også som generaliserede spændinger og tøjninger. Dette krav kan udtrykkes ved
følgende
ε αβ ⋅ σ αβ = ε i ⋅ σ i
(G.17)
Ved at skrive (G.17) ud og anvende den definerede sammenhæng mellem σi og σαβ i (G.16) kan det
ved at sammenligne venstre og højre side i (G.17) indses, at εi skal defineres på følgende måde
ε αβ
⎡ ε 1 ⎤ ⎡ ε 11 ⎤
⎡ε11 ε12 ⎤
⎢ ⎥ ⎢
⎥
=⎢
⎥ → ε i = ⎢ε 2 ⎥ ≡ ⎢ ε 22 ⎥
ε
ε
⎣ 21 22 ⎦
⎣⎢ε 3 ⎦⎥ ⎢⎣ 2ε12 ⎥⎦
(G.18)
Med definitionerne (G.16) og (G.18) kan de konstitutive betingelser i (G.15) omskrives til følgende
65
Del III Numeriske modeller
G CST- og LST-elementopbygning
E
(ε1 + νε 2 )
1 −ν 2
E
σ2 =
(νε1 + ε 2 )
1 −ν 2
E
1
(1 − ν ) ε 3
σ3 =
1 −ν 2
2
σ1 =
(G.19)
Den konstitutive betingelse, som angivet i (G.13), kan opskrives på matriceform ved at betragte
(G.19). Her skal der multipliceres med tykkelsen af skiven, da det er kraften og ikke spændingen, der
er angivet i (G.13). Dermed bliver [D] følgende
⎡1 ν
Et ⎢
[ D] =
ν 1
1 −ν 2 ⎢
⎢⎣0 0
0 ⎤
⎥
0 ⎥
1
1 − ν )⎥⎦
2(
(G.20)
hvor
t
er tykkelsen af skiven
[Byskov 2002, p192]
G.3
OPBYGNING AF LST-MODEL
Dette afsnit omhandler opbygningen af LST-modellen for projektkonsollen. Forskellen fra CSTmodellen er, at der for LST-modellen er indført yderligere 6 frihedsgrader som vist for delelementet
på figur 40.
x2
v12
v4
v11
v3
v8
v2
v6
v5
v10
v7
v9
v1
x1
Figur 40: Frihedsgrader for LST-modellen.
Flytningsfeltet skal være entydigt bestemt ud fra flytningerne af frihedsgraderne for at undgå overlapning af naboelementer. Dette betyder, at elementet får et flytningsfelt beskrevet ved 2. ordens
polynomier af x1 og x2. For en enhedsflytning af frihedsgrad 6, hvor alle andre frihedsgrader sættes
lig nul, vil udbøjningsfiguren på figur 41 være gældende. Det bemærkes, at den stiplede linie mellem
knude 5 og 6 er udeformeret.
66
G CST- og LST-elementopbygning
Del III Numeriske modeller
v12
v4
2
v11
6
v3
v8
v6
3
v5
v10
4
v7
5
v9
v2
1
v1
Figur 41: Enhedsflytning for frihedsgrad 6.
G.3.1
Flytningsinterpolationsmatrice [N]
Som for CST-modellen skal der ligeledes for LST-modellen bestemmes indgangene i [N] ved enkeltvis at opstille en forskrift for flytningstilstanden for hvert tilfælde som så bliver lig den respektive indgang i [N]. Flytningstilstanden er beskrevet på samme form som i (G.21), hvor flytningstilstanden i figur 41 er opskrevet.
N 26 ( x1 , x2 ) = C ⋅ l12 ⋅ l56
(G.21)
hvor
C
er en konstant [-]
l12
er ligningen for linien mellem knude 1 og 2 sat lig nul
l56
er ligningen for linien mellem knude 5 og 6 sat lig nul
Ved samme procedure som ved CST-modellen findes for LST-modellen forskriften for de to udeformerede liniestykker, hvorefter konstanten C beregnes. Ved beregning af N 26 ( x1 , x2 ) fås
N 26 ( x1 , x2 ) =
2 ⋅ N 2 ⋅ x1 ⋅ ( x1 − 2⋅bN )
b2
Ved at betragte figur 41 ses at følgende er gældende
N11 = N 22 , N13 = N 24 , N15 = N 26 , N17 = N 28 , N19 = N 2;10 , N1;11 = N 2;12
De resterende indgange er lig nul. Ved at beregne indgangene i flytningsinterpolationsmatricen [N]
ud fra samme form som (G.21) fås følgende
67
Del III Numeriske modeller
G CST- og LST-elementopbygning
a ⎞⎛
a ⎞
⎛
2 ⋅ N 2 ⋅ ⎜ x2 − ⎟⎜ x2 −
⎟
2
N
N⎠
⎝
⎠⎝
N11 = N 22 =
2
a
N ⋅ (2ax1 ⋅ N + ab - 2 ⋅ bx2 ⋅ N )(ax1 - bx2 )
N13 = N 24 =
(a 2 ⋅ b 2 )
N15 = N 26 =
2 ⋅ N 2 ⋅ x1 ⋅ ( x1 − 2bN )
N17 = N 28 = −
N19 = N 2;10 =
b2
4 ⋅ N ( ax1 − bx2 )( − x2 ⋅ N + a )
ba 2
4 ⋅ N ⋅ x1 ( − x2 ⋅ N + a )
ab
4 ⋅ N ⋅ x1 ⋅ ( ax1 − bx2 )
2
N1;11 = N 2;12 = −
G.3.2
ab 2
Tøjningsfordelingsmatricen [B]
Til beregningen af tøjningsfordelingsmatricen [B] anvendes (G.8). Ved differentiation af flytningsinterpolationsmatricen fås tøjningsfordelingsmatricen [B], jf. (G.22).
⎡
0
⎢
⎢
⎢
0
⎢
⎢
⎢ N ( 4 Nax1 − 4 Nbx2 + ab)
⎢
ab2
⎢
⎢
0
⎢
⎢
N ( 4 Nx1 − b)
⎢
⎢
b2
⎢
⎢
0
⎢
BT = ⎢
4 N ( Nx2 − a )
⎢
⎢
ab
⎢
⎢
0
⎢
⎢
4 N ( − Nx2 + a )
⎢
⎢
ab
⎢
⎢
0
⎢
⎢
⎢ 4 N 2 (−2ax1 + bx2 )
⎢
ab2
⎢
⎢
0
⎢
⎣
0
N ( 4 Nx2 − 3a )
a2
0
N ( −4 Nax1 + 4 Nbx2 − ab)
a2b
0
0
0
4 N ( −2 Nbx2 + ab + Nax1 )
a2b
0
−
4 N 2 x1
ab
0
4 N 2 x1
ab
N ( 4 Nx2 − 3a)
⎤
⎥
a
⎥
⎥
0
⎥
⎥
N ( −4 Nax1 + 4bNx2 − ab) ⎥
⎥
a2b
⎥
N ( 4 Nax1 − 4bNx2 + ab) ⎥
⎥
ab2
⎥
⎥
0
⎥
⎥
N ( 4 Nx1 − b)
⎥
⎥
2
b
⎥
4 N ( −2bNx2 + ab + Nax1 ) ⎥
⎥
a2b
⎥
4 N ( Nx2 − a )
⎥
⎥
ab
⎥
⎥
4 N 2 x1
−
⎥
ab
⎥
4 N (− Nx2 + a )
⎥
⎥
ab
⎥
⎥ (G.22)
4 N 2 x1
⎥
ab
⎥
⎥
4 N 2 (−2ax1 + bx2 )
⎥
2
⎦
ab
2
Som det ses, afhænger tøjningsfordelingsmatricen af x1 og x2 i 1. orden, hvorfor denne metode kaldes Linear Strain Triangle. Den lokale stivhedsmatrix for LST-elementet bestemmes ved (G.1).
68
G CST- og LST-elementopbygning
Del III Numeriske modeller
Da matricen [ B ]T [ D][ B ] afhænger af x1 og x2 skal fladeintegralet løses ved først at integrere over x2
som det indre integral og dernæst over x1 som det ydre integral. Ved at betragte figur 40 ses det, at
når der integreres over trekanten, skal den nedre grænse for x2 være en funktion af x1. Dermed kan
fladeintegralet (G.1) for LST-modellen opstilles som
[k ] = ∫ [ B]T [ D][ B]dA
A
b a
N N
(G.23)
[k ] = ∫ ∫ [ B] [ D][ B]dx2 dx1
xa
0
T
1
b
Udregning af [k] ved (G.23) er foretaget i maple, fil vedlagt på cd-rom som kned.mw. Da [k] er en
12x12 matrice er den i det følgende delt i 4 delmatricer for at lette overskueligheden, jf. (G.24) (G.28).
⎡[k ]
[k ] = ⎢ 1
⎣[k3 ]
3b
⎡
⎢ 4 (ν + 1) a
⎢
⎢
0
⎢
⎢
⎢
b
⎢
⎢
Et ⎢ 4 (ν + 1) a
k1 = ⎢
3
−1
⎢
⎢ 4 (ν + 1)
⎢
⎢
0
⎢
⎢
⎢
1
⎢
⎢⎣ 4 (ν + 1)
[k2 ]⎤
[k4 ]⎥⎦
(G.24)
0
b
4 (ν + 1) a
−1
4 (ν + 1)
0
−3b
2 ( −1 + ν 2 ) a
2 ( −1 + ν 2 )
ν
−b
2 ( −1 + ν 2 ) a
−ν
2 ( −1 + ν 2 )
3
4 ( −1 + ν )
−a
2b ( −1 + ν 2 )
ν
−3 ( 2 a2 + b2 − b2ν )
−b
2 ( −1 + ν 2 ) a
3
4 ( −1 + ν )
−ν
2 ( −1 + ν 2 )
−a
2b ( −1 + ν 2 )
0
−1
4 (ν + 1)
4b ( −1 + ν 2 ) a
2 ( −1 + ν 2 )
⎡ −b
⎢ (ν + 1) a
⎢
⎢ −2ν
⎢(
2
⎢ −1 + ν )
⎢ −b
⎢
Et ⎢ (ν + 1) a
k2 = ⎢
1
3 ⎢
⎢ (ν + 1)
⎢
⎢
0
⎢
⎢
⎢
⎢
0
⎢
⎣
⎡ −b
⎢ (ν + 1) a
⎢
⎢ 1
⎢ (
⎢ ν + 1)
⎢
⎢ 0
Et ⎢
k3 = ⎢
3 ⎢ −1
⎢ (ν + 1)
⎢
⎢ 0
⎢
⎢
⎢
⎢ 0
⎣
1
(ν + 1)
0
3 ( − a2 + a2ν − 2 b2 )
ν
4b ( −1 + ν 2 ) a
2 ( −1 + ν 2 )
2 ( −1 + ν 2 )
ν
−3a
2b ( −1 + ν 2 )
a
4 (ν + 1) b
0
−1
(ν + 1)
0
( −1 + ν 2 )
2ν
0
0
−2ν
0
0
2a
b ( −1 + ν 2 )
2b
(−1 + ν 2 ) a
0
0
−2ν
(−1 + ν 2 )
0
0
2ν
2a
0
−1
(ν + 1)
0
1
(ν + 1)
2
(−1 + ν 2 )
−2ν
( −1 + ν 2 )
2b
−b
(ν + 1) a
−2ν
(−1 + ν 2 ) b ( −1 + ν 2 )
1
(ν + 1)
2b
0
(−1 + ν ) a
(−1 + ν 2 )
(−1 + ν ) a
0
(−1 + ν 2 )
2ν
0
0
0
0
0
0
( −1 + ν 2 )
0
2a
b ( −1 + ν 2 )
( −1 + ν 2 )
0
1
(ν + 1)
(ν + 1) b
2
2
−2ν
−a
(G.25)
⎤
⎥
⎥
⎥
0
⎥
⎥
⎥
1
⎥
(ν + 1) ⎥
⎥
−a ⎥
(ν + 1) b ⎥⎥
−2ν ⎥
⎥
(−1 + ν 2 ) ⎥
⎥
−a ⎥
(ν + 1) b ⎥⎦
(G.26)
⎤
⎥
⎥
⎥
0 ⎥
⎥
−1 ⎥⎥
(ν + 1) ⎥
⎥
0 ⎥
⎥
1 ⎥⎥
(ν + 1) ⎥
⎥
−a ⎥
(ν + 1) b ⎦⎥
(G.27)
0
(−1 + ν ) a
2b
1
⎤
4 (ν + 1) ⎥
⎥
⎥
0
⎥
⎥
⎥
−1 ⎥
(
)
4 ν +1 ⎥
⎥
⎥
a
⎥
4 (ν + 1) b ⎥
⎥
⎥
0
⎥
⎥
⎥
3a
⎥
4 (ν + 1) b ⎦⎥
2ν
2a
b ( −1 + ν 2 )
−2ν
( −1 + ν 2 )
0
69
Del III Numeriske modeller
⎡ −2 ( 2 a 2 + b2 − b2ν )
⎢
2
⎢ b ( −1 + ν ) a
⎢
⎢
1
⎢
( −1 + ν )
⎢
⎢
4a
⎢
⎢
b ( −1 + ν 2 )
Et ⎢
k4 =
3 ⎢
−1
⎢
⎢
( −1 + ν )
⎢
⎢
⎢
0
⎢
⎢
⎢
1
⎢
( −1 + ν )
⎢⎣
G CST- og LST-elementopbygning
1
( −1 + ν )
2 ( −a2 + a2ν − 2 b2 )
b ( −1 + ν 2 ) a
−1
( −1 + ν )
4a
b ( −1 + ν 2 )
−1
( −1 + ν )
−1
( −1 + ν )
(ν + 1) b
1
( −1 + ν )
1
( −1 + ν )
(ν + 1) a
−2 ( 2 a 2 + b2 − b2ν )
b ( −1 + ν 2 ) a
−2a
2 ( −a 2 + a 2ν − 2 b2 )
−2 a
(ν + 1) b
1
( −1 + ν )
1
( −1 + ν )
(ν + 1) a
−1
( −1 + ν )
−1
( −1 + ν )
(−1 + ν ) a
−2b
0
b ( −1 + ν
2
)a
4b
2
0
−2b
−1
( −1 + ν )
−2 ( 2 a 2 + b2 − b2ν )
b ( −1 + ν 2 ) a
1
( −1 + ν )
⎤
⎥
⎥
⎥
⎥
0
⎥
⎥
⎥
−1
⎥
⎥
( −1 + ν )
⎥
⎥
4b
⎥
(−1 + ν 2 ) a ⎥⎥
⎥
1
⎥
( −1 + ν )
⎥
⎥
2
2
2
2 ( −a + a ν − 2 b ) ⎥
⎥
b ( −1 + ν 2 ) a
⎦⎥
1
( −1 + ν )
(G.28)
Analogt med CST-modellen, afsnit G.1.5, indføres de lokale stivhedsmatricer i den globale stivhedsmatrix, og matrixligningen løses herefter. Resultatet af elementmetodeanalysen med LSTelementer ses i hovedrapport afsnit 7.1.3.
G.4
LASTFORDELING PÅ KNUDER
I dette afsnit bestemmes, hvorledes linielasten skal fordeles ud på knuderne i elementmetoden for
CST-og LST-elementer. Problematikken er, at en linielast virker over en hel flade, og ved elementmetoden er det kun muligt at påføre konstruktionen laster i knudepunkter. Derfor bliver der i dette
afsnit opstillet en fordeling af de punktlaster, der er ækvivalente med den linielast, som konsollen er
belastet med.
Der anvendes virtuelle flytningers princip, hvor det virtuelle arbejde ved punktbelastning skal være
ækvivalent med den påførte linielast. Dette kan udtrykkes ved følgende relation
T
∫ {δ u} {w}dl = {δ v} { f }
T
l
(G.29)
hvor
{δu}
{w}
er flytningsvektoren for den virtuelle flytningstilstand [mm]
er lastvektoren for linielasten
{δv}
er flytningerne i knudepunkternes frihedsgrader [mm]
{f}
er lastvektoren for lasten i knudepunkternes frihedsgrader [N]
[ mmN ]
[Williams & Todd 2000, p298]
Det ses, at der i venstresiden af (G.29) skal integreres over liniestykket, da det her er en linielast,
mens der ved det ækvivalente arbejde fordelt på knuderne på højresiden af (G.29) er tale om punktarbejder. I (G.29) kan den virtuelle flytningsvektor udtrykkes ved sammenhængen beskrevet i (G.2)
{δ u} = [ N ]{δ v}
hvor
{δ u}
[N ]
er den virtuelle flytningsvektor i x1 og x2 retningen
er flytningsinterpolationsmatricen
{δ v} er den virtuelle frihedsgradsvektor
[Byskov 2002, p378]
Hermed kan (G.29) ved anvendelse af (G.30) omskrives til følgende
70
(G.30)
G CST- og LST-elementopbygning
Del III Numeriske modeller
T
∫ {δ u} {w}dl = {δ v} { f }
∫ ([ N ]{δ v}) {w}dl = {δ v} { f }
T
l
T
T
l
T
T
T
∫ {δ v} [ N ] {w}dl = {δ v} { f }
(G.31)
l
T
{δ v}T ∫ [ N ] {w}dl = {δ v}T { f }
l
Her er anvendt, at den transponerede af et matrixprodukt er lig de enkelte matricer transponeret og i
omvendt rækkefølge. Desuden sættes {δv} udenfor integralet, da det ikke afhænger af længdekoordinaten. Det ses af (G.31), at da {δv} skal kunne vælges vilkårligt, kan sammenhængen forkortes til
T
∫ [ N ] {w}dl = { f }
(G.32)
{f}
findes ved (G.32), hvor [N] og {w} er
l
Dermed kan lastvektoren for lasten i knudepunkterne
kendte.
G.4.1
Fordeling af linielast for CST-element
Både ved opbygning af konsollen med CST-elementer og LST-elementer er de eneste elementer, der
bliver belastet, de elementer som har vandret overside svarende til element nr. 2 og 4, jf. figur 42.
p
2
4
3
1
x2
x1
Figur 42: Konsol inddelt i elementer og påsat linielast.
Idet der ikke er last i x1-aksens retning kan lastvektoren for linielasten hermed opskrives til
⎡0⎤
{w} = ⎢ ⎥
⎣ p⎦
hvor
p
er linielasten som defineret i figur 42
[ mmN ]
Det ses af værdien af {w} , at de eneste indgange i [N] der giver værdier, er indgangene N2i, da indgangene N1i skal multipliceres med 0 fra {w} . Derudover er det for CST-elementerne kun indgangene svarende til frihedsgrad 4 og 6, der skal medtages, da der ikke er last ved frihedsgrad 2, jf. figur
36. Indgangene i [N] afhænger både af x1 og x2, og da det er flytninger i oversiden af elementet, der
71
Del III Numeriske modeller
skal anvendes, sættes x2 =
G CST- og LST-elementopbygning
a
N
. Dermed bliver [N] forkortet til følgende, når der betragtes CST-
elementer.
[ N ] = [ N24
N26 ]
a ⎞⎤ ⎡ − N ⋅ x1 + b
⎡ ⎛
⎢⎣ N ⎜⎝ x2 = N ⎟⎠⎥⎦ = ⎣⎢
b
N ⋅ x1 ⎤
b ⎦⎥
Ved indsættelse i (G.32) udregnes de indgange i lastvektoren for et CST-element forskellig fra nul til
b
N
⎡ − N ⋅ x1 + b
{f}= ∫ ⎢
⎣
b
0
T
N ⋅ x1 ⎤
⋅ p ⋅ dx1
b ⎥⎦
⎡1 b ⋅ p ⎤
⎢2 N ⎥
{f}= ⎢
⎥
⎢1 b ⋅ p ⎥
⎣⎢ 2 N ⎦⎥
Det ses, at lastfordelingen er ligeligt fordelt mellem knuderne for et CST-element. Lasten på knuderne er vist til venstre i figur 43.
f11 =
v4
1 b⋅ p
2 N
f 21 =
1 b⋅ p
2 N
v6
v3
f11 =
1 b⋅ p
6 N
v4
v5
v3
v8
1
f 21 =
v7
1
2 b⋅ p
3 N
f 31 =
1 b⋅ p
6 N
v12
v6
v11
v10
v5
v9
v2
v2
v1
v1
Figur 43: Lastfordeling på knuderne for et CST-element til venstre og for et LST-element til højre.
G.4.2
Fordeling af linielast for LST-element
Fremgangsmåden til bestemmelse af lastfordelingen på knuderne for et LST-element er helt analogt
til bestemmelsen for et CST-element. Da det af (G.32) ses, at fordelingen afhænger af [N], der er
forskellig for CST- og LST-elementerne, vil lastfordelingen være anderledes. For et LST-element
bliver indgangene, der skal anvendes i [N] ved betragtning af figur 40 dermed følgende
[ N ] = [ N24
N2 −12
N26 ]
a ⎞⎤ ⎡ (−2 x1 ⋅ N + b) ⋅ (− x1 ⋅ N + b)
⎡ ⎛
⎢⎣ N ⎜⎝ x2 = N ⎟⎠⎥⎦ = ⎣⎢
b2
4 N ⋅ x1 (− x1 ⋅ N + b)
b2
Ved indsættelse i (G.32) udregnes lastvektoren for et LST-element til
72
− N ⋅ x1 ⋅ (−2 x1 ⋅ N + b) ⎤
⎦⎥
b2
G CST- og LST-elementopbygning
b
N
Del III Numeriske modeller
⎡ (−2 x1 ⋅ N + b) ⋅ (− x1 ⋅ N + b)
{f}= ∫ ⎢
⎣
b2
0
4 N ⋅ x1 (− x1 ⋅ N + b)
b2
T
− N ⋅ x1 ⋅ (−2 x1 ⋅ N + b) ⎤
⎥⎦ ⋅ p ⋅ dx1
b2
⎡1 b ⋅ p ⎤
⎢6 N ⎥
⎢
⎥
2 b ⋅ p⎥
⎢
{f}=
⎢3 N ⎥
⎢
⎥
⎢1 b ⋅ p ⎥
⎢⎣ 6 N ⎥⎦
Det ses, at lastfordelingen ikke længere er ligeligt fordelt mellem knuderne for et LST-element. I
lastfordelingen skal der dog regnes med, at naboelementet ved samling af elementerne vil bidrage til
lasten i hjørneknuderne og ikke i midterknuden, hvilket gør, at midterknuderne kun skal modelleres
med dobbelt så stor punktlast, når elementet betragtes på konstruktionsniveau. Lasten på knuderne
for et LST-element er vist til højre i figur 43.
G.5
SPÆNDINGER I KONSOL BEREGNET UD FRA CST- OG LSTMODEL
For at kunne beregne spændingerne i konsollen ved både CST- og LST modellen, skal de globale
flytninger {V } for konsollens globale frihedsgrader beregnes ved at løse ligningssystem givet ved
(G.11), fil vedlagt på cd-rom som CSTmassivkonsol.m og LSTmassivkonsol.m. For at beregne tøjningerne i hver enkelt delelement skal tøjningsfordelingsmaticerne benyttes.
Ved at multiplicere tøjningsfordelingsmatricen med de globale flytninger der svarer til tøjningsfordelingsmaticens lokale frihedsgrader, findes tøjningerne til:
{ε ( xα )} = [ B]{V }
(G.33)
[Byskov 2002, p400]
Det bemærkes, at tøjningsfordelingsmatricen for CST-modellen ikke afhænger af x1 og x2 . For
LST-modellen er tøjningsfordelingsmatricen afhængig af x1 og x2 . For elastiske skiver er de lineære konstitutive relationer er givet ved Hookes lov. Ved indsættelse af tøjningerne for de enkelte delelementer beregnes spændingerne ud fra
{σ ( xα )} = 1t [ D ]{ε ( xα )}
(G.34)
[Byskov 2002, p191]
Beregningerne for spændingerne er vedlagt på cd-rom som CSTmassivkonsol.m og LSTmassivkonsol.m.
73
H Eliminering af indre frihedsgrader
H
Del III Numeriske modeller
ELIMINERING AF INDRE
FRIHEDSGRADER
Formålet med dette bilag er at gennemgå den teori, der ligger bag metoden omkring eliminering af
indre frihedsgrader, som er benyttet i dette projekt til at opstille stivhedsmatricer for superelementer.
Dette bilag er baseret på [Byskov 2002, pp401-406].
H.1
POTENTIEL ENERGI
Den potentielle energi af et system kan skrives som vist i afsnit F.1. For superelementet er den potentielle energi
Π P ({v}) = 12 {v} [k ]{v} − {v} {r }
T
T
(H.1)
hvor
ΠP
er den potentielle energi
{v}
[k ]
{r }
er flytningerne for superelementets frihedsgrader
er superelementets stivhedsmatrice
er lasten i superelementets frihedsgrader
Som det ses af (H.1) afhænger den potentielle energi for superelementet af samtlige frihedsgrader.
Ved at eliminere de indre frihedsgrader kan der opskrives et nyt udtryk for den potentielle energi,
udelukkende udtrykt ved de ydre frihedsgrader.
Antallet af frihedsgrader kaldes i det følgende ndof. Da der er to frihedsgrader i hver knude for CSTelementerne, jf. afsnit G.1, ses antallet af frihedsgrader i henhold til figur 44, at være:
ndoftotal = 2 ⋅ 47 = 94
ndofi = 2 ⋅18 = 36
ndofe = 2 ⋅ 29 = 58
hvor
ndoftotal er det totale antal frihedsgrader
ndofi er antallet af indre frihedsgrader
ndofe er antallet af ydre frihedsgrader
75
Del III Numeriske modeller
H Eliminering af indre frihedsgrader
54
55
56
57
33
11
13 12
37
58
50
51
32
31
9
10
8
49
30
7
4 3
60
61
14
15
39
16
40
41
42
17
2
1
46
26
45
25
24
44
23
18
19 20
21
43
47
27
5
6
48
28
29
38
59
62
34
35
36
52
53
22
63
Figur 44: Superelement B, med parametre
n=2, m=3, jf. hovedrapport afsnit 8.1.
H.2
ADSKILLELSE AF INDRE OG YDRE FRIHEDSGRADER
Der defineres to transformationsmatricer [Ti] og [Te], som skiller henholdsvis de indre og de ydre
frihedsgrader fra samtlige systemets frihedsgrader. Dermed kan flytningerne for de indre og ydre
frihedsgrader defineres som
{vi} = [Ti ]{v}
{ve} = [Te ]{v}
(H.2)
hvor
{vi}
[Ti ]
{ve}
[Te ]
er flytningerne for systemets indre frihedsgrader
er transformationsmatricen for de indre frihedsgrader
er flytningerne for systemets ydre frihedsgrader
er transformationsmatricen for de ydre frihedsgrader
Som det ses af (H.2) får [Ti] og [Te] henholdsvis dimensionen ndofi × ndoftotal og ndofe × ndoftotal .
De er nul-matricer med ét ettal i hver række, placeret efter hvilke frihedsgrader der er henholdsvis
indre og ydre. Dermed kan flytningerne i alle frihedsgraderne for superelementet skrives som
{v} = [Ti ] {vi} + [Te ] {ve}
T
T
(H.3)
Hvis (H.3) benyttes kan den potentielle energi (H.1) skrives som en funktion af de indre og ydre
frihedsgraders flytninger.
{
} [k ]{[T ] {v } + [T ] {v }}
− {[T ] {v } + [T ] {v }} {r }
Π P ({vi} ,{ve}) = 12 [Ti ] {vi} + [Te ] {ve}
T
T
T
i
e
T
i
T
i
T
T
i
T
e
De følgende regneregler for matricer anvendes til omskrivning af (H.4).
76
e
e
(H.4)
H Eliminering af indre frihedsgrader
Del III Numeriske modeller
([ A] + [B]) [C] = [ A][C] + [B][C]
T
([ A] + [B]) = [ A]T + [B]T
T
([ A][B]) = [B]T [ A]T
(H.5)
hvor
[A],[B],[C] er vilkårlig matricer
[Lay 1996, p106]
Ved hjælp af (H.5) omskrives (H.4).
{
} {
− {{v } [T ] + {v } [T ]}{r }
}
Π P ({vi} ,{ve}) = 12 {vi} [Ti ] + {ve} [Te ] [k ] [Ti ] {vi} + [Te ] {ve}
T
T
T
T
T
T
i
i
e
e
= 12 ({vi} [Ti ][k ][Ti ] {vi} + {vi} [Ti ][k ][Te ] {ve}
T
T
T
T
(H.6)
+ {ve} [Te ][k ][Ti ] {vi} + {ve} [Te ][k ][Te ] {ve})
T
T
T
T
− {vi} [Ti ]{r } − {ve} [Te ]{r }
T
T
Der defineres følgende hjælpestørrelser
[kii ] ≡ [Ti ][k ][Ti ]
T
[kie ] ≡ [Ti ][k ][Te ]
T
[kei ] ≡ [Te ][k ][Ti ]
T
[kee ] ≡ [Te ][k ][Te ]
[ri ] ≡ [Ti ]{r }
[re ] ≡ [Te ]{r }
T
(H.7)
Ved hjælp af (H.7) omskrives (H.6) til
Π P ({vi} ,{ve}) = 12 ({vi} [kii ]{vi} + {vi} [kie ]{ve} + {ve} [kei ]{vi} + {ve} [kee ]{ve})
T
T
T
T
− {vi} {ri} − {ve} {re}
T
H.3
T
(H.8)
VARIATION MED HENSYN TIL INDRE FRIHEDSGRADER
Den potentielle energi for systemet er dermed i (H.8) opskrevet som funktion af de indre og ydre
frihedsgraders flytninger. Da de indre frihedsgrader ikke deles med andre superelementer, kan disse
frihedsgrader varieres uafhængigt for dette superelement, og derved elimineres.
δ Π P {v } = 0 ∀ {δ vi}
i
(H.9)
Variationen kan som i bilag F beregnes ved først at substituere {vi } med {vi } + ε {δ vi } .
77
Del III Numeriske modeller
H Eliminering af indre frihedsgrader
T
T
Π P ({vi } + ε {δ vi }) = 12 ⎡({vi } + ε {δ vi }) [kii ] ({vi } + ε {δ vi }) + ({vi } + ε {δ vi }) [kie ]{ve}
⎣
T
T
+ {ve} [kei ] ({vi } + ε {δ vi }) + {ve} [kee ]{ve}⎤
⎦
− ({vi } + ε {δ vi }) {ri } − {ve} {re}
T
T
(H.10)
= ⎡{vi } [kii ]{vi } + {vi } [kii ] ε {δ vi } + ε {δ vi } [kii ]{vi } + ε {δ vi } [kii ] ε {δ vi }
⎣
T
1
2
T
T
T
+ {vi } [kie ]{ve} + ε {δ vi } [kie ]{ve} + {ve} [kei ]{vi } + {ve} [kei ] ε {δ vi }
T
T
T
T
T
T
T
T
+ {ve} [kee ]{ve}⎤ − {vi } {ri } − ε {δ vi } {ri } − {ve} {re}
⎦
Der differentieres med hensyn til ε .
∂Π P ({vi } + ε {δ vi })
∂ε
T
T
T
= 12 ⎡{vi } [kii ]{δ vi } + {δ vi } [kii ]{vi } + 2ε {δ vi } [kii ]{δ vi }
⎣
(H.11)
+ {δ vi } [kie ]{ve} + {ve} [kei ]{δ vi }⎤⎦ − {δ vi } {ri }
T
T
T
Da [k] er symmetrisk er [kii], [kie] og [kei] også symmetriske, og [kie] = [kei] jf. (H.7). Da hvert af
ledene i (H.11) er skalarer, er det tilladeligt at transformere disse uden at værdien ændres, hvilket er
gjort i (H.12).
∂Π P ({vi } + ε {δ vi })
∂ε
(
)
T
T
T
T
= 12 ⎡ {δ vi } [kii ]{vi } + {δ vi } [kii ]{vi } + 2ε {δ vi } [kii ]{δ vi }
⎢⎣
(
)
T
T
T
T
+ {δ vi } [kie ]{ve} + {δ vi } [kei ]{ve} ⎤ − {δ vi } {ri }
⎥⎦
(H.12)
= {δ vi } [kii ]{vi } + ε {δ vi } [kii ]{δ vi } + {δ vi } [kie ]{ve} − {δ vi } {ri }
T
T
T
T
Der evalueres for ε = 0 .
∂Π P ({vi } + ε {δ vi })
∂ε
= {δ vi } [kii ]{vi } + {δ vi } [kie ]{ve} − {δ vi } {ri }
T
T
T
(H.13)
ε =0
Variationen af den potentielle energi kræves, jf. (H.9) lig nul for alle værdier af δ vi .
{δ vi } [kii ]{vi } + {δ vi } [kie ]{ve} − {δ vi} {ri } = 0
[kii ]{vi} + [kie ]{ve} − {ri} = {0}
T
T
T
(H.14)
Da den potentielle energi Π P er negativ definit, det vil sige altid negativ, er [kii] også negativ definit
jf. (H.1). Determinanten af [kii] vil dermed altid være negativ, hvorfor [kii] er ikke-singulær og kan
inverteres, hvilket giver: [Lay 1996, p456]
{vi } = [kii ] ({ri } − [kie ]{ve})
−1
H.4
(H.15)
POTENTIEL ENERGI SOM FUNKTION AF YDRE FRIHEDSGRADER
Udtrykket (H.15) indsættes i (H.8), hvorved den potentielle energi udtrykkes udelukkende som en
funktion af {ve} . Symmetrien af skalarer, [kii], [kie] og [kei] udnyttes, og giver
78
H Eliminering af indre frihedsgrader
Del III Numeriske modeller
(
+ ([k ]
) [k ] ([k ] ({r } − [k ]{v }))
({r } − [k ]{v })) [k ]{v } + {v } [k ] ([k ] ({r } − [k ]{v }))
+ {v } [k ]{v }⎤ − ([k ] ({r } − [k ]{v })) {r } − {v } {r }
⎦
T
−1
Π P ({ve}) = 12 ⎡ [kii ] ({ri } − [kie ]{ve})
⎣
−1
ii
ii
i
T
−1
ii
i
ie
ie
−1
T
e
ie
e
e
ei
ee
e
ii
i
T
−1
T
e
e
ii
i
ie
ie
e
T
e
i
e
e
(H.16)
Π P ({ve}) =
1
2
(({r }
)
− {ve} [kie ] [kii ]
T
T
i
(
+ {ve} [kei ] [kii ]
T
−
(({r }
T
i
−1
−1
)[k ]([k ]
−1
ii
ii
({ri} − [kie ]{ve}))
)
− {ve} [kie ] [kii ]
T
−1
({ri} − [kie ]{ve}))
+ 12 {ve} [kee ]{ve}
T
){r } − {v } {r }
T
e
i
e
Der ganges nu ind i parenteserne.
⎛{ri }T [kii ]−1 [kii ][kii ]−1 {ri } − {ri }T [kii ]−1 [kii ][kii ]−1 [kie ]{ve}
⎞
⎟
Π P ({ve}) = ⎜
−
−
−
−
1
1
1
1
T
T
⎜ − {v } [k ][k ] [k ][k ] {r } + {v } [k ][k ] [k ][k ] [k ]{v }⎟
e
ie
ii
ii
ii
i
e
ie
ii
ii
ii
ie
e ⎠
⎝
1
2
+ {ve} [kei ][kii ]
{ri} − {ve} [kei ][kii ] [kie ]{ve} + 12 {ve} [kee ]{ve}
−1
−1
T
T
T
− {ri } [kii ] {ri } + {ve} [kie ][kii ] {ri } − {ve} {re}
−1
T
−1
T
T
(H.17)
Π P ({ve}) =
1
2
{ri} [kii ] {ri} − {ve} [kie ][kii ] {ri } + 12 {ve} [kie ][kii ] [kie ]{ve}
−1
T
−1
T
+ {ve} [kei ][kii ]
−1
T
{ri } − {ve} [kei ][kii ] [kie ]{ve} + 12 {ve} [kee ]{ve}
−1
−1
T
T
T
− {ri } [kii ] {ri } + {ve} [kie ][kii ] {ri } − {ve} {re}
−1
T
−1
T
T
Det udnyttes, at flere led er ens
Π P ({ve}) = − 12 {ri } [kii ]
−1
T
{ri } − 12 {ve} [kie ][kii ] [kie ]{ve}
−1
T
+ 12 {ve} [kee ]{ve} + {ve} [kie ][kii ]
T
−1
T
{ri } − {ve} {re}
T
(H.18)
Led, som ikke indeholder {ve} er konstante, når {ve} varieres, hvorfor de kan udelukkes af udtrykket
for den potentielle energi, (H.18), da det er variationen af den potentielle energi, der undersøges.
Π P ({ve}) =
1
2
{ve} ([kee ] − [kie ][kii ] [kie ]){ve} + {ve} ([kie ][kii ] {ri} − {re})
−1
T
T
−1
(H.19)
Der defineres følgende hjælpestørrelser
⎡⎣kee* ⎤⎦ ≡ [kee ] − [kie ][kii ]
−1
[kie ]
{r } ≡ {r } − [k ][k ] {r }
−1
*
e
e
ie
ii
(H.20)
i
Den potentielle energi, (H.19), kan dermed ved hjælp af (H.20) skrives som.
Π P ({ve}) =
1
2
{ve}
T
⎡⎣kee* ⎤⎦ {ve} − {ve} {re*}
T
(H.21)
Dette udtryk, (H.21), kan nu sammenlignes med det generelle udtryk for den potentielle energi,
(H.1). Som det ses, er udtrykkene på samme form, hvorfor den modificerede stivhedsmatrice ⎡⎣kee* ⎤⎦
og højresidevektor {re*} i (H.21), der kun afhænger af de ydre frihedsgrader, kan anvendes direkte i
79
Del III Numeriske modeller
H Eliminering af indre frihedsgrader
et elementmetodeprogram. Derved kan flytningerne i superelementets ydre frihedsgrader bestemmes. Superelementet er vist på figur 45.
Figur 45: Superelement uden indre knuder.
H.5
BESTEMMELSE AF INDRE FRIHEDSGRADERS FLYTNINGER
For at bestemme spændingen i de CST-elementer, som superelementet er opbygget af, er det imidlertid ikke nok at kende flytningerne i superelementets ydre frihedsgrader; også de indre frihedsgraders
flytninger må kendes.
Udtrykket for opdelingen af indre og ydre frihedsgrader, (H.3), kan ved hjælp af (H.15) omskrives til
{v} = [Ti ] [kii ] ({ri } − [kie ]{ve}) + [Te ] {ve}
T
−1
T
{v} = [Ti ] [kii ] {ri } − [Ti ] [kii ] [kie ]{ve} + [Te ] {ve}
T
−1
−1
T
{v} = [Ti ] [kii ] {ri } + ([Te ]
T
−1
T
T
− [Ti ] [kii ]
−1
T
(H.22)
[kie ]){ve}
Der defineres følgende hjælpestørrelser
⎡⎣Te* ⎤⎦ ≡ [Te ] − [Ti ] [kii ]
T
−1
T
[kie ]
{r } ≡ [T ] [k ] {r }
−1
T
*
i
ii
(H.23)
i
Derved kan flytningerne for de indre frihedsgrader ud fra (H.22) og (H.23) findes som
{v} = ⎡⎣Te* ⎤⎦ {ve} + {r *}
(H.24)
Da alle frihedsgraders flytninger nu er kendt, kan spændingerne nu beregnes konventionelt, som
gennemgået i bilag F.
80
I Topologioptimering
I
Del III Numeriske modeller
TOPOLOGIOPTIMERING
Det er valgt at foretage optimeringen af stivheden af en konstruktion ved hjælp af en optimalbetingelse. I dette afsnit er der redegjort for baggrunden og implementationen af optimalbetingelsen.
Først beskrives en målfunktion for konstruktionens stivhed og dens sidebetingelser. Dernæst bestemmes målfunktionens minimum svarende til maksimal stivhed. Afsnittet afsluttes med en redegørelse af filtret samt opdateringsrutinen.
I.1
OPRINDELIG MÅLFUNKTION
Der tages udgangspunkt i følgende målfunktion, som skal minimeres og beskriver den dobbelte tøjningsenergi af kontinuet som en funktion af den relative densitet af hvert element
n
(
g ({ ρ } ) = {V } ⎡⎣ K ({ ρ } ) ⎤⎦ {V } = ∑ {vi } ⎡⎣ ki ( ρi ) ⎤⎦ {vi }
T
i =1
T
)
(I.1)
hvor
er konstruktionens tøjningsenergi [Nm]
er den relative densitetsvektor af dimensionen n × 1 [-]
g
{ρ}
n
er antal elementer [-]
{V } , {vi } er den globale og lokale flytningsvektor af dimensionen m × 1 og l × 1 [m]
[ K ] , [ ki ]
er den globale og lokale stivhedsmatrix af dimensionen m × m og l × l [N/m]
m, l
er antal globale og lokale frihedsgrader [-]
[Sigmund 2001]
SIMP-modellen benyttes som estimat for et elements stivhed. Her antages stivheden af et element at
være proportional med stivheden af et element med udfyldt materiale, hvor proportionalitetsfaktoren
er elementets relative densitet opløftet i en potens
⎡⎣ ki ( ρi ) ⎤⎦ = ρi p [ ki 0 ] ,
p ≥ 1,
i = 1, 2,.., n
(I.2)
hvor
[ ki 0 ]
er den lokale stivhedsmatrix, svarende til udfyldt materiale, af dimensionen l × l [N/m]
p
er en potens, der afhænger af estimatets form [-]
Potensen kan varieres afhængigt af tilfredsheden med den optimale konstruktion, der fremkommer
ved beregning. Typisk vælges p = 3 . I specialtilfælde med p = 1 eller p = 2 , svarer (I.2) til henholdsvis Voigts eller det cellulære estimat.
81
Del III Numeriske modeller
I.2
I Topologioptimering
SIDEBETINGELSER
Ved minimering af målfunktionen, skal konstruktionens sidebetingelser være opfyldt:
•
•
•
statiske og kinematiske betingelser er opfyldt
anvendt volumen er lig volumen, som er til rådighed
partiklernes relative densitet er i det lukkede interval mellem 0 og 1
Sidebetingelserne er beskrevet matematisk ved følgende fire ligninger
[ K ]{V } = {F }
( statiske og kinematiske betingelser )
(I.3)
∑(ρ ) − f ⋅ n = 0
( volumenbetingelse )
(I.4)
ρi ≥ 0,
i = 1, 2,.., n
(I.5)
ρi ≤ 1,
i = 1, 2,.., n
( nedre grænse for rel. densitet )
( øvre grænse for rel. densitet )
n
i
i =1
(I.6)
hvor
{F}
er den globale lastvektor af dimensionen m × 1 [N]
f
er konstruktionens porøsitet [-]
Ved at benytte FEM til bestemmelse af (I.1), er sidebetingelsen (I.3) automatisk opfyldt. De øvrige
sidebetingelser kan opfyldes ved at gange vilkårlige skalarer på (I.4)-(I.6) og i stedet kræve de resulterende ligninger opfyldt
⎛
⎞
n
α ⎜ ∑ ( ρi ) − f ⋅ n ⎟ = 0,
⎝ i =1
βi ρi = 0,
for alle α
⎠
for alle βi ≥ 0, i = 1, 2,.., n
( volumenbetingelse )
(I.7)
(I.8)
γ i ( ρi − 1) = 0,
for alle γ i ≥ 0, i = 1, 2,.., n
( nedre grænse for rel.densitet )
( øvre grænse for rel. densitet )
(I.9)
hvor
α
{β }
er en vilkårlig skalar [Nm]
{γ }
er en vektor af vilkårlige skalarer af dimensionen n × 1 [Nm]
er en vektor af vilkårlige skalarer af dimensionen n × 1 [Nm]
De vilkårlige skalarer kaldes Lagrange-multiplikatorer, hvor α knytter sig til hele konstruktionen,
mens β og γ knytter sig til de enkelte elementer. Det ses af (I.8) og (I.9), at for elementer med relativ
densitet forskellig fra 0 og 1, skal β og γ være lig nul
⎧ βi = 0
,
⎩γ i = 0
ρi ∈]0;1[⇒ ⎨
I.3
i = 1, 2,.., n
(I.10)
ENDELIG MÅLFUNKTION
De resulterende ligninger (I.7)-(I.9) adderes til den oprindelige målfunktion (I.1), hvorved den endelige målfunktion bliver
82
I Topologioptimering
gM
Del III Numeriske modeller
({ ρ } )
⎛ n
⎞
= g ({ ρ } ) + α ⎜ ∑ ( ρ i ) − f ⋅ n ⎟
⎝ i =1
⎠
n
+ ∑ ( β i ρi
i =1
)
( nedre grænse for relativ densitet )
n
+ ∑ ( γ i ( ρi − 1) )
i =1
I.4
( volumenbetingelse )
(I.11)
( øvre grænse for relativ densitet )
OPTIMERINGSBETINGELSE
Den maksimale stivhed opnås ved den mindste tøjningsenergi, det vil sige minimum af målfunktionen. Denne findes ved differentiation af (I.11) med hensyn til den relative densitet af hvert element
∂g M
∂g
=
+ α + βi + γ i = 0,
∂ρi ∂ρi
i = 1, 2,.., n
(I.12)
Det er udnyttet, at de relative densiteter i summationsleddene for elementer forskellige fra det i'te
element er konstante, og derfor udgår ved differentiation. For elementer med relativ densitet forskellig fra 0 og 1, kan (I.10) udnyttes til at forenkle (I.12) til
∂g M
∂g
=
+ α = 0,
∂ρi ∂ρi
i = 1, 2,.., n
(I.13)
Multiplikation af (I.13) med ρi / α giver
−
ρi ∂g
⋅
= ρi , i = 1, 2,.., n
α ∂ρi
(I.14)
For elementer med relativ densitet forskellig fra 0 og 1, opnås den maksimale stivhed dermed når
(I.14) er opfyldt, og den relative densitet overholder (I.5) og (I.6), og α overholder (I.7).
I.5
Differentialet
OPRINDELIG FØLSOMHED
∂g
benævnes følsomheden af et element. Ved differentiation af (I.1) med hensyn til
∂ρi
den relative densitet for et element og udnyttelse af stivhedsestimatet (I.2) fås
∂g
T
= p ρi p −1 {vi } [ ki 0 ]{vi }
∂ρi
(I.15)
Følsomheden af et element kan således bestemmes uafhængigt af de øvrige elementer.
I.6
FILTRERET FØLSOMHED
For at undgå skaktern-mønsteret på grund af den valgte FEM-model, som forklaret i hovedrapport
afsnit 9.3, indføres et filter. Et elements filtrerede følsomhed bestemmes som et vægtet gennemsnit
af den oprindelige følsomhed af elementet og de nærmeste naboelementer ved følgende forskrift
83
Del III Numeriske modeller
I Topologioptimering
⎛
q
∂g
=
∂ρi
∂g
∑ ⎜⎜ ρ ( R − r ) ∂ρ
j =1
⎝
j
j
j
⎞
⎟⎟
⎠
(I.16)
q
∑ ( ρ ( R − r ))
j
j
j =1
hvor
R
er en given filterradius med centrum i elementets eget tyndepunkt [m]
r
er tyngdepunktsafstanden mellem elementet og det j'te element [m]
q
er antal elementer indenfor filterradien [-]
[Pedersen et al. 2006]
Parametrene er illustreret på følgende figur 46.
1
3
r
5
a
R
2
4
6
Figur 46: Ved beregning af et elements filtrerede følsomhed, indgår elementer med tyngdepunkter indenfor
en given filterradius R fra elementets eget tyngdepunkt.
Det ses af (I.16), at gennemsnittet er vægtet på en sådan måde, at naboelementer med lille relativ
densitet og med tyngdepunkter i stor afstand fra elementet får mindst betydning. Ved en filterradius
lig sidelængden af elementet, er den filtrerede og oprindelige følsomhed identisk, og filteret er uden
virkning.
Filtrets virkning er illustreret ved følgende eksempel med 6 elementer, der alle har ρ = 1 . Denne
model er understøttet og belastet svarende til figur 46, hvorfor det kan forventes, at elementerne 1 og
2 bidrager væsentligt til den samlede stivhed, mens element 6 har ringe betydning for den samlede
stivhed.
Med en filterradius R = 1.1⋅ a , svarende til 1.1 gange sidelængde af et element, beregnes følsomheden for element 3 af dets egen følsomhed samt naboelementerne 1, 4, 5. En sammenligning af den
oprindelige og filtrerede følsomhed er vist på figur 47.
84
I Topologioptimering
Del III Numeriske modeller
x 10
-14
2.5
Oprindelig
Filtreret
Følsomhed [Nm]
2
1.5
1
0.5
0
1
2
3
4
5
6
Elementnr i [-]
Figur 47: Oprindelig og filtreret følsomhed ved 6 elementer ved en filterradius R = 1.1 ⋅ a .
Ved en forøget filterradius R = 1.5 ⋅ a , indgår ydermere naboelementerne 2 og 6 ved beregning af
den filtrerede følsomhed af element 3. Idet element 2 har en relativ høj oprindelig følsomhed, fordi
elementet har stor betydning for konstruktionens stivhed, giver det en markant forøgelse af den filtrerede følsomhed for element 3, som vist på figur 48.
x 10
-14
2.5
Oprindelig
Filtreret
Følsomhed [Nm]
2
1.5
1
0.5
0
1
2
3
4
5
6
Elementnr i [-]
Figur 48: Oprindelig og filtreret følsomhed ved 6 elementer. ved en filterradius R = 1.5 ⋅ a
Den filtrerede følsomhed (I.16) indsættes i optimalbetingelsen (I.14).
I.7
OPDATERINGSRUTINE
Der er udformet en opdateringsrutine for elementer med relativ densitet forskellig fra 0 og 1. For at
hindre en uhensigtsmæssig hurtig ændring af den relative densitet, indføres en nedre og øvre grænse
for ændringen
85
Del III Numeriske modeller
I Topologioptimering
⎧(1 − d ) ρi
ai = max ⎨
0
⎩
⎧(1 + d ) ρi
bi = min ⎨
1
⎩
(I.17)
hvor
a
er en nedre grænse for hvert element per iterationskørsel [-]
b
er en øvre grænse for hvert element per iterationskørsel [-]
d
er en given maksimal skridtlængde per iterationskørsel [-]
[Sigmund 2001]
For at sikre konvergens mod den optimale toplogi, justeres højresiden af (I.14) efter værdien af venstresiden af (I.14) på følgende måde
ρi , k +1
η
⎧
⎛ 1 ∂g ⎞
⎪
a,
ρi , k ⋅ ⎜ − ⋅
⎟ <a
⎪
⎝ α ∂ρi ⎠
⎪
η
η
⎛ 1 ∂g ⎞
⎛ 1 ∂g ⎞
⎪
= ⎨ ρi , k ⋅ ⎜ − ⋅
⎟ , a ≤ ρi , k ⋅ ⎜ − ⋅
⎟ ≤b
⎝ α ∂ρi ⎠
⎝ α ∂ρi ⎠
⎪
⎪
η
⎛ 1 ∂g ⎞
⎪
b,
b < ρi , k ⋅ ⎜ − ⋅
⎟
⎪
⎝ α ∂ρi ⎠
⎩
(I.18)
hvor
k
er antal iterationskørsler
η
er en given parameter, der bestemmer opdateringsrutinens effektivitet [-]
[Sigmund 2001]
Der er valgt d = 0.2 og η = 0.5 , som typiske anvendes [Sigmund 2001].
Fordi volumenbetingelsen (I.7) er en monotont aftagende funktion af den relative densitet, kan justeringen af de enkelte elementers relative densiteter ved hver iterationskørsel foretages i en indre bisektionsløkke, også kaldet vinkelhalveringsløkke, således:
1. Gæt en værdi af α
2. Udregn venstresiden af (I.14) for hvert element
3. Opdater de relative densiteter ved (I.18)
4. Gentag punkt 1 hvis volumenbetingelsen (I.7) ikke er opfyldt for denne værdi af α
[Sigmund 2001]
86
J Singularitet af spændinger
J
Del III Numeriske modeller
SINGULARITET AF SPÆNDINGER
Der er i analyserne af konsollen ved hjælp af elementmetode observeret, at spændingen omkring den
fri rand ved indspændingen stadig øges, for stigende finhed af diskretisering. Det kan vises, at når
finheden af diskretisering går mod uendelig, går også randspændingerne mod uendelig; spændingerne har en singularitet. I det følgende beskrives, hvorledes denne kan beregnes.
J.1
PROBLEMATIK
Spændingerne ved konsollens indspænding er et specialtilfælde af det på figur 49 viste tilfælde, hvor
to forskellige materialer er sat sammen, og udsat for mekanisk last.
σ∞
σ∞
E2 ,ν 2
E1 ,ν 1
r
θ
Figur 49: Skitse af problem. σ ∞ er den ydre last, svarende til r = ∞ .
Spændingernes variation i det på figur 49 viste polære koordinatsystem kan beskrives som
σ ij ( r ,θ ) =
K
ω
⎛r⎞
⎜ ⎟
⎝ L⎠
⋅ fij (θ ) + σ ij 0 (θ )
(J.1)
hvor
σ ij
er spændingstensoren
K
r, θ
L
ω
er spændingsintensitetsfaktoren
er de polære koordinater, radius og vinkel
er en karakteristisk længde for den aktuelle geometri
er spændingssingularitetsparameteren
fij
er vinkelleddet
σ ij 0
er det konstante spændingsled
[Munz & Yang 1992, p857]
87
Del III Numeriske modeller
J Singularitet af spændinger
Det er en omfattende opgave at bestemme alle parametrene, og det er ikke beskrevet i det følgende. I
stedet beskrives hvilke afhængigheder der er mellem parametrene, og variationen af spændingerne
ved overgangen mellem de to materialer beskrives.
J.2
AFHÆNGIGHEDER
ω afhænger af de elastiske konstanter for de to materialer. Det gør fij , K og σ ij 0 også, og de afhænger derudover også af den påtrykte last, σ ∞ , vist på figur 49, og geometrien af systemet. Dermed
vil de eneste variable for et system med kendte materialer, geometri og lastsituation, være koordinaterne r og θ . ω er større end nul for langt de fleste tilfælde.
Som det ses af (J.1) vil spændingstensoren σ ij dermed gå mod uendelig for radius r gående mod nul.
Forholdet mellem de elastiske konstanter for de to tilstødende materialer har afgørende betydning,
og en effektiv modul for hvert materiale kan defineres som
E* =
E
ν (1 + ν )
(J.2)
hvor
E*
E
er den effektive modul
er elasticitetsmodulen
ν
er Poissons forhold
[Munz & Yang 1992, p858]
J.3
SPÆNDINGSVARIATION
Spændingernes afhængighed af forholdet mellem de effektive moduler for de to tilstødende materialer er vist på figur 50, hvor normalspændingen σθ normeret med den ydre spænding σ ∞ , vist på
figur 49, er plottet som funktion af radius r normeret med den karakteristiske længe L, for θ = 0 .
Dette er optegnet for forskellige forhold af de effektive parametre, defineret som
E1* − E2*
E1* + E2*
hvor
( )1 , ( )2
88
er indeks der angiver materiale 1 eller 2, som vist på figur 49
(J.3)
J Singularitet af spændinger
Del III Numeriske modeller
Figur 50: Normeret spænding som funktion af normeret afstand for forskellige forhold af effektive moduler[Munz & Yang 1992, p861].
J.4
SINGULARITET FOR KONSOL
Indspændingen af konsollen kan, som nævnt ovenfor, betragtes som et specialtilfælde af dette problem, idet indspændingen kan betragtes som en overgang mellem konsollen, der har de elastiske
koefficienter E og v, og et uendeligt stift materiale. Dermed bliver forholdet defineret ved (J.3) meget stort, og ifølge figur 50 vil der opstå en væsentlig spændingssingularitet, som det er vist i de
foregående afsnit.
89