Foglalkoztatáspolitikai és Munkaügyi Minisztérium Humánerőforrás-fejlesztés Operatív Program
Dr. Páczelt István – Dr. Nándori Frigyes - Dr. Sárközi László Dr. Szabó Tamás - Dr. Baksa Attila - Dluhi Kornél
A végeselemes modellezés kontinuummechanikai alapjai Szakmérnöki jegyzet
Készült
„A felsőoktatás szerkezeti és tartalmi fejlesztése” CAD/CAM/FEM kompetencia kurzusok projekt keretében Miskolci Egyetem
Miskolc 2006.
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
2
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
TARTALOMJEGYZÉK 1.
A VÉGESELEM-MÓDSZER ALAPJAI 1.1. Alapfogalmak 1.2. Rugalmasságtani összefoglaló 1.3. Kompatibilis elmozdulási elemmodell 1.4. Elemek illesztése, a rendszer merevségi mátrixa, redukált terhelési vektora, kinematikai peremfeltétel figyelembevétele 1.5. Egyváltozós feladatok 1.6. Kétváltozós rugalmasságtani feladatok vizsgálata izoparametrikus elemekkel 1.7. Térbeli elemek 2. A MEGOLDÁS ÉS ANNAK HIBÁJA 2.1. Az együttható Mátrix szerkezete 2.2. Hibaanalízis 3. MODELLEZÉSI KÉRDÉSEK 3.1. Alszerkezettechnika 3.2. Adott elmozdulások figyelembevétele 3.3. Adott elmozdulásmezőben fennálló szakadás, kezdeti hézag figyelembevétele 3.4. Excentrikus csatlakozás 3.5. Ferdehatásvonalú támasz figyelembevétele 3.6. Periódikus szerkezet 3.7. Rugalmas ágyazás 4. I-DEAS RENDSZER HASZNÁLATA 4.1. Általános jellemzők: 4.2. Rajzolás az I-DEAS-ban 4.3. Végeselemes analízis 5. C ÁLLVÁNY VIZSGÁLATA 5.1. Geometria létrehozása 5.2. Végeselemes modell 5.3. A feladat megoldása 5.4. Számítási eredmények 6. TENGELYSZIMMETRIKUS FELADAT ELEMZÉSE 6.1. A peremértékfeladat 6.2. A geometria létrehozása 6.3. Végeselemes hálózás 6.4. peremfeltételek előírása 6.5. A feladat megoldása 6.6. Az eredmények kiértékelése 7. TENGELYSZIMMETRIKUS FELADAT ELEMZÉSE 7.1. A peremértékfeladat ismertetése 7.2. A geometria felépítése 7.3. Végeselemes hálózás 7.4. Peremfeltételek előírása 7.5. A feladat megoldása 7.6. Az eredmények szemléltetése 8. TÉRBELI FELADAT ELEMZÉSE 8.1. A geometria előállítása 8.2. Peremfeltételek megadása 8.3. Végeselemes háló generálása 8.4. Feladat megoldása, eredmények kiértékelése 8.5. A háló finomítása 9. TÉRBELI FELADAT ELEMZÉSE 9.1. A geometria módosítása 9.2. Peremfeltételek 9.3. Végeselemes háló 9.4. Eredmények 10. IRODALOMJEGYZÉK
3
7 7 9 18 22 25 37 58 65 65 67 71 71 73 73 75 76 79 79 81 81 83 85 87 88 88 90 91 93 93 93 94 95 96 96 99 99 100 101 101 102 103 105 105 106 107 108 109 111 111 112 113 113 117
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
4
ELŐSZÓ A jegyzet felsőfokú végzettségű, elsősorban mérnökök számára készült, akik tanulmányaik során megismerkedtek a szilárdságtan és a lineáris rugalmasságtan alapfogalmaival. Jártasságot szereztek a vektor és tenzorszámítás, illetve a mátrix algebra területén. A jegyzet célul tűzi ki bonyolult geometriával és terheléssel rendelkező gépészeti szerkezetek szilárdságtani analízisét. A megoldás a végeselem-módszer alkalmazásával történik. Ennek megfelelően a jegyzet rövid összefoglalókra építve áttekintést ad a módszer általános kérdéseiről, majd sorra veszi a mérnöki gyakorlatban leginkább elterjedt végeselemes eljárásokat: egyváltozós (rúd) feladatokat, kétváltozós (síkbeli, forgásszimmetrikus, illetve hajlított lemez) feladatokat és térbeli feladatokat. A módszer tárgyalása során külön hangsúlyt kapnak az izoparametrikus elemcsalád alkalmazása, illetve a modellezés speciális problémái (alszerkezettechnika, ferde görgő, excentrikus csatlakozás, elmozdulás-mezőbeli szakadás stb.). Végül bemutatjuk az I-DEAS programrendszer végeselemes alkalmazását néhány alapfeladat megoldásán keresztül. Érdemes itt is kiemelni, hogy a vektor és tenzor jellegű mennyiségek kétféle jelöléssel is előfordulnak. Ennek megfelelően jelöléskor a vektor vastagon szedett dőlt kisbetű (a, b), míg a tenzor vastagon szedett dőlt nagybetű (A, T). Ezen mennyiségeknek megfelelő mátrixos jelölésben vastagon szedett álló kis, illetve nagy betű szerepel (a, b, A, T). A szerzők az alábbiakban osztották szét a munkát a fejezetek megírásával kapcsolatosan. 1.1-1.4. fejezet Nándori Frigyes, 1.5, 3. fejezet: Páczelt István, 1.6 fejezet: Sárközi László, Nándori Frigyes, 1.7, 2, 6, 7. fejezet: Szabó Tamás, 4, 5. fejezet Baksa Attila, 8, 9. fejezet Dluhi Kornél. Miskolc, 2006. november Sikeres tanulmányokat kivánnak
A szerzők
1. A VÉGESELEM-MÓDSZER ALAPJAI 1.1. Alapfogalmak A végeselem-módszer (továbbiakban VEM) több mint 30 éve szerepel a gépészmérnökképzés tananyagában [1]. A módszer kialakulását a Courant által a csavarási feladat közelítő megoldásánál használatos szakaszonkénti (háromszögletű tartományok feletti) csavarási feszültségfüggvény approximációja jelentette 1943-ban. 1956-ban Turner és társai síkrugalmasságtani feladatot oldottak meg az elmozdulásmező négyszögletű altartományok feletti közelítésével, a hagyományos Ritz-féle módszer lokális közelítő függvényeken keresztüli alkalmazásával. Clough 1960-ban ennek az eljárásnak a végeselem-módszer nevet adta. Az elmúlt ötven évben a módszer látványos fejlődésének vagyunk szemtanúi. A 60-as évekre a rugalmasságtani feladatainak megoldását szolgáló elemcsaládok kifejlesztése, sokoldalú modellezési lehetőséget nyújtó végeselem-programok (ASKA, NASTRAN, SAP) megjelenése a jellemző. A 70-es években elkezdődik a számítási hibák analízisének kutatása Babuska cikkének megjelenésével. Elindul 1973-ben Szabó Barna javaslatára, az un. pverziójú számítás a hozzátartozó elemek kidolgozásával. Sorra kerülnek a nemlineáris feladatok vizsgálatára alkalmas módszerek kidolgozása, számítógépi programok alkalmazásba vétele (NONSAP, ABAQUS, ADINA, ANSYS, COSMOS/M, FEAP, MARC, SYSTUS stb.). Megjelennek a p-verziójú elemeket hordozó programok, PROBE, StressCheck, RASNA. A CAD rendszerekkel összekapcsolt rendszerek alakulnak ki az 1980-as években, amelyeknek a fejlődése mind a mai napig tart (CATIA, I-DEAS, MSC/NASTRAN, Patran, Pro/Engineer, SolidWorks, stb.). A kapcsolt feladatokra (szilárdságtani, hőtani, áramlástani, villamosságtani stb.) programok nyernek kidolgozást a 1990-es évek óta (FLUENT, ProCast, SysWeld, DEFORM stb.) A módszernek a gépészeti tervezés folyamatában való elhelyezéséről tájékoztat a következő 1.1 ábra. Itt a legutolsó rublikába eső különböző variációs elvek és a számítástechnika fejlődése alapozta meg a módszer elterjedését. Equation Chapter 1 Section 1
1.1. ábra. Végeselem-módszer elhelyezése A továbbiak röviden összefoglalják a későbbiekben felhasználásra kerülő matematikai fogalmakat és összefüggéseket [4].
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
Függvény: u ( x ) ahol u a függő, x a független változó. Funkcionál: F ( x, u , u′) egy olyan függvény, melynek a független változói között függvények esetleg deriváltjaik is szerepelnek. Pl. F ( x, u , u`) = c1 ∫u dx + c2 ∫u `dx Variáció: δ u az u függvény variációja, lásd 1.2. ábra, értelmezés szerint δu = αv( x) ahol α: paraméter, amely a különböző variációknál más és más, v(x): egy másik függvény.
1.2. ábra. Variáció értelmezése d δ u = α dv = δ u ` , azaz a dx ( ) dx
deriválás és a variáció képzés sorrendje felcserélhető.
A funkcionál első variációja a fenti F funkcionál esetén definició szerint: ∂F ∂F δu + δ u′ δ F = ′ ∂u ∂u ∂F ∂F ∂F dx + du + du ′ dF = ∂x ∂u ∂u ′
F (u + δ u ) = F (u ) + δ F + δ 2 F + … A variációképzés azonosságai közül érdemes az alábbiakra felhívni a figyelmet:
δ ( F1 ⋅ F2 ) = δ F1 ⋅ F2 + F1 ⋅ δ F2 δ ( FF ) = δ F ⋅F( F−δ) F ⋅F 1
1
2
2
1
2
2
2
δ ( F n ) = n ⋅ δ ( F n −1 )
δ ∫ u ( x)dx = ∫ δ u ( x)dx ( L)
( L)
Integrál átalakítási tételek közül az alábbiakra lesz szükség: →
∫B ⊗ ∇dV = ∫B ⊗ n dA
V
L
Gauss-Osztograszkij tétel
A
L
L
L
dw d dv L ∫0v dx dx = ∫0 dx ( v ⋅ w ) dx − ∫0 dx w dx = [v ⋅ w]0 − ∫0v′w dx
A vektor és tenzor jellegű mechanikai jellemzők jelölésére egyaránt használjuk a szimbólikus és a mátrixos (általában (x,y,z) koordinátarendszerben kapott koordinátákból alkotott) jelölést. A vektor és tenzor jele vastagon szedett dőlt kis– illetve nagybetű. Mátrixos jelölésnél pedig álló helyzetű kis- illetve nagybetű. Legyen A és B az (x,y,z) koordinátarendszerben adott 3 sorú és 3 oszlopú, azaz (3,3)-as tenzor:
8
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
a xx Α = a yx azx
axy a yy a zy
axz a yz , azz
bxx B = byx bzx
bxy byy bzy
bxz byz bzz
A két tenzor kétszeres skaláris szorzása az alábbi módon értelmezett: A ⋅⋅B = axx ⋅ bxx + axy ⋅ bxy + … + azz ⋅ bzz
Ha A = AT akkor a tenzor szimmetrikus, ha pedig B = − BT akkor a tenzor aszimmetrikus, ha A = AT és B = c d , azaz a c és d vektorok általános szorzata cx d x B = c y d x cz d x
cx d y . .
cx d z . .
akkor A ⋅ ⋅B = A ⋅⋅ ( cod ) = c ⋅ A ⋅ d .
1.2. Rugalmasságtani összefoglaló Tekintsünk egy V térfogatú lineárisan rugalmas térbeli testet, melyet A határoló felület vesz körül. Az A felület felbontható egy Au - melyen az elmozdulás elő van írva - és egy Ap felületrészre, ahol pedig a terhelés előírt.
(A = A
u
∪ Ap ) . A test egy tetszőleges pontjába
mutató helyvektor r = xe x + ye y + ze z . Lásd 1.3 ábrát.
1.3. ábra. Rugalmas test
A lineáris rugalmasságtan ismeretlen változói rendre az elmozdulás vektor , az alakváltozási-, illetve feszültségi tenzor. Skaláris koordinátái az alábbi módon értelmezhetők: elmozdulás: (3 db ismeretlen)
u ( r ) = ue x + ve y + we z
(1.1)
alakváltozási tenzor: (6 db ismeretlen)
9
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
εx A ( r ) = 12 γ yx 12 γ zx
γ xy εy 1 2 γ zy
1 2
1 2 1 2
γ xz ∂u ∂u ∂v γ yz , ε x = ,… , γ xy = + ,… ∂x ∂y ∂x ε z
(1.2)
melyben az ε i , γ ij : fajlagos nyúlások és szögtorzulások dimenziótlan mennyiségek. A feszültségi tenzor: (6 db ismeretlen) σ x τ xy τ xz (1.3) T ( r ) = τ yx σ y τ yz τ zx τ zy σ z a σ és τ normál és nyíró feszültségeket jelöl. A rugalmasságtan linearizált elmélete kis elmozdulás (az elmozdulás lényegesen kissebb mint a test mérete), illetve kis alakváltozás (az alakváltozás jellemzői lényegesen kissebbek mint 1) esetén alkalmazható. Ekkor az ismeretlenek meghatározására szolgáló egyenletek az alábbiak: A=
1 ( u ∇ + ∇ u) 2
T = D ⋅⋅ A r ∈ V
r ∈V
- kinematikai egyenlet
(1.4)
- anyagegyenlet
illetve a ν T = 2G A + 1 − 2ν AI ⋅ I
r ∈V
- Hook-törvény
T ⋅ ∇ + ρ ⋅ k = 0 - egyensúlyi egyenlet
(1.5)
Peremfeltételek:
u = u0
r ∈ Au
T ⋅ n = p r ∈ Ap
(1.6)
- kinematikai peremfeltétel
(1.7)
- dinamikai peremfeltétel
Itt G – az anyag csúsztató rugalmassági modulusa, ν - a Poisson tényező, ρ a tömegsűrűség, k pedig a tömegen megoszló terhelés intenzitás vektora. AI az alakváltozási tenzor első skalár invariánsa, D az anyagjellemzők később részletezett mátrixa. A rugalmasságtan feladata a fenti differenciálegyenlet-rendszer megoldása. Közelítő számítás esetén: elmozdulásra alapozva a számítást a fenti egyenleteket kell egymást követve kielégíteni. Ha pedig közelítő feszültségekből indulunk ki, akkor a kinematikai egyenletek helyett a vele egyenértékű kompatibilitási egyenletet kell biztosítani, amely a következő alakban írható: ∇ × A×∇ = 0
(1.8)
1.2.1. Energia elvek, variációs módszerek A továbbiak először a fenti differenciálegyenlet-rendszer néhány egyenletét kielégítő, speciális mezőket értelmeznek. Kinematikailag lehetséges elmozdulásmező alatt olyan elmozdulást értünk, amely kielégíti a *
kinematikai peremfeltételt, folytonos és deriválható. Jele: u . ∗
Származtatható belőle egy kinematikailag lehetséges alakváltozás: A
10
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
∗ ∗ 1 ∗ A = u ∇ + ∇ u, 2
∗
∇ × A× ∇ = 0 r ∈ Au
(1.9)
∗
azaz A kinematikailag lehetséges alakváltozás a geometriai egyenletet kielégíti és egyúttal kompatibilis. Statikailag lehetséges feszültségmező alatt olyan T jelű feszültséget értünk, amely kielégíti az egyensúlyi egyenletet és a dinamikai peremfeltételt, azaz:
T ⋅∇ + ρ ⋅ k = 0 r ∈ V
T ⋅ n = p r ∈ Ap ,
(1.10)
Származtatható belőle egy statikailag lehetséges alakváltozás: A = D −1 ⋅⋅T de ez általában nem kompatibilis, és a kinematikai peremfeltétel sem biztosított. Azaz: ∇ × A × ∇ ≠ 0 r ∈V
(1.11)
u ≠ u0
(1.12)
1.2.1.1.
r ∈ Au .
Virtuális munka elve:
A fentiekben bevezetett statikailag illetve kinematikailag lehetséges feszültség- és elmozdulásmezőből kiindulva, az alábbi átalakítások tehetők: Az egyensúlyi egyenletet a kinematikailag lehetséges elmozdulásmezővel megszorozva, majd a kapott skaláris mennyiséget térfogaton integrálva ∗
T ⋅∇ + ρ ⋅ k = 0
⋅ u,
∫
V
∗
(
)
∗
∫ u⋅ T ⋅∇ dV + ∫ u⋅ ρ ⋅ k dV = 0
V
V
a következő szorzat deriválási szabály és az elmozdulásmező derivált tenzora felbontási lehetőségének figyelembevételével
(
)
∗ ∗ ∗ u⋅ T ⋅ ∇ = u⋅ T ⋅∇ − u ∇ ⋅⋅T ,
∗ ∗ ∗ 1 ∗ 1 ∗ u ∇ = u ∇ + ∇ u + u ∇ − ∇ u , 2 2 ∗
(
∗
∗
Ψ =−Ψ T
A
)
∗ ∗ ∗ u⋅ T ⋅∇ = u⋅ T ⋅ ∇ − A⋅⋅T
∗
, Ψ ⋅⋅T = 0
a Gauss-Osztograszkij tétel alkalmazásával, nyerjük a virtuális munka elve legáltalánosabb alakját: ∗
∫ A⋅⋅T dV =
V
∫
A = Au + Ap
∗
∗
u⋅ T ⋅ n dA + ∫ u⋅ ρ ⋅ k dV
(1.13)
V
A virtális munka elv variációs alakjának felírásához először értelmezzük az elmozdulás mező variációját (virtuális elmozdulást), melynek jele: δ u
11
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI ∗
δ u = u− u,
(1.14)
u : egzakt elmozdulás
Nyilvánvalóan teljesül: δ u = 0 , ha r ∈ Au . Hasonlóan értelmezhető az alakváltozás variációja: ∗ ∗ 1 A = A u = A ( u + δ u ) = u + δ u ∇ + ∇ 2 1 1 = ( u ∇ + ∇ u )+ (δ u ∇ + ∇ δ u ) 2 2
( u + δ u ) =
δA
A
Azaz: ∗
A = A+δ A A virtuális munka elv egzakt u elmozdulással a következő alakú:
(1.15)
∫u ⋅ T ⋅ n dA + ∫u ⋅ ρ ⋅ k dV = ∫ A ⋅⋅T dV A
V
V
Ezt kivonva az (1.13) egyenletből kiadódik a virtuális munka elv variációs alakja:
∫δ A ⋅⋅T dV = ∫ δ u ⋅ ρ ⋅ kdV + ∫δ u ⋅ p dA .
V
V
(1.16)
Ap
A baloldal, az ú.n. belsőerők alakváltozási munkájának variációja megegyezik a jobboldalon felírt ú.n. külső erők virtuális munkájával. A T mező helyett szabad (nem föltétlen statikailag lehetséges) mezőt feltételezve az elv az alábbi alakban is megfogalmazható: Ha az (1.16) tetszőleges kinematikailag lehetséges elmozdulásmező variációja mellett fennáll, akkor T feszültségi tenzormezőmező statikailag lehetséges. 1.2.1.2.
Potenciális energia minimuma elv származtatása:
∗
Egy u kinematikailag lehetséges elmozdulásmezőből kiindulva a megfelelő egyenletek ∗
∗
segítségével származtatható az A alakváltozási és a T feszültségi tenzormező. A potenciális energia értelmezése a következő:
∗ ∗ Π p u = Π p =
1 ∗ ∗ A⋅⋅T dV 2 V∫ ∗
U : alakváltozási energia
∗
∫ u⋅ ρ ⋅ k dV −
−
V
∗
∫ u⋅ p dA
(1.17)
Ap
∗
W k : külső erők (virtuális ) munkája
A potenciális energia minimuma elv szerint a kinematikailag lehetséges elmozdulás mezőkhöz tartozó potenciális energiák az egzakt megoldásnál minimummal rendelkeznek, azaz: ∗ Π p u ≥ Π p ( u)
(1.18) ∗
egyenlőség akkor áll fenn, ha u ≡ u . Az elv az alábbi lépések és átalakítások alapján bizonyítható.
12
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
∗ ∗ ∗ ∗ ∗ Π p = Π p u = Π p ( u + δ u ) , ahol A = A + δ A , T = D ⋅⋅ A
1 ( A + δ A) ⋅⋅D ⋅⋅ ( A + δ A ) dV − ∫ ( u + δ u ) ⋅ ρ ⋅ k dV − ∫ ( u + δ u ) ⋅ p dA = 2 V∫ V Ap
Π p ( u + δ u) =
T
1 A ⋅⋅ D ⋅⋅ A dV − ∫u ⋅ ρ ⋅ k dV − ∫ u ⋅ p dA+ 2 V∫ V Ap
=
Π ( u ) : egzakt értékhez tartozó T
+ ∫δ A ⋅⋅ D ⋅⋅ A dV − ∫δ u ⋅ ρ ⋅ k dV − ∫ δ u ⋅ p dA+ V
V
Ap
δΠ : a potenciális energia első variációja δT
1 + ∫δ A ⋅⋅ D ⋅⋅δ A dV 2V δ 2 Π ≥ 0, =0
Π p ( u + δ u ) = Π p ( u ) + δ Π p + δ 2Π p
(1.19)
Az (1.16) virtuális munka elv szerint: δ Π p = 0 , vagyis valóban fennáll az (1.18) alatti egyenlőtlenség. A variáció zérus előírását a variációszámítás stacionaritási feltételének nevezzük. 1.2.1.3.
Lagrange-féle variációs elv
A Lagrange-féle variációs elv szerint a teljes potenciális energia variációja a megoldásnál zérus. Ennek igazolására induljunk ki a variáció értelmezéséből!
δ Π p = ∫δ A ⋅⋅T dV − ∫δ u ⋅ ρ ⋅ k dV − ∫ δ u ⋅ p dA = 0 V
∇ δu =
V
(1.20)
Ap
1 1 ( ∇ δ u + δ u ∇ )+ ( ∇ δ u − δ u ∇ ) Miután δΨ = −δ Ψ T 2 2 δA
δΨ
( ∇ δ u ) ⋅⋅T = ( ∇ ⋅ T ⋅ δ u ) − δ u ⋅ (T ⋅∇ ) = δ A ⋅⋅T
∫
A = Ap + Au
n ⋅ T ⋅ δ u dA − ∫δ u ⋅ (T ⋅∇ ) dV − ∫δ uρ ⋅ k dV − ∫ δ u ⋅ p dA = 0 V
V
Ap
továbbá az Au felületen δ u = 0 így az előző egyenlet formailag a következő alakra hozható: − ∫δ u ⋅ V
(T ⋅∇ + ρ ⋅ k ) = 0 : egyensúlyi egyenlet
dV + ∫ δ u ⋅ Ap
( n ⋅T − p)
dA = 0
(1.21)
= 0 : dinamikai peremfeltétel
Az elv tehát tetszőleges variáció esetén biztosítja az egyensúlyi egyenlet és a dinamikai peremfeltétel teljesülését. Általában u -t közelítjük úgy, hogy adott függvények és ismeretlen paraméterek kombinációját vesszük, majd így írjuk elő a variáció eltűnését. Ez a feltétel mindig ad egy lehetséges legjobb megoldást. A módszer célja olyan mező kiválasztása, amely az egyensúlyi egyenletet, illetve a dinamikai peremfeltételt a lehető legjobban biztosítja [1]. 13
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.2.1.4.
Ritz-féle módszer
A potenciáis energia minimuma elvet alkalmazva, legyen a kinematikailag lehetséges elmozdulás a következő: ∗
∗
N
N
N
i =1
j =1
k =1
u = u0 + ∑ci ⋅ ϕi ( r ) e x + ∑cN + jψ j ( r ) e y + ∑c2 N + k χ k ( r ) ez
(1.22)
∗
ahol u0 kielégíti a kinematikai peremfeltételt, a választott folytonos koordinátafüggvényekre az egyszeres véges-mértékű deriválhatóságon kívül áll az alábbi megkötés
ϕi = ψ i = χ i = 0 , r ∈ Au , i = 1,… , N
(1.23)
Ezzel potenciális energia kifejezés egy többváltozós valós függvény. Π p = Π p ( c1 , c2 , … , c3 N )
(1.24)
A variáció eltűnésének feltétele az alábbi alakban írható
δΠ p =
∂Π p ∂c1
δ c1 +
∂Π p ∂c2
δ c2 + … +
∂Π p ∂c3 N
δ c3 N = 0
(1.25)
Itt δ ci tetszőleges kell, hogy legyen! Ezáltal kapunk egy lineáris egyenletrendszert a paraméterekre nézve. Bevezetve az ismeretlenek oszlopvektorát cT = [ c1 c2
c3 N ]
(1.26)
a stacionaritási feltétel tömör formában is írható:
δΠ = δ cT ∂Π p
= ∂c
∂Π ∂c
=0
∂Π p
∂Π p
∂c1
∂c2
…
∂Π p ∂c3 N
T ∂Π p =0⇔ = 0, i = 1,… 3 N ∂ci
(1.27)
Példa 1.1 Ritz-módszer bemutatására egyváltozós feladatot vizsgálunk.
P1.1-1. ábra. Húzott rúd
A P1.1-1. ábrán látható húzott-nyomott rúd terhelése p hosszmentén megoszló súlyterhelés, és FL koncentrált
erő. A rúd pontjainak rúdirányú elmozdulása u ( x ) , csak a keresztmetszet x koordinátájától függ. A rugalmasságtan általános egyenletei ebben az egyváltozós esetben az alábbi alakra egyszerűsödnek:
εx =
∂u du = = u ′ - kinematikai egyenlet ∂x dx
(a)
σ x = E ⋅ ε x = E ⋅ u ′ - Hooke-törvény
(b)
14
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
dN + p = 0 - egyensúlyi egyenlet dx
(c)
Kinematikai peremfeltétel: u ( x = 0 ) = 0 . A dinamikai peremfeltétel a rúdvég koncentrált terhelésének és rugós
megtámasztásának megfelelően: N = A ⋅ σ x azaz,
|
N ( x = L ) = A ⋅ E ⋅ u ′ L = FL + Fc = FL − c ⋅ uL
(d)
A potenciális energia ebben az esetben: A rúd belső alakváltozási energiája A külső terhelések munkája
1 ε x ⋅ σ x dV 2 V∫ A dx
Π p (u ) = L
rugóenergia
L
−
∫u p dx − u
L
⋅ FL
+
0
1 c ⋅ uL2 2
L
1 1 2 A E ( u ′ ) dx − ∫u p dx − u L ⋅ FL + c ⋅ uL2 2 ∫0 2 0
Πp =
(e)
Ennek képezve a variációját kapjuk, hogy L
L
δΠ p = ∫ A E u ′ δ u ′ dx − ∫δ u ⋅ p dx + c ⋅ u L ⋅ δ uL − δ uL ⋅ FL = 0 0
0
(δ ( u ′ )
2
= 2u ′δ u ′
)
NL ′ δΠ p = − ∫ ( A ⋅ E ⋅ u ′ ) + p δ u dx + ( A ⋅ E ⋅ u ′ )L − ( FL − c ⋅ uL ) δ uL = 0 0 = 0 : egyensúlyi egy. L
= 0:
(f)
dinamikai peremfeltétel
Közelítő megoldás kereséséhez induljunk ki az alábbi másodfokú polinomból: u ( x ) = c0 + c1 x + c2 x 2
(g)
Kinematikailag lehetséges, ha folytonos és deriválható, továbbá ha kielégíti a kinematikai peremfeltételt, azaz u(0)=0, amiből következik, hogy c0 = 0 azaz ∗
u ( x ) = c1 x + c2 x 2
(h)
Behelyettesítés után a potenciális energia 2
2 L L∗ ∗ ∗′ 1 ∗ ∗ ∗ 1 Π p u = ∫A⋅ E ⋅ u dx − ∫ u⋅ pdx + c ⋅ uL − uL ⋅ FL = 2 20 0 L
=
L
2 1 1 2 A⋅ E ( c1 + 2c2 x) dx − ∫ ( c1x + c2 x2 ) ⋅ pdx + c ⋅ ( c1L + c2 L2 ) − ( c1L + c2 L2 ) ⋅ FL = 2 ∫0 2 0
L2 1 L2 L3 L3 = A⋅ E c12 ⋅ L + 4c1c2 + 4c22 − p ⋅ c1 + c2 + 2 2 3 3 2 1 + c c12 L2 + 2c1c2 L3 + c22 L4 − ( c1L + c2 L2 ) ⋅ FL 2 A variációképzéssel:
δ Πp =
∂Π p ∂c1
δ c1 +
∂Π p ∂c2
δ c2 = 0
⇒
∂Π p ∂c1
= 0,
∂Π p ∂c2
=0
azaz kifejtve:
∂Π p ∂c1 ∂Π p ∂c2
= AELc1 + AEL2 c2 − p = AEL2 c1 +
L2 + cL2 c1 + cL3 c2 − LFL = 0 2
4 L3 AEL3 c2 − p + cc1 L3 + cL4 c2 − L2 FL = 0 3 3
ugyanezt mátrixos formában felírva:
15
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
L AE 2 L
L2 L2 c1 c + 3 4 3 L c2 L 3
L3 c1 p L2 + LFL = 3 L4 c2 p L3 + L2 FL 2
(i)
Speciális esetként először vizsgáljuk meg azt az esetet, amikor csak koncentrált erővel terhelt a rúd:
P1.1-2. ábra. Koncetrált terhelés
Ekkor az (i) egyenlet megoldása:
FL . A⋅ E Ebben az esetben a lehetséges elmozdulás a következő alakban adódik c2 = 0 , és c1 =
∗
FL x AE melyből látható, hogy a Π minimum elv érvényben van hiszen teljesül: u=
∗
(j)
∗
N = AE u ′ = FL dinamikai peremfeltétel dN = 0 → N = áll. N L = FL egyensúlyi egyenlet dx Másik speciális esetként vizsgáljuk meg azt az esetet, amikor csak megoszló terheléssel terhelt a rúd
P1.1-3. ábra. Megoszló terhelés
Ekkor a (i) egyenlet az alábbi alakban írható
pL 2 AE 4 pL c1 + Lc2 = 3 3 AE melyből: c1 + Lc2 =
p⋅L p , c2 = − A⋅ E 2A⋅ E Ebben az esetben a lehetséges elmozdulás a következő alakban adódik c1 =
∗
pL p ⋅x− ⋅ x2 AE 2A⋅ E A kapott megoldás egzakt, mivel az egyensúlyi egyenletet, illetve a dinamikai peremfeltételt kielégíti. u=
∗
(k)
∗
N = AE u ′ = pL − px = p ( L − x ) Hasonló módon belátható, hogy a legáltalánosabb esetben, amikor rugalmas megtámasztású a rúd vége szintén egzakt megoldást kapunk. Ha viszont a megoszló terhelés jellege megváltozik –nem konstans, vagy szakaszonként különböző állandó- akkor a fenti másodfokú közelítés már nem vezet egzakt megoldásra.
16
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
Több testből álló rendszer
1.2.1.5.
1.4. ábra. Kételemű rendszer
Tekintsük az 1.4. ábrán vázolt 1 és 2 jelű testekből álló szerkezetet. Az egyes elemekhez tartozó mennyiségeket felső indexbe tett sorszám jelöli. Mindkét elemre egyidejűleg az e felső indexszel hivatkozunk. A V e térfogatú e elemet az Ae felület határolja. A V e térfogaton a ρ e k e sűrűségű megoszló terhelés, az Aep felületen a p sűrűségű felületi terhelés működik, míg az Aue felületen ismert az u0 elmozdulás. Az Ae felület megmaradó Ace részén az elem a szomszédos elemmel érintkezik. Mindkét testre érvényesek a rugalmasságtan egyenletei, azaz Ae =
1 e ( u ∇ + ∇ ue ) 2
T e = D e ⋅⋅ Ae
r ∈V e
r ∈V
- kinematikai egyenlet
- anyagegyenlet
(1.28) (1.29)
T e∇ + ρ k e = 0 r ∈ V e − egyensúlyi egyenlet
(1.30)
mint mezőegyenlet, továbbá érvényesek a peremfeltételek [1]: ue = u0
r ∈ Aue
T ⋅u = p e
e
− kinematikai peremfeltétel
r ∈ Ape
e
− dinamikai peremfeltétel
(1.31)
illetve az Ace felületen az illesztési feltételek:
u1 = u 2
kinematikai
(1.32)
T ⋅ n = − T ⋅ n2 dinamikai 1
1
2
A dinamikai illesztési feltétel a közös érintkezési felület pontjainak kölcsönhatását fejezi ki. Ennek értelmében az egymásra átadódó feszültségek egymással ellentétesek. Itt n1 és n 2 a testekből kifelé mutató normálvektorok. Vizsgáljuk meg hogy ebben az esetben hogyan használható a potenciális energia minimuma elv, illetve a megfelelő variációs elv. ∗
Tételezzük fel, hogy az ue elmozdulásmező kinematikailag lehetséges, illetve teljesíti a kinematikai illesztési feltételt 17
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
Feltesszük továbbá, hogy
δ u1 = δ u2
r ∈ Ac12
(1.33)
A két testre vonatkozó variációs elv felépítéséhez induljunk ki a teljes potenciális energiák variációjából [3]:
∑ {δ U 2
e
e =1
− δ Wke } = 0
2
∑{∫δ A ..T dV − ∫δ u e
e =1
V
e
e
V
e
e
}
⋅ ρ ⋅ k dV − ∫δ ue ⋅ p dA = 0 e
(1.34)
A
δU e
δ Wke
Az első integrál az alábbi átalakítások alapján átírható:
∫δ u ⋅ (T ⋅∇ )
e
dV =
Ve
∫ (δ u ⋅ T ⋅∇ )
Ve
e
dV − ∫ (δ u ∇ )e ..T e dV = Ve
δ A +δ Ψ
= ∫δ u ⋅ T ⋅ n dA − ∫δ A ..T e dV e
e
e
e
Ae
Ve
innen
∫δ A
e
Ve
⋅⋅T e dV = ∫δ ue ⋅ T e ⋅ ne dA − ∫δ u ⋅ (T ⋅ ∇ ) dV e
Ae
Ve
ahol Ae = Aue + Ape + Ace δ u=0
A kapott formulákat visszaírva és a tagokat átrendezve:
∑ ∫δ u ⋅ [T ⋅ ∇ + ρ k ] 2
e =1
V e
e
e
e dV − ∫ δ ue ⋅ (T ⋅ n − p ) dA − Aep
− ∫ δ u ⋅ (T ⋅ n + T ⋅ n 1
1
2
2
) dA = 0
(1.35)
Ac12
Ennek alapján a variációs elv biztosítja az egyensúlyi egyenlet, a dinamikai perem-, és illesztési feltétel teljesülését
1.3. Kompatibilis elmozdulási elemmodell 1.3.1. Az elmozdulásmező közelítésének felépítése, csomóponti elmozdulás vektor Végeselem-módszer alkalmazásakor első lépésben a tartományt véges kiterjedésű részekre u.n. elemekre bontjuk. Az elemeket sorszámozzuk. Elemenként külön – külön közelítjük az elmozdulást úgy, hogy az a teljes testre kinematikailag lehetséges elmozdulássá legyen illeszthető. Ez azt jelenti, hogy kinematikailag lehetséges elmozdulás közelítésből indulunk ki, azaz teljesül az elemek határán az illesztési vagy kinematikai peremfeltétel, továbbá a deriváltak szakaszonként (elemenként) folytonosak. Ebben a szakaszban feltételezzük, hogy az elmozdulásmező elemenként a hely (x,y,z)=x vonatkoztatási koordinátarendszer függvényeként áll elő. Később látni fogjuk, hogy gyakran érdemes az elemhez kötött helyi koordináta – rendszert alkalmazni. Illesztés céljából az elemek határán jelöljünk ki pontokat
18
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
u.n. csomópontokat. Ezeket is sorszámozzuk. Az egymáshoz kapcsolódó elemek egybeeső csomópontjainak elmozdulásait fogjuk azonossá tenni. A csomópontok számát és a közelítés típusát úgy kell megválasztani, hogy a csomóponti elmozdulások azonossága biztosítsa a teljes érintkezési tartományon a folytonosságot. Az így felépített elemet kompatibilis elmozdulási elemnek fogjuk nevezni [1]. Az előzőkkel összhangban egy e – jelű elem x rendszerbeli u elmozdulás-mezőjét u e = u e (x ) = Φ e (x )c e = Φ ec e
(1.36)
alakban közelítjük, ahol ϕT Φ e = 0T 0T
0T ϕT 0T
0T 0T ϕT
(1.37)
az alapfüggvénymátrix, ce az ismeretlen paraméterek vektora, ϕT (x) pedig rendszerint hatványfüggvényeket tartalmazó sorvektor: ϕT (x) = 1 x
y 2 … (1.38) A továbbiakban az a cél hogy az állandók helyett áttérjünk a csomópontok elmozdulásaira. Legyen az e jelű i, j, k, csomópontokkal rendelkező elem csomóponti elmozdulásainak vektora a következő: y
z
x2
qie e ui e q q e = ej (1.39) , qie = vi q k w i Feltéve, hogy a csomóponti elmozdulások koordinátáinak száma megegyezik a közelítésben felvett paraméterek számával a ce vektor kifejezhető q e vektorral az alábbiak szerint: -1
Φ e ( x i ) e e -1 e e e c = (1.40) q = G q ≡ V q itt a felsőindexben szereplő -1 az invertálásra utal. Természetesen lehetőség van arra is, hogy az ismeretlen paraméterek száma nagyobb legyen mint a csomóponti elmozdulás koordináták száma. Ekkor más technikára (u.n. pótlólagos állandók kijelölése, vagy belső csomópontok értelmezése) van szükség. Ebben a fejezetben ettől eltekintünk. (1.40) alapján a (1.36) közelítés az alábbi alakban írható: e
u e = u e (x) = Φ e (x)V e qe ≡ N e (x)q e = N eq e
(1.41)
ahol az N e mátrixot az elem approximációs mátrixának nevezzük. Nyilvánvaló , hogy N e a csomópontok szerint particionálható: N e = N i
Nj
N k …
e
(1.42)
Vizsgáljuk meg, hogy az N e mátrixnak milyen feltételeket kell kielégíteni. Nyilvánvaló, hogy az elmozdulásnak a csomópontban meg kell egyeznie a csomóponti elmozdulásvektorral, azaz N ie (xi ) = E , N ie (x j ) = 0
(1.43)
19
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
itt E az egységmátrix, 0 pedig a nullmátrix. Természetes elvárás, hogy az elmozdulás közelítése tartalmazza az un. merevtestszerű mozgást, azaz az eltolásra és a forgásra zérus alakváltozási energiát adjon. További elvárás, hogy az elem térfogati és torzulási energiáját külön-külön lehessen számolni. Mivel a közelítő elmozdulásnak kinematikailag lehetségesnek kell lenni ez azt jelenti, hogy az adott előírásokat a csomóponti paraméterek megválasztásával kell biztosítani (a deriválhatóság a fenti polinomokra mindig teljesül). Az elmozdulásmező közelítéséből kiindulva származtathatók az elem további szilárdsági jellemzői.
1.3.2. Alakváltozás-, feszültségi vektor, terhelési vektorok Az alakváltozási tenzor szimmetriáját felhasználva a tenzor elemeiből 6 méretű vektor állítható elő: e ∂u ∂ ε x ∂x ∂x ε ∂∂yv 0 y ε z ∂∂wz 0 e A ⇒ ε = = ∂u ∂v = ∂ + γ xy ∂y ∂x ∂y γ yz ∂v + ∂w 0 ∂z ∂y γ xz ∂w + ∂u ∂ ∂x ∂z ∂x e
0 ∂ ∂y
0 ∂ ∂x ∂ ∂z
0
0 0 e ∂ u ∂z ⋅ v = ∂ ⋅ ue 0 w ∂ ∂y ∂ ∂z
∂
ε e = ∂ ⋅ N e ⋅ q e = B e ⋅ qe
(1.44)
Be
ahol a B e mátrix is felbontható a csomópontok szerint: B e = Bie B ej B ek … Feszültségmező leírására hasonlóan értelmezhető a 6 méretű feszültségvektor:
T ⇒ σT = σ x σ y σ z τ xy τ yz τ xz illetve a Hooke törvény alapján felírható csomóponti elmozdulásvektorral is:
σ e = D ⋅ ε ε = D ⋅ Be ⋅ qe
(1.45)
(1.46)
ahol D az anyagjellemzők mátrixa. Térbeli izotróp rugalmas anyagra: c1 c 2 c D= 2 0 0 0
c2 c1 c2 0 0 0
c2 c2 c1 0 0 0
0 0 0 c3 0 0
0 0 0 0 c3 0
0 0 0 0 0 c3
(1.47)
szerkezetű, ahol c1 = E
1 −ν (1 +ν )(1 − 2ν )
c2 = E
ν
(1 +ν )(1 − 2ν )
20
c3 =
E =G 2 (1 + ν )
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A későbbiek miatt érdemes a terhelési vektorokat is oszlopvektorba rendezni. A peremen és a térfogaton megoszló terhelések oszlopvektorai a következők: px ρk x e e p ⇒ p = p y , ρk ⇒ ρk = ρk y pz ρk z e
e
(1.48)
1.3.3. Az elem potenciális energiája A végeselem-módszer másik fontos lépése egy hibaelv megválasztása, amely lehetővé teszi az állandók meghatározását. Hibaelvként a potenciális energia minimuma elvhez tartozó variációs elvet választjuk. A továbbiakban előállítjuk az elem potenciális energiáját a közelítő mezők felhasználásával. Az előző fejezetben bevezetett vektorokkal az e jelű elemre: Π ep =
1 eT e e ε D ε dV − ∫ u eT pe dA − ∫ u eT ρk dV ∫ 2 Ve Ae Ve
(1.49)
p
majd a csomóponti elmozdulásokat bevezetve a végesdimenzióban felírt potenciális energia: Π ep =
1 q eT B eT De B e q e dV − ∫ q eT B eT p dA − ∫ q eT B eT ρk dV 2 V∫e Ae Ve p
végül a csomóponti elmozdulásvektort kiemelve az alábbi tömör alak írható [1]:
1 Π ep = q eT K eq e - q eT f e 2
(1.50)
Ebben a kifejezésben előforduló K e elemi merevségi mátrix és f e elemi terhelési vektor az alábbi szerkezetű: BieT K e = ∫ B eT De B e dV = ∫ B eTj De Bie Ve Ve
B
e i
K iie …dV = K eji
K ije K
e jj
… …
(1.51)
f e = f pe + fρek f pe =
∫N
eT
p e dA , fρek =
Aep
∫N V
eT
ρk e dA
(1.52)
e
Nyilvánvaló, hogy a merevségi mátrix szimmetrikus, továbbá a merevségi és terhelési mátrix csomópontok szerint particionálható. Az elem merevségi mátrixa az elem alakváltozási energiájával kapcsolatos. Ezért ha a q e az elem merevtestszerű mozgását írja le, akkor az alakváltozási energia zérus, egyébként pedig pozitív: q eT K e q e ≥ 0
(1.53)
amiből az következik, hogy a mátrix elfajuló.
21
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.4. Elemek illesztése, a rendszer merevségi mátrixa, redukált terhelési vektora, kinematikai peremfeltétel figyelembevétele Elemek illesztését a közös csomópontba eső elmozdulások azonosságával biztosítjuk. Legyen példaként az i jelű csomópontokba befutó elemek jele rendre e, e+1 illetve k, s. Ekkor az illesztés szerint [1]: q ie = q ie +1 = q ik = q is = qi
(1.54)
vagyis azt is mondhatjuk, hogy a megkülönböztető felső index elhagyható. A vizsgált rendszer teljes potenciális energiája az elemek energiáinak összegéből ( N e .az elemek száma) illetve a W k koncentrált csomóponti erők munkájából áll elő: e 1 1 Π p = ∑ q eT Kq e − q eT f e − W k = qT Kq − qT f (1.55) 2 e =1 2 Ez a kifejezés értelmezi a szerkezet q csomóponti elmozdulás vektorát, terhelési vektorát, valamint a szerkezet K merevségi mátrixát. A szerkezet csomóponti terhelési és elmozdulási vektora az illesztési feltétel alapján
N
Ne
∑q
f = … + qieT fie + qie +1T fie +1T + qisT fis + qTi fik
eT e
=
e =1
= … + qTi fie + fie +1 + fik + fis … = q Tf fi
azaz fi = ∑ fie
(1.56)
e∈i
tehát az összegzést mindazon elemekre el kell végezni, amelyek az i- jelű csomópontot tartalmazzák. Itt qT = q1T
qT2
… qTncs
(1.57)
T f T = f1T f2T … f ncs ncs pedig a szerkezet csomópontjainak száma. A szerkezet merevségi mátrixa az energiával kapcsolatos: Ne
∑q
eT
(1.58)
K e qe = q e Kq e
e=1
ahol K = [K ij ] i, j = 1, …, ncs K ij =
∑K
e ij
(1.59)
e∈i , j
Az összegzés alapján ij indexű blokk mindazon elemeknél szerepel, amelyek tartalmazzák egyidejűleg az i és a j jelű csomópontot. A csomóponti elmozdulások számítása a potenciális energia variációjának eltűnése alapján történik. Ennek megfelelően
22
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
δΠ p = δ qT
∂Π p ∂q
=0
azaz
δ q T ( Kq − f ) = 0
(1.60)
ahol δ q a kinematikai peremfeltételt kielégítő csomóponti elmozdulásvektor variációja. Legyen q j = q ju adott csomóponti elmozdulás, mely azt jelenti, hogy δ q j = 0 Ekkor a (1.60) egyenletben a j-edik blokksor 0-val szorzódik. Cseréljük ki a j-edik csomópont elmozdulását az ismeretlen R j külső terheléssel, amely nyilván reakció erő:
K j1 q1 K jj qju f j R j ismert
ismeretlen
(1.61)
Az ismert és az ismeretlen mennyiségek cseréjével kapjuk, hogy K11 … 0 … q1 f1 − K1 j q ju = (1.62) K i1 … E … -R j f j − K jj q ju ahol E egységmátrix. Az egyenletrendszer mérete megmarad, viszont a szimmetria elromlik. Ha a támasztóerő nem érdekes, akkor a j-edik blokksor főátlón kívüli elemei is nullázhatók, s így az együtthatómátrix szimmetriája helyreáll. Hasonló eredményre vezet a kinematikai előírás rugalmas megtámasztással való biztosítása. Ekkor megfelelően nagy rugóállandót kell alkalmazni. A Κ rangját megkapjuk ha képezzük az összes ismeretlen és a merevtestszerű mozgás szabadságfokának a különbségét. A csomóponti paraméterek megkötése révén elvégezve a kinematikai előírásokat, megszűnik a merevtestszerű mozgás lehetősége. A megoldandó egyenletrendszert ekkor is Kq = f
(1.63)
alakban szokás írni. Az együttható mátrix és terhelési vektor azonban már tartalmazza a kinematikai előírásokat is. Az egyenletrendszer jellemzői közül a nagy méretek miatt fontos az együttható mátrix zérustól különböző elemeinek elhelyezkedése. Az egyenletrendszer szalagszerkezetű, amelyet az alábbi ábra szemléltet:
1.5. ábra. K szerkezete
23
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
Ez azt jelenti, hogy a zérustól különböző blokkok egy adott sávba esnek, amelyet az elemen lévő sorszámkülönbség maximális értéke határoz meg. Nagysága alapvetően a sorszámozástól függ. Példaként tekintsük a következő végeselem felosztást és konstruáljuk meg a hozzá tartozó sematikus merevségi mátrixot
1.6. ábra. Téglalap felosztása 3 elemre
1.7. ábra. K mátrix szerkezete Itt látható, hogy a sávszélesség: főátló +5 elem. Ennél van kedvezőbb számozás is, mely egyúttal optimális számozást jelent.
1.8. ábra. Optimális számozás
1.9. ábra. K optimális szerkezete Az egyenletrendszer megoldása az un. direkt vagy iterációs eljárás alapján történik.
24
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.5. Egyváltozós feladatok Sok esetben a 3 dimenziós térben elhelyezkedő testet, a mechanikai viselkedés szempontjából, kinematikai és feszültségi hipotézisek felhasználásával egy, ill. kétváltozós feladatokká lehet áttranszformálni, vagyis redukálni. Ezzel a feladatok megoldhatósági körülményei egyszerűsödnek. A felvett hipotézisek alkalmazását a gyakorlat visszaigazolja. Igen lényeges osztályát jelenti az ilyen típusú szerkezeteknek, a rúdként tárgyalható modellek. A vizsgálatainkat síkbeli szerkezetekkel kezdjük.
1.5.1. Síkbeli rúdszerkezetek Rúdnak nevezzük azokat a testeket, amelyeknél a test egy kitüntetett térgörbére merőleges geometriai méretei lényegesen kisebbek a térgörbe irányában mérthez képest. Ha a térgörbe egyenes, akkor egyenes rudakról beszélünk. A test „merőleges metszetei” a rúd keresztmetszetét jelölik ki. Feltételezésünk szerint a keresztmetszet súlypontja a térgörbén helyezkedik el, amit tömören középvonalnak nevezünk. Elöljáróban síkbeli egyenes középvonalú és állandó keresztmetszetű (prizmatikus) húzottnyomott, hajlított-nyírt rudakat fogunk vizsgálni. A nyírási energia elhanyagolásával az. u. n. Bernoulli hipotézisű rudakhoz jutunk [1]. A rúdszerkezetek modelljeivel számos gyakorlati probléma szilárdságtani elemzése kényelmesen és nagy megbízhatósággal megoldható. A gépészetben az erőátviteli hajtóművek tengelyeinek méretezése, a gépállványok első durva méreteinek meghatározása, csarnokszerkezetek tervezése stb. feladatorientáltan elkészített végeselemes programok révén a mindennapos tervezői analízis eszköze. A bemutatott elmélet az elmélyültebb munkát, a mechanikai szemléletmód erősítését szolgálja.
1.5.1.1.
Bernoulli – féle hipotézis, variációs egyenletek
A vizsgált x–z síkban fekvő rúdszerkezet egy tetszőleges e jelű elemét a végein elhelyezkedő i és j csomópontokkal jellemezzük. A rúdhoz kötött helyi koordinátarendszert az i, j csomópontokon átmenő ξ tengely és a rúdkeresztmetszetben elhelyezkedő η, ζ főtengelyek alkotják. A rúd tengelye az x tengellyel β szöget zár be, hossza Le . A rúd elmozdulásánál feltételezzük, hogy a rúd keresztmetszete merőleges marad a meggörbült középvonalra, azaz érvényes a Bernoulli-féle hipotézis. Ekkor, szilárdságtani ismereteink alapján mondhatjuk, hogy a rúdban alakváltozási energiát csak rúdirányú feszültségek adnak. A keresztmetszet mentén húzás-nyomásból állandó, hajlításból lineárisan megoszló lefutású feszültség keletkezik. Ehhez tartozóan a rúdirányú elmozdulás a keresztmetszet egy tetszőleges P pontjában (lásd az 1.10.ábrát) a következőképp írható fel:
uP = uP (ξ ,η , ζ ) = u (ξ ) − w′ (ξ ) ζ ,
(1.64)
( ) , u (ξ ) , w (ξ ) a ξ ill. ζ irányú elmozdulás. dξ Az (1.64) összefüggés felhasználásával a tengelyirányú fajlagos nyúlás
ahol ( )′ =
d
ε ξ (ξ ) = ε (ξ ) = u′ (ξ ) − w′′ (ξ ) ζ ,
(1.65)
míg az egyszerű Hooke-féle anyagegyenlet alapján a normál feszültség
σ ( ξ ) = E ε (ξ ) ,
(1.66)
25
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
ahol E a Young modulus.
1.10. ábra. Síkbeli rúd elmozdulása, terhelése Jelölje a rúdon ható megoszló terhelést, hosszirányban pξ , keresztirányban pζ , melyeknek mértékegysége [N/mm]. A rúd végein − Fξ 0 , Fξ L rúderő, − Fζ 0 , Fζ L nyíróerő és − M η 0 , M η L hajlítónyomaték hat. A fenti terheléseket figyelembevéve, továbbá tekintettel, hogy az alakváltozási energia csak a σ (ξ ) feszültségből származik, a rúd teljes potenciális energiája két integrálon keresztül és a rúdvégeken ható koncentrált erők és nyomatékok terhelési munkájából áll össze. Πp =
1 2
∫ ∫ εξ E εξ dA dξ − ∫ ( u pξ + w pζ ) dξ L
A
L
− ( u ( L ) Fξ L − u ( 0 ) Fξ 0 + w ( L ) Fζ L − w ( 0 ) Fζ 0 − ( − w′ ( L ) M η L + w′ ( 0 ) M
)
(1.67)
)
η0
A minimum feltételt kijelölő δ Π p = 0 stacionaritási feltételből a δ u és δ w mezők függetlensége miatt - az u rúdirányú elmozdulás vonatkozásában részletezett módon -, az alábbi mezőegyenletek és peremfeltételek vezethetők le:
δu Π p = 0 = ∫ L
∫ δ u′ E (u′ − w′′ζ ) dA dξ − ∫ A
− (δ u ( L ) Fξ L − δ u ( 0 ) Fξ 0 − (δ u ( L ) Fξ L − δ u ( 0 ) Fξ 0 = ( AEu ′ − Fξ i
δ u pξ d ξ
L
) = ∫ δ u′ AE u′ dξ − ∫ δ u pξ dξ L
L
(1.68)
)
′ ) δ u 0 − ∫ δ u ( AEu′ ) + pξ dξ ⇒ L
L
A E u ′′ + pξ = 0,
(1.69)
26
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
(N − F )δu
L
ξi
0
= 0, N = AEu ′
(1.70)
A részletek mellőzésével a keresztirányú w elmozdulás vonatkozásában, az egyensúlyt kifejező alapegyenlet
Iη E wIV − pζ = 0,
(1.71)
és a dinamikai peremfeltételt adó variációs egyenletek
(F
ζ
)
− Fζ i δ w
L 0
= 0,
(M
η
)
− Mη i δ w′
L 0
= 0,
(1.72)
azaz
Fζ = − Iη E w′′′ , és Mη = − Iη E w′′
(1.73)
Látható, hogy megoszló terhelés hiányában az (1.69) és (1.71) mezőegyenletek lineáris u , ill. harmadfokú w polinommal elégíthetők ki. Amennyiben a rúd hossza mentén megoszló terhelések lineárisan változnak, akkor a mezőegyenletek partikuláris megoldását harmadfokú ill. ötödfokú polinom szolgáltatja. Ebben az esetben a végeselem közelítő mezői egyúttal pontos megoldások, vagyis jelen esetben a teljes potenciális energia minimuma elv egzakt megoldást szolgáltat a rúdszerkezet vonatkozásában. Nagy előnye a módszernek, hogy kis elmozdulások és alakváltozások feltételezése mellett, a statikailag többszörösen határozatlan szerkezetek minden nehézség nélkül vizsgálhatók. Figyelmet a kinematikai perem- és illesztési feltételek kielégítésére kell csak összpontosítani. A fenti (1.69)-(1.73) alatti variációs egyenletekből, peremfeltételekből következik, hogy az elemek közötti mezők folytonossági feltételek u mezőnél C 0 osztályú – azaz a függvény folytonos –, a w mezőnél C1 osztályú folytonosságot – azaz a derivált folytonosságát is megköveteljük. A síkban elhelyezkedő különböző irányítottságú elemek miatt a ξ − ζ helyi koordinátarendszerben értelmezett, csomópontonként megjelenő u, w, ϕη = − w′ elmozdulási paraméterek transzformációjára lesz majd szükség.
Elmozdulásmező közelítése
1.5.1.2.
A fenti levezetésből következik, hogy az elemen belüli elmozdulásmező u, w . Ezeket polinomok segítségével közelítjük. A polinomok tagjainak egy részénél az együtthatókat csomópontonként felvett két elmozdulási és egy szögelfordulási értékkel tudjuk kifejezni, ill. az inhomogén differenciálegyenletek partikuláris megoldásaihoz tartozó tagokat pótlólagos állandóként, paraméterként fogjuk a továbbiakban szerepeltetni. Definiálva az elem helyi koordinátarendszerben értelmezett q e általánosított csomóponti vektorát, az a e pótlólagos állandók vektorát, a felsorolt műveletek végrehajtása után az alábbi approximációhoz jutunk. Vagyis az elemen belüli elmozdulásvektor e
u u (ξ ) = = N e (ξ ) q e + N e (ξ ) a e , w ahol a csomóponti elmozdulásvektorhoz tartozó approximációs mátrix e
N e (ξ ) = N i (ξ ) N j (ξ ) , e
1 − ξ 0 ahol N i e (ξ ) = 1 − 3ξ 2 + 2ξ 3 0
(1.74)
(1.75) e
, − L (ξ − 2 ξ 2 + ξ 3 ) 0
27
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
e
ξ 0 0 N j (ξ ) = , ξ = ξ / Le . 2 3 2 3 3ξ − 2ξ L ( ξ − ξ ) 0 A pótlólagos állandókkal megszorzott approximációs mátrix e
L2 (ξ 2 − ξ e N (ξ ) = 0
) L (ξ 3
3
0
e
q e,T = qTi qTj , a e ,T = aTu
1.5.1.3.
)
−ξ
L4 (ξ 4 − 2 ξ 3 + ξ 2 ) L5 (ξ 5 − 3 ξ 3 + 2ξ 2 ) 0
e
0
(1.76)
aTw , qi e ,T = [u , w, − w′]i . e
e
Merevségi mátrix, redukált terhelési vektorok
Az (1.67) diszkretizálása után véges dimenziójú feladatot kapunk, azaz a diszkretizált teljes potenciális energia
ahol fqe(p) =
∫ e
L
e
pξ Ne,T (ξ ) p d ξ , fae(p) = ζ
fqe( p ) K eqa q e −2 e fa ( p ) K eaa a e pξ Ne,T (ξ ) p d ξ ζ
K eqq e K aq
1 Π p = Π p ( q , a ) = q e ,T a e ,T ( 2 e
∫ e
L
e
e qq
K =
∫B
e ,T
Le
AE 0 Be (ξ ) d ξ , Ke = (ξ ) aa 0 I η E
)
(1.77)
(1.78) e
∫ Le
AE 0 Be (ξ ) d ξ . Be,T (ξ ) 0 I η E
Itt 1 − L e B (ξ ) = 0
1 L
0
0
12ξ − 6 L2
4 − 6ξ L
0 0
6 − 12ξ L2
2 − 6ξ L
e
0
(1.79)
e
2ξ − L 3ξ 2 − L2 0 0 . (1.80) 2 2 3 2 3 0 12 (ξ − ξL ) + 2L 20ξ − 18ξL + 4L 0 A potenciális energiában szereplő, a helyi koordinátarendszerben értelmezett vegyes indexű merevségi mátrix, jelen esetben K eaq = K eqa,T = 0 . Be (ξ ) =
Az a e paraméterekhez tartozó K eaa merevségi mátrix és annak inverze zárt alakban felírható és így a ∂Π ep / ∂ a e = K aa a e − f ae( p ) = 0 minimum feltételből az a e kiszámolható. A számítások elvégzése után a pótlólagos állandók vektora a két mező vonatkozásában au
e ,T
aw
e ,T
e
pξ i = − 2 AE
1 pξ i − pξ 6 AEL
pζ i = 24 Iη E
1 pζ j − pζ i 120 Iη EL
(
(
j
)
,
)
(1.81) e
.
28
(1.82)
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A szimmetrikus K eqq merevségi mátrix az alábbi 0 AE / L 12 Iη E / L3 0 0 − 6 Iη E / L2 e K qq = − AE / L 0 − 12 Iη E / L3 0 − 6 Iη E / L2 0
− AE / L
0
0
− 6 Iη E / L2
0
− 12 Iη E / L3
4 Iη E / L
0
6 Iη E / L2
0
AE / L 2
2 − 6 Iη E / L 2 Iη E / L 0 2 6 Iη E / L 4 Iη E / L
0
0
6 Iη E / L
0
12 Iη E / L3
2 Iη E / L
0
6 Iη E / L2
e
A redukált csomóponti terhelési vektor az (1.78) alatti integrál kiszámítása után a következő összefüggések révén számolható, (az áttekinthetőség érdekében a mátrix elemeket vesszővel választjuk el): fq ( p ) e ,T = [ L [
pξ i 2
pζ i
(
1 pξ j − pξ i 6
+
(
) ] , L[
pζ i
+
2 pξ i
(
3 pζ j − pζ i 20
(
)
(
)
) ],
1 pξ j − pξ i ] , 12 2 3 pζ i pζ i 7 1 e L[ pζ j − pζ i ] , L2 [ pζ j − pζ i ] ] + + 2 20 12 20 A helyi koordinátarendszerben felírt diszkretizált potenciális energiát az elemek közötti elmozdulásmező folytonosságának biztosítása érdekében az x−z globális koordinátarendszerben értelmezett U, W elmozdulásokkal és a síkra merőleges ϕ y szögelforduláson keresztül lehet kifejezni. A helyi, rúdhoz kötött koordinátarendszerben lévő csomóponti általánosított elmozdulást a globálbeli értékeken keresztül az alábbi összefüggés révén fejezhetjük ki: − L2 [
+
1 pζ j − pζ i 30
) ] , L[
(
+
)
e
u cos β e qi = w = − sin β ϕη = − w′ 0 i
sin β cos β 0
0 0 1
e
e
U W ≡ Te q e , 0 i ϕY i
(1.83)
vagyis az elem csomóponti általánosított elmozdulásvektora e
qi T0 0 q = = 0 T0 q j e
e
e
qi e e q ≡ T q , j
(1.84)
ahol Te az elem transzformációs mátrixa. Ezek után az elem teljes potenciális energiája 1 Π p = Π p ( q e , ae ) = Π p ( q e , a e ) = qe ,T ( Te ,T K eqq Te q e − 2Te ,T fqe( p ) 2 1 = q e ,T K eqq q e − q e ,T f qe( p ) + ..., 2 ahol K eqq = Te ,T K eqq Te a globális rendszerbeli merevségi mátrix,
) + ... (1.85)
(1.86)
29
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
f qe( p ) = Te ,T fqe( p )
(1.87)
a globális rendszerbeli redukált csomóponti általánosított terhelési vektor. Ezek ismeretében az elemek csatolása az ismert szabályok alapján már könnyen elvégezhető. A hőhatás figyelembevétele
1.5.1.4.
A rúdban a hőmérséklet-megoszlást a θ = θ (ξ , ζ ) függvényen keresztül adjuk meg, α -a fajlagos hőtágulási együttható. Feltételezzük, hogy a rúdban a hőmérsékletmező az alábbi összefüggés alapján lineárisan változik
θ = θ (ξ , ζ ) = θi + (θ j − θi ) ξ / L + ζ ( ∆θi + ( ∆θ j − ∆θ i ) ξ / L
)
(1.88)
ahol θi az i-edik keresztmetszet súlypontjának hőmérsékletét, ∆θi pedig a hőmérséklet i-edik keresztmetszetbeli ζ menti lineáris változását jellemzi. Hőhatás esetén a keletkező normálfeszültség
σ (ξ ) = E {ε ξ (ξ ) − αθ (ξ , ζ )} .
(1.89)
A képlet szerint látható, ha pl. egy rúd meg van akadályozva a megnyúlásában (két vége merev lapokra támaszkodik), akkor egyenletesen melegítve a rudat θ (ξ , ζ ) = θ áll , alakváltozás nem lép fel, de a keresztmetszet menti állandó nyomó feszültség jön létre σ (ξ ) = − Eαθ áll . Ha a rúd egyik végét mozgásában nem akadályozzuk, akkor a hőmérséklet
emelkedésből származó fajlagos nyúlás ε ξ (ξ ) = αθ áll L / L azonos lesz az αθ áll értékkel, vagyis a rúdban nem lép fel hőfeszültség. Általánosan mondható, ha a test homogén, izotróp, és a testben a hőmérséklet lineárisan változik, továbbá a test szabadon tud terjeszkedni, akkor a testben a hőhatásból nem származnak feszültségek, annak ellenére, hogy a testben elmozdulások felléptek. A hőhatásból adódóan a teljes potenciális energia módosul. Szimmetrikusan felírva Πp =
1 2
∫ ∫( L
ε ξ − αθ
A
) E ( ε ξ − αθ ) dA dξ − ∫ ( u pξ
+ w pζ ) dξ − ... ,
(1.90)
L
majd a diszkretizálást elvégezve a redukált hőterhelési vektor: f q (θ ) e , T =
[
−
I Eα A Eα θi + θ j ) , η ( ( ∆ θ i − ∆ θ j ) , − Iη E α ∆ θ i , 2 L I Eα A Eα e ∆ θ j − ∆ θ i ) , Iη E α ∆ θ j ] θi + θ j ) , η ( ( 2 L
(1.91)
ill. a hőmérséklethez tartozó pótlólagos állandók vektora e
a(θ )
1.5.1.5.
e ,T
θ − θi , 0, 0, 0 . = α j 2L
(1.92)
Rugalmas ágyazás hatása
Érintkezésben álló testek egyikét gyakran un. Winkler típusú közeggel szokásos helyettesíteni. A vasúti sín, az utak betonburkolata, szerszámgépek szánrendszerének
30
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
vezetékei stb. ezzel a modellel jól megközelíthetők. A szóban forgó modellnél a testet egymástól függetlenül álló rugókkal helyettesítjük. A rugókban keletkező erő arányos az elmozdulással. Esetünkben a megoszló terhelés intenzitása pw = cζ w , ahol cζ un. rugóállandó, ágyazási tényező. Az egységnyi hosszra eső fajlagos alakváltozási energia 0.5 w pw = 0.5 cζ w2 . Ily módon a tartó menti rugalmas ágyazásból származó energia e U rug =
1 w cζ w dξ . 2 ∫L
(1.93)
A ζ tengely irányú w elmozdulás (1.74) alatti közelítésével, a behelyettesítés és integrálás után, L = Le figyelembevételével U
e rug
=U
e rug
K eqq 1 e ,T e ,T ( q , a ) = 2 qw aw K e aq e
e K qa e K aa rug
e
q ew e , a w
(1.94)
ahol q ew,T = q ew,,Ti q ew,,Tj , q ew,,Ti = [ w, − w′]i , a ew,T = [ a1 , a2 ]w , e
K eqq ,rug
156 − 22 L 54 cζ L 4 L2 13L = 156 420 szimm.
K
e aq , rug
1/ 60 = cζ L 101L 2520
K
e aa , rug
1/ 630 = cζ L L / 252
5
9
e
13L − 3L2 , 22 L 4 L2
L / 280 11 2 L 1260
(1.95)
1/ 60 − L / 280 109 23 2 , L − L 2520 2520
L / 252 23 2 L 2310
.
(1.96)
(1.97)
A felírásból következik, hogy ebben az esetben mivel K eaq kapcsoló mátrix nem zérus, a pótlólagos állandók hatásának eliminálása a 3.1 részben ismertetésre kerülő alszerkezettechnikánál ismertetett eljárás révén oldható meg. Ennek értelmében a redukált merevségi mátrix e e K eqq ,rug , red = ( K qq − K qa ( K eaa ) K eaq −1
)rug .
A rugalmas ágyazású tartóhoz rendelhető teljes potenciális energia
Πp =
1 2
∫
Iη E ( w′′ ) dξ +
L
2
1 2
∫ L
cζ w2 dξ − ∫ w pζ dξ − .... .
(1.98)
L
Képezve a Π p első variációját, a differenciálegyenletére a következőt kapjuk:
rugalmas
Iη Ew IV + cζ w = p ζ .
ágyazású
prizmatikus (1.99)
Ennek megoldása a
31
tartók
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
4α 4 =
cζ
(1.100)
Iη E
tag bevezetésével ξ
4α w = ∑ aiVi (ξ ) + V4 (ξ − τ ) pζ (τ ) dτ , cζ ∫0 i =1 4
(1.101)
ahol a Vi (ξ ) Krülov függvények
V1 = ch (αξ ) cos (αξ ) , V2 =
1 ch (αξ ) sin (αξ ) + sh (αξ ) cos (αξ ) , 2
1 1 sh (αξ ) sin (αξ ) , V4 = ch (αξ ) sin (αξ ) − sh (αξ ) cos (αξ ) , 2 2 alakban számolhatók. Ily módon a rugalmas ágyazású tartó w elmozdulásának polinomos közelítése már nem ad pontos megoldást. A végeselemes megoldás pontosítható egyrészt a pótlólagos állandók, másrészt a végeselemek számának emelésével. A fenti végeselemnél a felvett pótlólagos állandók száma kettő.
V3 =
1.5.2. Térbeli rúdszerkezetek Az előző fejezetben bemutatott elvek segítségével a térbeli rudak elmozdulásmezőn alapuló közelítése könnyen elvégezhető [1]. Mζ Fζ
Mη
ζ
η i
Fη
e
ξ Fξ
j
Mξ
1.11. ábra. Térbeli rúd igénybevétele A 1.11. ábra szerint húzás-nyomás, két tengelykörüli hajlítás-nyírás, továbbá csavarás jelenti az igénybevételeket. A Bernoulli-féle hipotézis és a St. Venant-féle szabad csavarás feltételezése mellett, az alakváltozási energia számításánál a nyírást és a keresztmetszet öblösödését megakadályozó gátlások hatását elhanyagoljuk. Nyírás figyelembevételekor az un. Timoshenko-féle rúdelmélet alapján szokás számolni, míg gátolt csavarás esetén a
32
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
vékonyfalú szelvényekre kidolgozott elméletek jönnek számításba. Ezek végeselemes tárgyalása az [1] irodalomban részletesen megtalálható. pζ z
ϕζ w v
ζ
pη
ϕη
η i
y
x e
ξ u
j
ϕξ pξ
1.12. ábra. Térbeli rúd terhelése, elmozdulása A fent említett elmozdulásra vonatkozó hipotézisek alapján a test tetszőleges P pontjának elmozdulása
uP = uP (ξ ,η , ζ ) = u (ξ ) − v′ (ξ ) η − w′ (ξ ) ζ + ϕξ′ (ξ ) Φ (η , ζ ) ,
(1.102)
vP = vP (ξ ,η , ζ ) = v (ξ ) − ϕξ (ξ ) ζ ,
(1.103)
wP = wP (ξ ,η , ζ ) = w (ξ ) + ϕξ (ξ ) η ,
(1.104)
ahol Φ = Φ (η , ζ ) deplanációs függvény, amelynek értéke végeselem-módszer révén is meghatározható [Páczelt I. – Szabó T.: Estimation of torsional rigidity by means of the finite element method, Acta Technica Acad. Sci. Hung., 104(1-3), 1991/92, p. 211-236.] A keletkező alakváltozások, az (1.102)alapján
ε ξ = ε ξ (ξ ,η , ζ ) = u′ (ξ ) − v′′ (ξ ) η − w′′ (ξ ) ζ + ϕξ′′ (ξ ) Φ (η , ζ ) , ∂Φ −ζ ∂η
(1.105)
ϕξ′ (ξ ) ,
(1.106)
∂Φ + η ϕξ′ (ξ ) . ∂ζ
(1.107)
γ ηξ = γ ηξ (ξ ,η , ζ ) = γ ζξ = γ ζξ (ξ ,η , ζ ) =
A csavarásnál értelmezést nyer a keresztmetszet I c csavarási keresztmetszet jellemző, amely az 2
∂Φ ∂Φ +η + −ζ Ic = ∫ { ∂ζ ∂η A integrál segítségével számolható ki.
2
} dA
(1.108)
33
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A csavarást jellemző ϕξ (ξ ) keresztmetszet merevtestszerű szögelfordulás mezőjét, hasonlóan 0
a rúdirányú elmozduláshoz C osztályú függvénnyel kell közelíteni. Az η − ξ síkbeli hajlításnál a v elmozdulásnak és annak ξ szerinti deriváltjának is folytonosnak kell lennie. Ily módon az elem belüli elmozdulásmező e
u ϕ ξ e u (ξ ) = = N e (ξ ) q e + N e (ξ ) a e , v w ahol a csomóponti elmozdulásvektor e
(1.109)
e
q e,T = qTi qTj = u , v, w, ϕξ , − w′, v′ i ,
(1.110)
továbbá a lineárisan megoszló terhelés hatását pontosan figyelembevevő pótlólagos állandók vektora e
(1.111) a e ,T = aTu aTv aTw . A síkbeli esetben bemutatott approximációs mátrix eredményeit felhasználva a térbeli esetre a közelítés függvényei könnyen felírhatók. A potenciális energia prizmatikus rudat feltételezve Πp =
1 2
∫ ∫ εξ E εξ dA dξ + 2 ∫ I G ϕξ dξ − ∫ ( u pξ +vpη + w pζ ) dξ 1
2
c
L
A
L
L
− ( u ( L ) Fξ L − u ( 0 ) Fξ 0 + v ( L ) Fη L − v ( 0 ) Fη 0 + w ( L ) Fζ L − w ( 0 ) Fζ 0
)
− ( ϕξ ( L ) M ξ L − ϕξ ( 0 ) M ξ 0 + v′ ( L ) M ζ L −v′ ( 0 ) M ζ 0 − w′ ( L ) M η L + w′ ( 0 ) M
η0
)
Diszkretizálás után ismét azt kapjuk, hogy K eaq = 0 , továbbá a pótlólagos állandók zárt alakban kiszámolhatók. A (12,12) méretű K eqq mátrix és csomóponti terhelési vektor az alábbiak szerint számítható ki: e
K
e qq
K ii K ij = , K ji K jj qq
(K )
e
ji qq
= ( K Tij ) , e
(1.114)
qq
34
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
ahol
(K )
e
ii qq
(K )
e
ij qq
( K jj )
e qq
0 0 AE / L 3 12 Iζ E / L 0 0 0 0 12 Iη E / L3 = 0 0 0 − 6 Iη E / L2 0 0 6 Iζ E / L2 0 0 0 0 − AE / L 3 − 12 Iζ E / L 0 0 0 0 − 12 Iη E / L3 = 0 0 0 0 6 Iη E / L2 0 0 − 6 Iζ E / L2 0
0
0
0
0
0
− 6 Iη E / L2
I cG / L 0
0 4 Iη E / L
0
0
0
0 0 AE / L 3 12 Iζ E / L 0 0 0 0 12 Iη E / L3 = 0 0 0 0 6 Iη E / L2 0 0 − 6 Iζ E / L2 0
0 0
0
− 6 Iη E / L2
− I cG / L
0
0
2 Iη E / L 0
0
0
0
0
0
6 Iη E / L2
I cG / L
0
e
2 6 Iζ E / L 0 , 0 0 2 Iζ E / L qq 0
0
0
e
2 6 Iζ E / L 0 , 0 0 4 Iζ E / L qq 0
0 − 6 I ζ E / L2 0 0
0
4 Iη E / L
0
0
0
4 Iζ E / L
e
. qq
Ezek integrálásánál felhasználtuk a Be (ξ ) = Bie (ξ ) B j e (ξ ) mátrixot, ahol e 1 0 0 0 0 0 − L 1 0 0 0 0 0 − L e Bi (ξ ) = , 1 0 2 (12ξ − 6) 0 0 0 (6ξ − 4) / L L 1 0 0 6 12 0 6 4 / 0 ξ ξ L − − − − ) ( ) 2 ( L e
1 0 0 0 0 0 − L 1 0 0 0 0 0 − L e . B j (ξ ) = 1 0 2 (6 − 12ξ ) 0 0 0 6ξ − 2) / L ( L 1 0 6 − 12ξ ) 0 − (6ξ − 2) / L 0 0 2 ( L
A redukált csomóponti terhelési vektor:
35
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
(
)
(
)
(
)
(
fqe( p )
e
pξ i 1 + pξ j − pξ i ] L[ 2 6 pη i 3 L [ 2 + 20 pη j − pη i ] pζ i 3 L [ 2 + 20 pζ j − pζ i ] 0 p − L2 [ ζ i + 1 pζ j − pζ i ] 12 30 p 1 L2 [ η i + pη j − pη i ] 12 30 = pξ i 1 L[ + pξ j − pξ i ] 2 3 pη i 7 L[ + pη j − pη i ] 2 20 pζ i 7 + pζ j − pζ i ] L[ 2 20 0 L2 pζ i + 1 p − p . [ ] j i ζ ζ 12 20 − L2 pη i + 1 p − p ] [ η j ηi 12 20
)
(
)
(
)
(
)
(
)
(
)
(
(1.115)
)
Az au e és a we pótlólagos állandók vektora azonos a síkbeli esetnél kapottal, az η irányú v elmozduláshoz tartozóan e
pη i 1 pη j − pη i . av = 24 Iζ E 120 Iζ EL Hőhatás esetén hasonló eredmények írhatók fel.
(
e ,T
)
(1.116)
A lokális rendszerbeli mennyiségek az alábbi ortogonális transzformáció segítségével nyernek átszámítást. Jelölje a rúdhoz kötött koordinátarendszer tengelyeinek irányába mutató egységvektorokat eξ , eη , eζ . A globális rendszerbeli Descartesi tengelyek irányába mutató egységvektorok legyenek e x , e y , e z . Ekkor áll e U e x ⋅ eξ V = e ⋅ e y ξ W i e z ⋅ eξ
e x ⋅ eη e y ⋅ eη e z ⋅ eη
e x ⋅ eζ u e e e y ⋅ eζ v = T0e ,T qelm ,i . e z ⋅ eζ w i
(1.117)
Hasonló összefüggés igaz a szögelfordulási paraméterekre is. Tehát az i-edik csomópont vonatkozásában fennáll, hogy
36
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
e
u U v V e W T0T 0 w e ,T e = ϕ = T0 qi . T ϕ ξ T 0 x 0 ϕ ϕ y η ϕζ ϕ z i i e
(1.118)
Az egész elem vonatkozásában a transzformáció e
e
e
qi T0T 0 qi e ,T e q = = =T q , T q q j 0 T0 j i e
(1.119)
azaz e
e
e
qi T0 0 qi e e q = = =T q . q q j 0 T0 j i Vagyis a globális rendszerbeli merevségi mátrix e
K eqq = Te ,T K eqq Te ,
(1.120)
(1.121)
míg a redukált csomóponti terhelési vektor f qe( p ) = Te ,T fqe( p )
(1.122)
egyszerű szorzással a lokális rendszerbeli értékekből kiszámolható.
1.6. Kétváltozós rugalmasságtani feladatok vizsgálata izoparametrikus elemekkel 1.6.1. Feladattípusok 1.6.1.1.
Síkalakváltozás (SA)
Amennyiben a vizsgált test geometriája és terhelése következtében létezik egy olyan irány, amely mentén a test pontjai nem mozdulnak el, valamint ezen kitüntetett irányhoz tartozó helykoordinátától a reá merőleges síkban fellépő elmozdulásvektor koordinátái függetlenek, síkalakváltozásról szokás beszélni. Ez az eset áll fenn például egy hosszú, nem feltétlenűl körgyűrű keresztmetszetű, nagynyomású cső esetén, mikor is a csőtest pontjai csak a tengelyre merőleges metszetben mozdulnak el. Legyen a kitüntetett irány z. Ekkor a szóbanforgó állapot kialakulásához az szükséges, hogy a térfogaton megoszló ρ k terhelésnek és az Ap felületen megoszló p terhelésnek ne legyen z irányú összetevője. Igy azután az elmozdulásmező és a terhelési függvények u = u( xy ) = ue x + ve y
ρ k = ρ k ( x, y ) = ρ ( k x e x + k y e y )
(1.123)
p = p ( x, y ) = p x e x + p y e y
37
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
alakban írhatók. Az A alakváltozási tenzor fentiekből adódóan εx 1 A = A( x, y ) = γ yx 2 0
1 γ xy 2
εy 0
0 εx 0 ⇒ ε = εy γ xy 0
(1.124)
míg a T feszültségi tenzor σ x τ xy 0 T = T ( x, y ) = τ yx σ y 0 0 0 σ z ahol σ z = ν (σ x +σ y ) Az ε z ≡ 0 miatt az alakváltozási energia számításánál csak a T tenzor síkbeli részével kell dolgozni, tehát σ z τ xy T = (1.125) τ yx σ y Így végül is síkalakváltozás esetén, homogén izotróp anyagot feltételve az anyagtörvény 1 −ν σ x ν E ν 1 −ν T ⇒ σ = σ y = (1 +ν )(1 − 2ν ) τ xy 0 0
1.6.1.2.
0 εx 0 ε y ≡ Dε 1 − 2ν γ xy 2
(1.126)
Síkfeszültségi állapot (SF)
A síkfeszültségi állapotot az jellemzi, hogy most a kitüntetett z irányra merőleges síkokon nem keletkezik σ z = τ xz = τ yz = 0 feszültség. Ehhez az szükséges, hogy a ρ k és p terhelési függvényeknek ne legyen z irányú összetevője. A vékony tárcsa középfelületére, ahol is z = 0 , a terhelésnek, melyet az oldalperemen írunk elő, négyzetes függvényként kell változni. Fentiek alapján a feszültségi tenzornak csak a síkbeli része lehet zérustól különböző σ x τ xy 0 T = τ yx σ y 0 = T ( x, y ) 0 0 0 Ismét homogén, izotróp anyagot tételezünk fel, így a z irányú fajlagos nyúlás
εz = −
ν
(ε x + ε y ) 1 −ν míg az A alakváltozási tenzor
38
(1.127)
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
εx 1 A = γ yx 2 0
0 εx (1.128) 0 = A( x, y ) ⇒ ε = ε y εy γ xy 0 ε z Tekintettel megint az alakváltozási energia kiszámítási módjára elegendő csak a tenzorok síkbeli részét megtartani: 1 γ xy 2
σ x E T ⇒ σ = σ y = 1 −ν 2 τ xy
1.6.1.3.
1 ν 0 εx 0 ε y = Dε ν 1 1 −ν γ xy 0 0 2
(1.129)
Általános Síkfeszültségi állapot (ÁSF)
A SF állapot szigorú kiindulási feltételeinek enyhítése céljából ez esetben feltételezzük, hogy a σ z mindenhol zérus, a σ x , σ y , τ xy a z –nek páros függvényei, a τ xz és a τ zy pedig a z – nek páratlan függvényei úgy, hogy közben a tárcsa alsó és felsőlapjain zérus értékűek. Itt jegyezzük meg, hogy egyes munkákban az SF illetve ÁSF feladatokat saját síkjukban terhelt lemezfeladatoknak is nevezik. Fenti feszültségi koordinátákra vonatkozó feltételek teljesüléséhez a terhelési függvények a térfogaton egyrészt
ρ k ( x, y , z ) = ρ k ( x, y , − z ) ≡ 0
(1.130)
másrészt a paláston p = px e x + p y e y + pz e z p x ( x, y , z ) = p x ( x, y , − z )
,
p y (x, y, z ) = p y (x, y, − z )
(1.131)
pz (x, y, z ) = − pz (x, y, − z ) alakúak kell, hogy legyenek. Az így megoldott terhelési feltételek mellett az egyes mechanikai mennyiségeket a b vastagság mentén integrálva átlagértékeket kapunk. Így értelmezhető az átlagos feszültségi és alakváltozási tenzor T=
1 Tdz , b (∫b )
(T ↔ A)
(1.132)
valamint az „átlagos” elmozdulás és terhelési vektor u=
1 u dz b (∫b )
(1.133)
p=
1 pdz = p x e x + p y e y b (∫b )
(1.134)
Így az integrálás elvégzésével az ÁSF állapotot is kétváltozósként lehet kezelni. A későbbiekben az átlagolásra utaló felülvonást elhagyjuk.
39
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.6.1.4.
Tengelyszimmetrikus alakváltozás (TSZ)
Ez esetben a 1.13.ábrán látható z tengelyű forgástest terhelése és megfogása független a kerületi irányban mért ϕ koordinátától
1.13. ábra. Egy forgásszimmetrikus test geometriája és egy tetszőleges meridiánmetszet mentén jelentkező elmozdulás koordináták Így az alkalmasan választott henger-koordináta-rendszerben a test tetszőleges pontjának elmozdulás vektora u = u ( R, Z )e R + w( R, Z )e z
(1.135)
alakú, azaz a vizsgálatokat egy tetszőleges meridiánmetszet mentén az Rz síkban kétdimenziós feladatként lehet elvégezni. Az alakváltozási és feszültségi vektorok ∂u ∂R ε R u σ R ε σ R ϕ ε= = , σ= ϕ εz σ z ∂w γ Rz ∂z τ Rz ∂u ∂w + ∂z ∂R között homogén izotróp anyagra az anyagállandók mátrixa
ν ν 1 −ν ν 1 −ν ν E D= ν ν 1 −ν (1 +ν )(1 − 2ν ) 0 0 0
0 0 0 1 − 2ν 2
(1.136)
(1.137)
teremt kapcsolatot. σ = Dε
(1.138)
40
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.6.2. Síkbeli elemek A valóságos mérnöki, szilárdságtani feladatok mindig a térbeli, háromdimenziós euklideszi térhez köthetők. Mégis, számos esetben a vizsgált test geometriai alakja, az anyagjellemzők és a testre működő külső erőrendszer tulajdonságai lehetővé teszik, hogy matematikailag a problémát kétváltozósként lehessen kezelni. A szilárdságtan tipikus kétváltozós feladattípusait a következő fejezet tárgyalja és az ott bemutatott formalizmusból látható majd, hogy a különböző feladattípusokhoz alapvetően azonos elemtípusokat lehet alkalmazni. A gyakorlati feladatok megoldásában különösen jól használhatók a lineáris és kvadratikus, három illetve négyszög geometriájú izoparametrikus elemek. Ez utóbbi, általános definició szerint azt jelenti, hogy az elem geometriai pontjait és az elem menti elmozdulás mezőt ugyanolyan, természetes koordináta-rendszerben adott interpolációs függvényekkel közelítjük [1], [2].
1.6.2.1.
Négycsomópontú elem
A 1.14.a.ábra egy konvex egyenesoldalú, négycsomópontú elemet mutat az xy globális koordináta-rendszerben, amelyet egy két egység élű négyzet-tartományra kívánunk leképezni. Ennek érdekében az x ( ξ , η ) és y (ξ ,η ) leképező függvényeket bilineáris alakban írjuk fel: x(ξ ,η ) = a1 + a2ξ + a3η + a4ξη = [1 ξ η ξη ] a = ϕ T (ξ ,η )a y (ξ ,η ) = b1 + b2ξ + b3η + b4ξη = ϕT (ξ ,η )b
(1.139)
1.14. ábra. Négyszög alakú végeselem a/ leképezés két egység élű négyzettartományra b, c/ nem konvex elemek amelyek nem biztosítják az egyértelmű leképezést ( det J ≤ 0 ) (1.139) -ben az ai és bi állandókat a csomópontok, azaz a sarokpontok x(ξi ,ηi ) = xi
,
y(ξi ,ηi ) = yi
koordinátái alapján lehet meghatározni. Az 1.14 b. ábrából kiolvashatóan a négy pont ξ ,η koordinátájának behelyettesítésével
41
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
x1 1 −1 −1 1 a1 x2 = 1 1 −1 −1 a2 x3 1 1 1 1 a3 x4 1 −1 1 −1 a4
(1.140)
ugyanez tömörebben x = Ga,
a = G −1x
Hasonlóképpen y-ra y = Gb
b = G −1y
Az állandókat visszaírva (1.139)-be a leképezés 4
x(ξ ,η ) = ϕT (ξ ,η )G −1x = ∑ N i (ξ ,η ) xi
( x ↔ y)
(1.141)
i =1
alakban áll elő, ahol az N i (ξ ,η ) un. alakfüggvények felépítése a következő: 1 1 (1 − ξ )(1 − η ) N 3 = (1 + ξ )(1 + η ) 4 4 1 1 N 2 = (1 + ξ )(1 − η ) N 4 = (1 − ξ )(1 + η ) 4 4 Könnyen ellenőrizhető, hogy a (1.142) függvények összege:
N1 =
(1.142)
4
∑ N (ξ ,η ) = 1 i =1
(1.143)
i
Mármost tekintettel az izoparametrikus elemek definíciójára a (1.142) függvények birtokában egyértelmű, hogy az elemek mentén az x irányú u és az y irányú v elmozduláskoordinátákat az 4
u = ∑ N i (ξ ,η )ui i =1
4
és
v = ∑ Ni (ξ ,η )vi
(1.144)
i =1
formulákkal közelítjük, ahol ui , vi konkréten az i-edik csomópontbeli x és y irányú elmozduláskoordinátákat jelenti. Végül e ponton belül néhány megjegyzés: A (1.142) formulákkal adott alakfüggvények elsőrendű folytonos deriváltakkal rendelkeznek. Ahhoz azonban, hogy ez biztosítsa a (1.144) elmozdulásmező folytonosságát az elem mentén, nyílván egy-egyértelmű leképzés szükséges az x,y és a ξ ,η koordináta-rendszerek között, melyhez teljesülni kell a ∂x ∂y ∂ξ ∂ξ > 0, det J = det ∂x ∂y ∂η ∂η azaz a J –Jacobi mátrix determinánsnak pozitívnak kell lennie. Ez azonban csak akkor lehetséges, ha a négyszög konvex, vagyis valamennyi belső szög 180o –nál kisebb (ld. 1.14..ábrák). Egyszerűen bizonyítható az elem teljessége
42
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
4
4
i =1
i =1
u = ∑ N i (ξ ,η )ui = ∑ N i (c0 + c1 + c2 yi ) =
(u ↔ v )
= ∑ N i c0 + ∑ N i yi ) c1 + ∑ N i yi c2 = c0 + c1 x + c2 y i i i ami fizikailag azt jelenti, hogy a közös oldallal rendelkező szomszédos elemek határai mentén-és így a teljes végeselemekkel behálózott kétdimenziós tartományon az elmozdulásmező folytonos. Az eddig ismertetett gondolatok valamennyi izoparametrikus elemtípusra érvényesek, azaz az elemháló sűrítésével a konvergencia kritériumok (a közelítő mezőhöz és az egzakt mezőhöz tartozó potenciális energiák különbségre a zérushoz tart) automatikusan teljesülnek. Ez is magyarázata az izoparametrikus elemek széleskörű alkalmazásának.
1.6.2.2.
Nyolc-csomópontú elem
A 1.15.ábra egy nyolc-csomópontú, görbeperemű izoparametrikus elemet mutatat, melynek nyolc alakfüggvényét [1], [2] alapján az alábbiak sorolják fel: 1 (1 − ξ )(1 − η )(−ξ − η − 1) 4 1 N3 = (1 + ξ )(1 + η )(ξ + η − 1) 4 1 N5 = (1 − ξ 2 )(1 − η ) 2 1 N 7 = (1 − ξ 2 )(1 + η ) 2 N1 =
1 (1 + ξ )(1 − η )(ξ − η − 1) 4 1 N 4 = (1 − ξ )(1 + η )(−ξ + η − 1) 4 1 N 6 = (1 − η 2 )(1 + ξ ) 2 1 N8 = (1 − η 2 )(1 − ξ ) 2 N2 =
(1.145)
1.15. ábra. Nyolc-csomópontú, izoparametrikus elem Természetesen a (1.145) alakfüggvények is teljesítik a (1.143) követelményt, továbbá itt is érvényes, hogy egy adott alakfüggvény az adott csomópontbeli helyettesítési értéke egy, míg minden más csomópontbeli helyettesítési érték zérus.
1.6.2.3.
A háromcsomópontú és a hatcsomópontú elemek
A 1.16. ábra a háromcsomópontú és a görbeperemű, hatcsomópontú elemeket mutatja:
43
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.16. ábra. a/ a háromcsomópontú, b/ a hatcsomópontú, görbeperemű izoparametrikus elem Az előbbiek alakfüggvényei: N1 = 1 − ξ − η
N2 = ξ
N3 = η
(1.146)
míg a hatcsomópontú elemé: 1 1 N 4 − N6 2 2 1 1 N3 = η − N5 − N 6 2 2 N5 = 4ξη
N1 = 1 − ξ − η −
N2 = ξ −
1 1 N 4 − N5 2 2
N 4 = 4ξ( 1 − ξ − η)
(1.147)
N 6 = 4η( 1 − ξ − η)
1.6.3. Az elem mechanikai jellemzői 1.6.3.1.
Az alakváltozási vektor előállítása
A korábbi fejezetekben mondottak értelmében az x, y globális koordináta-rendszerbeli mezők deriváltjait is approximálni kell. A két koordináta-rendszer közötti leképezést az előzőekkel összhangban ncs
x = ∑ N i (ξ ,η ) xi i =1
ncs
y = ∑ N i (ξ ,η ) yi
(1.148)
i =1
formulák adják, ahol ncs az adott elemmodell csomópontjainak száma (Előző példáinkban 3, 6 illetve 4 vagy 8). Kétváltozós feladatok esetén az elmozdulás koordináták közelítésére az ncs
u = ∑ N i (ξ ,η )ui i =1
ncs
v = ∑ N i (ξ , n)vi
(1.149)
i =1
összefüggések szolgálnak (1.149) szokás szerint tömörebben is írható! u N ue = = 1 v 0
0 N2
0
N1 0
N2
...
N ncs 0
e
0 e q ≡ N e (ξ ,η )q e N ncs
ahol
44
(1.150)
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
q eT = u1 , v1 ,...ui , vi ,...uncs , vncs az elem 2 ⋅ ncs méretű elmozdulásvektora. Ezek után felírandó az alakváltozási vektor
∂u ∂N i ui ∑ ∂x i ε x ∂x ∂N i ∂v (1.151) vi εe = ε y = ∑ = ∂y ∂y i γ xy ∂u ∂v ∂N i ∂Ni + + u v ∑ i i ∂x ∂y ∂x i ∂y alakban, amiből látszik, hogy a ehhez szükséges az alakfüggvények globál koordináta∂N i ∂N i rendszerbeli parciális deriváltjainak , ∂x ∂y ∂N i ∂N i ∂ξ ∂N i ∂η ∂x ∂ξ ⋅ ∂x + ∂η ⋅ ∂x = ∂N i ∂N i ∂ξ + ∂N i ⋅ ∂η ∂y ∂ξ ∂y ∂η ∂y
∂ξ ∂x = ∂ξ ∂y
∂η ∂N i ∂x ∂ξ ∂η ∂N i ∂y ∂η
(1.152)
meghatározása. (1.152) -et is célszerű tömörebben felírni
J11−1 J12 −1 ( ∂ G N i ) = J ( ∂ L N i ) = −1 ( ∂ L Ni ) J 22 −1 J 21 J-1 a Jacobi mátrix inverze ahol ∂ ( G ) a globálrendszerbeli deriváltak vektora −1
(∂ L )
(1.153)
a lokálrendszerbeli deriváltak vektora.
A ξ ,η rendszerbeli parciális deriváltakra analóg módon írható, hogy ∂N i ∂x ∂ξ ∂ξ = ∂N i ∂x ∂η ∂η
∂y ∂N i ∂ξ ∂x ∂y ∂N i ∂η ∂y
(1.154)
vagyis
J11 J 21
( ∂ L Ni ) = J ( ∂ G Ni ) =
J12 (∂ N ) J 22 G i
(1.155)
így azután a Jacobi mátrix (1.148) felhasználásával
∂ ∑ Ni i xi ∂ξ J = 2,2 ∂ ∑ Ni i xi ∂η
∂∑ Ni
∂N i yi ... ∂ξ ... x1 ∂ξ xi = ∂ ∑ N i ∂N i ... ... x i yi ∂η ncs ∂η i
45
y1 yi yncs
(1.156)
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
módon számítható. Ezek után (1.153) figyelembevételével az x és y szerinti parciális deriváltakat ξ és η szerinti parciális deriváltakból az alábbi módon lehet előállítani: ∂N i ∂N ∂Ni = J11−1 + J12 −1 i ∂x ∂ξ ∂η ∂N i ∂N ∂Ni = J 21−1 + J 22 −1 i ∂y ∂ξ ∂η
(1.157)
Ez azt jelenti, hogy így előállítható a következő formula által definiált B e el- mozdulásalakváltozás transzformációs mátrix, amelyet szorozva a q e elem csomóponti elmozdulásvektorral, közvetlenül számítható az
∂N1 ∂x εe = 0 ∂N1 ∂y
∂Ni ∂x
0 ∂N1 ∂y ∂N1 ∂x
...
0 ∂Ni ∂y
0 ∂Ni ∂y ∂N i ∂x
u1 v 1 ... = B e (ξ ,η )qe ui vi
(1.158)
elem alakváltozási vektora.
1.6.3.2.
Az elemi merevségi mátrix és a redukált terhelési vektor számítása
Az előző pontokban ismertetett kétdimenziós feladattípusok sajátosságait figyelembe véve a potenciális energia kifejezése mátrixos formában a (1.150) elmozdulás mező közelítéssel és az ε alakváltozási vektort értelmező (1.151) formula felhasználásával 1 Π ep = q eT ∫ BT DB b dA q e − qeT (fεe + f pe + fqke ) 2 ( Ae )
(1.159)
ahol
Ke =
∫B
T
DB b dA
- elemi merevségi mátrix, az
Ae
∫B
fεe =
T
Dε 0 bdA
- kezdeti alakváltozásból,
Ae
f pe =
∫ N pb d Γ T
Γ
f ρek =
- felületi terhelésből, és
e
∫N
T
ρ kb dA
- a térfogati terhelésből
e
A
számítandó terhelési vektor. (1.159) felírásakor kihasználtuk, hogy ÁSF esetén az elemi térfogat, b vastagságú tárcsa esetén dV = b dA
(1.160)
Ugyanez SA esetén b = 1 egységnyi szeletre, míg TSZ állapotváltozáskor (1.159)-be b = 2π R
helyettesítésével kapjuk hogy
46
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
dV = 2π RdA .
1.6.3.3.
(1.161)
Numerikus integrálás
Az (1.159)-es kifejezésben szereplő felületi integrálok a ξ ,η változók függvényei, így az integrálást ξ és η szerint -1 és +1 intervallumra vonatkozóan kell elvégezni. Mivel az elemi felület dA = dx dy = det Jd ξ dη
formában írható. Ezért az elemi merevségi mátrix
K = e
+1 +1
∫ ∫ B DBb det Jdξ dη T
(1.162)
−1 −1
illetőleg a térfogati terhelés csomóponti vektor (és hasonlóan a kezdeti alakváltozásból adódó csomóponti vektor) 1 1
f ρek =
∫ ∫N
T
ρ k b det J dξ dη
(1.163)
−1 −1
szerint állítható elő. A peremen működő terhelés számításához és a Jacobi mátrix determinánsára van szükség. Nem részletezve a levezetést csak utalva a [1]-es hivatkozásra, például egy elem η = 1 izoparametrikus koordinátával jellemezhető oldala mentén müködő x irányú px megoszló terhelésből
p f = ∫ N pb d Γ = ∫ NT (ξ ,η = 1) x b det J Γ dξ 0 Γe −1 1
e p
T
(1.164)
elemi redukált csomóponti tehervektor származtatható. Kérdés ezek után a (1.162), (1.163) típusú integrálok előállítása, melynek széles körben alkalmazott módszere a Gauss-féle numerikus integrálási technika. Ennek értelmében a (8.17) formula integranduszban szereplő mátrix-szorzatot F e ( x, y ) -al jelölve a Ke =
∫B
T
( x, y )D( x, y )B( x, y )b( x, y )dA ≡
Ae
=
∫ F ( x, y )dA = e
Ae
1 1
NG NG
−1 −1
i =1 j =1
∫
(1.165)
e e e e i j det J (ξ i ,η j )F (ξi ,η j ) míg például a ∫ F ( x(ξ ,η ), y (ξ ,η ) det J (ξ ,η )dξ dη = ∑∑WW
(1.164) szerinti perem menti integrál 1
f pe = ∫ NT pbd Γ = ∫ NT (ξ ,η = 1)p(ξ )b(ξ ,η = 1) det J Γ (ξ ,η = 1)d ξ Γe
(1.166)
−1
formában számolható. Fenti képletekben NG a ξ illetve η irányban felvett Gauss integrációs pontok számát Wi , W j pedig az integrációs súlyfokokat jelenti, amelyeket számszerűen a 1.6.3.1 és 1.6.3.2 táblázat is bemutat.
47
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.6.3.1. táblázat: Négyszögelem Gauss pontok és súlyok [2]
48
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.6.3.2. táblázat: Háromszögelem Gauss pontok és súlyok [3]
1.6.4. Lemezek 1.6.4.1.
Reissner-Mindlin-féle lemez modell
Minden olyan vékony, síkfelületekkel határolt testet melyeknél egyértelműen kijelölhető egy sík középfelület lemeznek nevezünk. Ha a testet terhelő erőrendszer olyan, hogy a test elmozdulásának jelentős részét a középfelületre merőleges irányban okozza, akkor mechanikai értelemben is lemezfeladatról van szó, egyébként síkfeladatról szokás beszélni. A lemezelmélet olyan hipotézisekkel él, amellyel az eredeti háromváltozós feladatot kétváltozósra lehet viszavezetni. Természetesen az ilyen szerkezetek megtámasztásánál, csatlakozásánál, illetve az erőbevezetés helyein a tényleges feszültségállapot térbeli, amely az alkalmazott hipotézisekkel nem írható le. Ilyenkor átmeneti és térbeli elemeket szokás alkalmazni. A továbbiak a klasszikus lemezmodell vizsgálatával foglalkoznak [3].
1.17. ábra. Lemez jellemző adatai Legyen a b vastagságú lemez középfelülete az xy síkban (1.17. ábra). Kis lehajlású lemezról beszélünk, ha a z irányú w elmozdulás kisebb mint a lemezvastagság kb. 20 százaléka. Hajlított lemezek, nyírási energiáját is figyelembevevő elméletét, Reissner-Mindlin
49
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
elméletnek is szokás nevezni [3]. Ezek szerint egy térbeli pont elmozdulása az elmodulási hipotézis alapján származtatható.
1.18. ábra Elmozdulás hipotézis Elmozdulási hipotézis szerint a középfelület csak z irányba mozdul el, továbbá a középfelület normálisa merevtestszerűen elfordul (egyenes marad és hossza nem változik, azonban nem marad merőleges a középfelületre). Feltételezzük, hogy az elmozdulás és alakváltozás kicsi. A normális x és y tengely körüli szögelfordulása ϕ x , és ϕ y szintén kicsi. Ekkor az elmozdulás egy tetszőleges P pontban az alábbi alakban írható (lásd 1.19. ábrát, melyen a szögelfordulást felnagyítva szemléltetjük): u = ϕ y ( x, y ) ⋅ z v = −ϕ x ( x, y ) ⋅ z
(1.167)
w = w ( x, y )
A középfelülethez kötött ϕx , ϕ y , w ( x, y ) meghatározhatjuk az alakváltozás jellemzőit is.
segítségével bármely térbeli pontban
∂ϕ y ∂u = zϕ′y = z ⋅ ∂x ∂x ∂ϕ x ∂v εy = = −z ⋅ ∂y ∂y ∂w εz = =0 ∂z
εx =
γ xy =
∂ϕ ∂ϕ ∂u ∂v + = z ⋅ y − x ∂y ∂x ∂x ∂y ∂w ∂u ∂w γ xz = + = + ϕy ∂x ∂z ∂x ∂w ∂v ∂w γ yz = + = − ϕx ∂y ∂z ∂y
Érdemes a z-től függő, és a z-től független jellemzőket külön mátrixba rendezni:
50
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
∂ϕ y εx ∂x ∂ϕ ε = ε y = z ∂yx = z ⋅ κ γ xy ∂ϕ y ∂ϕx − ∂x ∂y
(1.168)
κ
ahol κ a görbületek oszlopvektora: 0 0 κ = 0 − ∂∂y 0 − ∂∂x ∂ε
∂ ∂x
w 0 ϕx ∂ ϕ y ∂y
(1.169)
u
A z-től független szögtorzulások pedig a következő kételemű vektorba rendezhetők: ∂
γ xz ∂x γ= =∂ γ yz ∂y
w 1 ϕx = ∂ γu −1 0 ϕ y 0
∂γ
(1.170)
u
Látható, hogy az elmozdulási hipotézishez tartozó szögtorzulás két tagból áll. Az első tag a normális középfelülettel együtt történő mozgása miatt, a második pedig az ahhoz képest jelentkező szögelfordulás miatt van. Ezt szemlélteti az alábbi ábra példaként γ yz esetén:
1.19. ábra. Szögelfordulás és szögtorzulás Természetesen a valóságban a perem síkbeli erőkkel való terheletlensége miatt a szögtorzulás a peremen zérus kell, hogy legyen. A normális közelítő és egzakt alakját szemlélteti a következő ábra:
1.20. ábra. Normális elfordulása
51
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A klasszikus lemezelmélet másik hipotézise a feszültségi állapottal kapcsolatos. Ennek alapján a középfelületre merőlegesen ébredő normálfeszültség elhanyagolható: σz ≅ 0 Nyilvánvaló, hogy az ε z = 0, σ z = 0 egyidejű feltételezése ellentmond a Hooke-féle anyagtörvénynek. A Hook-törvény alapján származtathatók az elmozdulás hipotézishez tartozó feszültség koordináták is. A z-től függő koordináták: 1 E σ x −νσ y ) σ x = ( (ε x +νε y ) E 1−ν 2 1 E ε y = (σ y −νσ x ) σ y = (νε x + ε y ) E 1−ν 2 τ E E γ xy = xy τ xy = Gγ xy = γ xy mivel G = G 2 (1 +ν ) 2 (1 +ν )
εx =
Mátrixos alakban pedig: σx 1 ν 0 ε x ∼ E ∼ z 1 0 σ = σy = ν ε = D ⋅ ε = D ⋅κ y 1 − ν2 1 −ν 0 0 2 γ xy τ xy A z-től független nyírófeszültségek pedig:
(1.171)
τ xz = kG γ xz τ yz = kG γ xz mátrixosan írva:
τ xz (1.172) τ = = kGγ τ yz A nyírófeszültségekben szereplő k az un. nyírási tényező,mely abból a feltételezésből határozható meg, hogy a közelítő konstans nyírófeszültséghez és az egzakt nyírófeszültséghez tartozó alakváltozási energia megegyezik. Értéke: 5/6.
1.21. ábra. Nyírási tényező származtatása Az előbbiekben bevezetett alakváltozási és feszültségi jellemzőkből már számítható a lemez alakváltozási energiája.
52
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A továbbiakban a külső erők virtuális munkájának felírása céljából tekintsük a lemez peremének egy darabját. A peremhez kötött helyi (s,n) koordinátarenszerben ébredő feszültségek közül σn , τsn a vastagság mentén lineáris, τ zn pedig konstans. Külső/belső erőrendszer redukálása a középfelületre
1.22. ábra. Feszültségek, élerők, élnyomatékok A ϕs , ϕn , w elfordulás és elmozdulás koordináták munkájának felírásához redukálni kell a feszültségeket is a középfelület pontjaiba. A σ n -ből és τsn -ből származó nyomaték M n és M sn , míg a τnz -ből származó eredő Qn az alábbi formulákkal számíthatók: M n = ∫ zσ n dz
M sn = ∫ zτsn dz
Qn = ∫ τ zn dz
b
b
b
Meg kell jegyezni, hogy a feszültségkoordinátákhoz hasonlóan ezen élnyomatékok és élerő az x,y koordinátarenszerben értelmezett megfelelő tenzorokból is származtathatók. A külső erők virtuális munkája ennek alapján a következő: Wk = ∫ wpdA − ∫ M n ϕs ds + ∫ M sn ϕn ds + ∫ Qn wds A
Γ
Γ
Γ
ahol p ( x, y ) a felső illetve alsó lapon működő normális irányú terhelés eredője, A a középfelület Г pedig annak a pereme. A lemez megtámasztásától függően az alábbi peremfeltételeket szokás megkülönböztetni (lásd 1.23.ábrát ): Befogás esetén:
w = 0 , ϕs = ϕn = 0
szabad perem esetén:
M n = M sn = Qn = 0
egyszerű alátámasztás esetén:
w = 0 , M n = M sn = 0 w = ϕn = 0 , M n = 0
vagy pedig:
53
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.23. ábra. Megtámasztási módok 1.6.4.2.
Izoparametrikus lemezelem
1.6.4.2.1.
Potenciális energia
Az előző pontban megismert összefüggések fehasználásával felépíthető a REISSNER – MINDLIN féle lemez teljes potenciális energiája [4]: Πp =
1 1 εT ⋅ σ dzdA + ∫∫γ T ⋅ τ dzdA − Wk ∫∫ 2 Ab 2 Ab
Πp =
1 0 1 T 2 ∼ 1 κ ∫z D dz κ dA + ∫γ T k G b ⋅ γ dA − Wk ∫ 2A b 2A 0 1
∫
(1.173)
Dny
Dh = D z 2 dz b
1 1 Π p = ∫κT ⋅ Dh ⋅ κ dA + ∫γ T ⋅ Dny ⋅ γ dA − Wk 2A 2A ahol
1 ν b3 E b3 Dh = D = ν 1 12 12 (1 − ν 2 ) 0 0
, 1−ν 2 0 0
1 0 Dny = k G b 0 1
(1.174)
továbbá Wk = ∫ wpdA − A
∫M
0 n
ϕs ds +
Γp
∫M
Γp
0 sn
ϕn ds +
∫ Q wds 0 n
(1.175)
Γp
Itt a 0-val indexelt élnyomatékok és élerő előírt értékek a peremen. A perem normálisának és érintőjének irányvektora legyen:
nx en = ny
,
−n es = y nx
54
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
ekkor
ϕn = ϕ x nx + ϕ y n y , ϕsn = −ϕ x n y + ϕ y nx 1.6.4.2.2.
Izoparametrikus elem
Izoparametrikus elemet választva formálisan az elem lehajlási és szögelfordulási mezőit kell a megszokott módon alakfüggvények és csomóponti paraméterek szorzataként approximálni a C 0 osztályú folytonosság biztosításához [1]: e
w e u = ϕ x = N eqe ϕ y
(1.176)
ahol
q eT = q1T
qTne cs
e
, qieT = wi
ϕ xi
ϕ yi
e
Ebből kiindulva képezhetők az alakváltozás jellemző vektorai:
κ e = ∂ ε ⋅ u e = ( ∂ ε ⋅ N e ) ⋅ q e = B eh ⋅ q e
γ e = ∂ γ ⋅ u e = ( ∂ γ N e ) ⋅ q e = B eny ⋅ q e
ahol az elmozdulás-alakváltozás mátrixok felépítése a következő:
B eh = [… Bi …] h e
(B )
e i h
0 0 ∂N = 0 − i ∂y ∂N i 0 − ∂x
∂N i ∂x ∂N i ∂y
e
B eny = [… Bi …] ny e
(B )
e i ny
∂N i ∂x = ∂N i ∂y
0 − Ni
Ni 0
e
Az alakfüggvények deriválása a síkfeladatoknál leírtakhoz hasonlóan történik Az elem potenciális energiája: det J d ξd η
Πp e
↑ 1 1 e = qeT ⋅ ∫ B eT ⋅ D ⋅ B dA ⋅ q e + q eT h h h 2 2 Ae K eh
K e = K eh + K eny
55
∫
Ae
e e B eT ny ⋅ Dny ⋅ B ny dA⋅ q − Wk
K eny
(1.177)
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1 Π p = q eT K e q e − Wk 2 Jól látszik, hogy a hajlításból és a nyírásból származó energiákhoz különböző merevségi mátrixok rendelhetők. A külső erők virtuális munkája is kifejezhető a csomóponti elmozdulás és szögelfordulás paramétereivel:
e eT Wk = q ∫ N eT Ae
p 0 dA − N eT ∫ Γep 0
Qn0 0 ds + 0
0 0 + q eT ∫ N eT −n y M n0 − nx M sn0 ds = q eT f e Γep n ny x Végül az elem potenciális energiája a szokásos tömör alakban írható: 1 Π p = q eT K e q e − q eT f e 2 A gyakorlatban az un. Lagrange-típusú approximációt felhasználó 4, 8, illetve 9 csomópontú izoparametrikus négyszög alakú elem a legelterjedtebb.A négycsomópontú elem egy bilineáris approximációt tartalmaz, alakja az alábbi:
1.24. ábra. 4 csomópontú elem ahol az alakfüggvények a következők: 1 (1 + ξξi )(1 + ηηi ) i = 1,..., 4 4 A 9 csomópontú négyszögelem oldalélei már parabolák is lehetnek. Ekkor az approximáció bikvadratikus, teljes másodfokú polinom. A 9 csomópontú Lagrange-féle elem alakfüggvényei a következők. A sarokpontokban: Ni =
1 ξηξi ηi (1 + ξξi )(1 + ηηi ) i = 1,..., 4 4 az oldalfelező pontokban: Ni =
Ni =
1 ξξi (1 + ξξi ) (1 − η2 ) i = 5, 7 2
56
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1 ηηi (1 + ηηi ) (1 − ξ2 ) i = 6,8 2 a kilencedik, középső pontban: Ni =
N9 = (1 − ξ2 )(1 − η2 )
Ezen kívül használatos a szűkített approximációt tartalmazó 8 csomópontú elem is (un. serendipity elem).
1.25. ábra. 8 csomópontú elem Ekkor az alakfüggvények a következők: 1 (1 + ξξi )(1 + ηηi )( ξξi + ηηi − 1) , i = 1, 2, 3, 4 4 1 , i = 5, 7 Ni ( ξ, η) = (1 − ξ2 ) (1 + ηηi ) 2 1 , i = 6,8 Ni ( ξ, η) = (1 − η2 ) (1 + ξξi ) 2 Ezen belűl szokás csak a szögelfordulás közelítését szűkíteni (un. Heterosis elem). Ezek a verziók az elem b vastagságának csökkenésekor fellépő numerikus nehézségeket kívánják leküzdeni. Ekkor ugyanis a nyírási energia tart a zérushoz ami a megfelelő merevség mátrix kondicióját rontja. Számítási tapasztalatok azt mutatják, hogy a különböző elemekhez a numerikus integrálást is megfelelően kell megválasztani. Ni ( ξ, η) =
1.6.4.2.3.
Kirchhoff-féle vékony lemez
A korábbi elmozdulási feltételek itt is érvényben vannak, igaz, hogy γ xz = γ yz = 0 . Azaz a normálisok normálisok maradnak a deformáció során, nem úgy mint a vastag lemez elméletnél feltételeztük. A normális szögelfordulása ezért most kifejezhető a lehajlás deriváltjaival: ϕx =
∂w ∂w , ϕy = − ∂y ∂x
Így az elem alakváltozása már másodrendű deriváltakat is tartalmaz, ami miatt C1 -rendben folytonos elemeket kell alkalmazni. A vékony lemezekre számos elemfajta került kidolgozásra. Ezek ismertetésétől most eltekintünk.
57
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.7. Térbeli elemek A térbeli elemek közül a legelterjedtebbek az izoparametrikus elemek, amelyek sokféle geometriai alakot ölthetnek. Ezek közül a leggyakrabban alkalmazottak az egyenes illetve görbült élekkel rendelkező hatlapú (hexahedron) „tégla”-, négyoldalú (tetraéder) „gúla”-, vagy ötoldalú (pentaéder) „ék alakú” elemek. Ha minden él csak egyenes vonalú lehet, akkor az elemen belül az approximáció lineáris, ha mindegyik görbülhet is, akkor, pedig legalább kvadratikus. Természetesen görbült- és egyenes él, azaz lineáris- és kvadratikus közelítés (approximáció) vegyesen is előfordulhat egy-egy elemen belül. A görbült elem élein, a végpontokon kívül, a felező pontokban is van csomópont. Megjegyezzük, hogy magasabb approximáció alkalmazása esetén, az éleken akár kettő vagy három közbenső csomópont is előfordulhat, valamint további csomópontok lehetnek az oldallapokon és az elem belsejében is [1]. Itt részletesen ismertetjük a 8 csomópontú lineáris-, a 20 csomópontú kvadratikus hexahedron, valamint a 10 csomópontú kvadratikus tetraéder elem alakfüggvényeit. Ezt követően a merevségi mátrix és terhelési vektor előállítása a különböző elemekre egységes formalizmus szerint történik. Az elem csomópontjai három - x,y,z irányú elmozdulási - szabadságfokkal rendelkeznek, ezért terhelés a csomópontokban is csak ugyan ilyen irányú erő lehet.
1.7.1. Az alakfüggvények és a közelítés A 1.26. ábra egy térbeli nyolc-csomópontú izoparametrikus elemet szemléltet. Az elem alakfüggvényei kvázi lineáris alakban adhatók meg [3]: N1 =
1 1 (1 − ξ )(1 −η )(1 − ζ ) N 2 = (1 + ξ )(1 − η )(1 − ζ ) 8 8
(1.178)
N3 =
1 1 (1 + ξ )(1 + η )(1 − ζ ) N 4 = (1 − ξ )(1 + η )(1 − ζ ) 8 8
(1.179)
N5 =
1 1 (1 − ξ )(1 − η )(1 + ζ ) N 6 = (1 + ξ )(1 − η )(1 + ζ ) 8 8
(1.180)
1 1 (1 + ξ )(1 + η )(1 + ζ ) N8 = (1 − ξ )(1 + η )(1 + ζ ) 8 8 De tömören is felírhatjuk: N7 =
Ni =
1 1 (1 + ξξi )(1 + ηηi )(1 + ζζ i ) ahol ξi , ηi , ζ i = 8 −1
(1.181)
a sarok pontok koordinátái.
1.26. ábra A nyolc-csomópontú izoparametrikus térbeli elem
58
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.27. ábra. A húsz-csomópontú hexahedron elem Az 1.27. ábrán látható húsz-csomópontú elem alakfüggvényeit tömör formában adjuk meg a sarok és élfelező pontokban: • sarokpontokban i = (1, 2, 3, 4, 9, 10, 11, 12 ) Ni = •
1 (1 + ξξi )(1 + ηηi )(1 + ζζ i )(ξξi + ηηi + ζζ i ) 8
az élek felezőpontjában
1 (1 + ηηi )(1 + ζζ i ) (1 − ξ 2 ) 4 1 • ahol η = 0 és i = ( 6, 14, 16, 8 ) Ni = (1 + ξξi )(1 + ζζ i ) (1 − η 2 ) 4 1 • ahol ζ = 0 és i = (18, 19, 20, 17 ) Ni = (1 + ξξi )(1 + ηηi ) (1 − ζ 2 ) 4 A kereskedelmi szoftverekben a térbeli elemek közül leggyakrabban a tetraéder elemek fordulnak elő. Az irodalom a lineáris 4 csomópontú változat alkalmazását nem ajánlja, mert az eredmények esetenként túl nagy hibával terheltek. Ezért itt a 10 csomópontú kvadratikus változatot mutatjuk be (1.28. ábra), amely tulajdonképpen a 4 csomópontú változat alakfüggvényeinek bővítésével, éspedig az adott sarokpontba befutó közbenső csomópontok alakfüggvényei 0.5 szörösének levonásával állnak elő. •
ahol ξ = 0 és i = ( 7, 15, 13, 5 )
Ni =
1.28. ábra. Kvadratikus tetraéder elem Az alakfüggvények: N1 = 1 − ξ − η − ζ −
1 ( N5 + N7 + N10 ) 2
N2 = ξ −
59
1 ( N 5 + N 6 + N8 ) 2
(1.182)
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
N3 = η −
1 ( N6 + N7 + N9 ) 2
1 ( N8 + N9 + N10 ) 2
(1.183)
N 7 = 4η (1 − ξ − η − ζ )
(1.184)
N4 = ζ −
N 5 = 4ξ (1 − ξ − η − ζ )
N 6 = 4ξη
N10 = 4ζ (1 − ξ − η − ζ ) N8 = 4ξζ N9 = 4ηζ (1.185) Az elem geometriájának leképzése a lokális a ξ ,η , ξ természetes koordináták a globális x, y, z koordináták között az alakfüggvények segítségével minden fent bemutatott elemre egységesen írható fel: ncs
x ( ξ ,η , ζ ) = ∑N i (ξ ,η , ζ ) xi
x↔ y↔z
(1.186)
i =1
ahol ncs a vizsgált elem csomópontjainak a száma. Az izoparametrikus elemekre jellemző módon az elmozdulás közelítését approximációját ugyan ezekkel az alakfüggvényekkel írjuk le. Ezért az alakfüggvényeket szokás interpolációs függvényeknek is nevezni. ncs
u = ∑N i ui u ↔ v ↔ w
(1.187)
i =1
Mátrixos formában átírva kapjuk, hogy
u N1 u = v = 0 w 0
0
0
… N ncs
N1
0
…
0
N1 … N
0
0
N ncs
0
0
e
u1 v 1 0 w1 0 = Nq N ncs uncs vncs w ncs
(1.188)
qe
1.7.2. Alakváltozások és feszültségek A szimmetrikus alakváltozási- és feszültségi tenzorok független elemeiből alakváltozási- és feszültségi oszlopvektorokat képezünk. Az alakváltozási vektor elemeit az elmozdulási mező deriváltjaiból állítjuk elő: ∂u ∂ ε x ∂x ∂x ε ∂∂yv 0 y ∂w 0 ε ∂z A ⇒ ε = z = ∂u ∂v = ∂ + γ xy ∂y ∂x ∂y γ yz ∂v + ∂w 0 ∂z ∂y ∂w ∂u ∂ γ xz ∂x + ∂z ∂x
0 ∂ ∂y
0 ∂ ∂x ∂ ∂z
0
0 0 ∂ u ∂z v = ∂ u = ∂ Nq = Bq 0 B w ∂ ∂y ∂ ∂z
∂
60
(1.189)
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
ahol a B alakfüggvények deriváltjait, tartalmazó mátrixot szokás alakváltozási-elmozdulási mátrixnak is nevezni. A számítások elvégzéséhez be kell vezetni a globális x, y, z koordinátarendszerben értelmezett deriváló operátort
∂ ∂∂ξx ∂x ∂ G = ∂∂y = ∂∂ξy ∂ ∂ξ ∂z ∂z
∂η ∂x ∂η ∂y ∂η ∂z
∂ ∂ξ ∂∂η ∂ζ ∂ ∂z ∂ζ ∂ζ ∂x ∂ζ ∂y
(1.190)
∂L
J −1
amelyet a lokális ξ ,η , ζ koordináták szerinti deriváltakkal állíthatunk elő. A differenciáló operátorok közötti kapcsolat megfordítva is érvényes, azaz ∂ G = J −1 ∂ L
∂ L = J ∂G ,
(1.191) ncs
ahol J Jacobi mátrix a következő módon számítható, mivel x = ∑ N i xi i =1
y1 z1 ∂∂Nξ1 ∂∂Nξ2 … ∂N∂ξncs x1 y2 z2 ∂N 2 ∂y ∂N ncs x2 ∂z = ∂N1 (1.192) … ∂η ∂η ∂η ∂η ∂η ∂y ∂N 2 ∂z ∂N1 … ∂N∂ζncs xncs yncs zncs ∂ζ ∂ζ ∂ζ ∂ζ A feszültségmezőt az alakváltozások ismeretében lineárisan rugalmas anyagot feltételezve a Hooke-féle anyag-egyenlettel adjuk meg. A független feszültségi mennyiségeket tartalmazó feszültségi vektor
∂x ∂ξ J = ∂∂ηx ∂x ∂ζ
∂y ∂ξ
∂z ∂ξ
T ⇒ σT = σ x σ y σ z τ xy τ yz τ xz
(1.193)
az alábbi mátrixegyenlettel írható fel σ = Dε = DBq
(1.194)
ahol D az anyagjellemzők mátrixa
c1 c 2 c D= 2 0 0 0
c2 c1 c2 0 0 0
c2 c2 c1 0 0 0
0 0 0 c3 0 0
0 0 0 0 c3 0
0 0 0 0 0 c3
(1.195)
az egyes mátrix elemek a következő képletekkel adottak. c1 = E
1 −ν (1 +ν )(1 − 2ν )
c2 = E
ν
(1 +ν )(1 − 2ν )
61
c3 =
E =G 2 (1 + ν )
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.7.3. A Merevségi mátrix A csomóponti elmozdulásokkal kifejezett alakváltozások és feszültségek ismeretében az elem merevségi mátrixa az alakváltozási energia sűrűség térfogaton vett integráljából képezhető. Az izoparametrikus elem egyik előnye éppen abban mutatkozik, meg, hogy ezt az integrálást a ξ ,η , ζ természetes koordinátarendszerben nagyon könnyen elvégezhetjük. A K e merevségi mátrixot végül is numerikus integrálással a Gauss féle kvadratúra felhasználásával állítjuk elő. K e = ∫ B eT DB e Ve
1 1 1
=
dV det J d ξ dη d ζ
∫∫∫B
D Β e det J d ξ dη d ζ = F (ξi , η j , ζ k )
eT
−1−1−1
= … ∑∑∑wi w j wk F (ξi , η j , ζ k ) = Κ ij NG NG NG
e
(1.196)
i, j = 1, … , ncs
i =1 j =1 k =1
ahol egy-egy csomóponthoz, illetve csomópontpárhoz tartozó mátrixblokk külön is megmutatható:
K ij = ∫∫∫BieT DB ej det J d ξ dη d ζ
(1.197)
( 3,3)
1.7.4. Tehervektorok Az elem terhelése lehet térfogaton- és felületen megoszló erő. Ezen elemtípus esetén a nyomatéki terhelés nem értelmezett. Csak két különböző csomóponton működő erőpárral fejthetünk ki nyomatéki terhelést.
1.7.4.1.
Térfogati terhelés
A térfogaton megoszló terhelés munkájának integráljából származtatható a tehervektor: f ρe = ∫ eN eT ρ k dV = ∫ V
1
1
∫ ∫
1
−1 −1 −1
N eT ρ k J d ξ dη d ζ
(1.198)
F
− gyorsulás ahol térfogati terhelés intenzitása lehet ρ k = − forgás − súly
. Az integrálást a térfogaton
ismét a Gauss-féle kvadratúra segítségével hajtjuk végre. A tehervektorban az elem csomópontjainak a számával megegyező számú blokk található. Vagyis a megoszló erőrendszerből a csomópontokba redukált tehervektort kaptunk. f ρe = [ fie ] i = 1, … , ncs
(1.199)
(3, 1)
1.7.4.2.
Felületi terhelés
A másik gyakori terhelési típus a felületen megoszló terhelés. Az ebből származó tehervektor előállítása kissé komplikáltabb, mert a felületen képeznünk kell a területvektort is:
62
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1.29. ábra. Felületi terhelés az elem oldalán Az 1.29. ábrán az elem egyik oldalát a felületi normálissal ellentétes irányban p nagyságú megoszló erőrendszer terheli, a helyi rendszerben kifejezve: p = − pn . Az 1.30. ábrán látható felületelemvektor a ξ ,η koordinátairányú deriváltak segítségével írható fel:
a1 =
∂r ∂ξ
a2 =
∂r ∂η
(1.200)
1.30. ábra. Felületelem vektor Az (1.200) és az 1.30. ábra alapján, a felületelemen működő erő az alábbi módon írható
p dA = − pn dA = − p dA = − p ( a1 × a2 ) d ξ dη = − p a 3 d ξ dη
(1.201)
melyből a csomóponti redukált tehervektor a kijelölt integrál végrehajtásával áll elő: 1 1
f pe = − ∫ ∫ N eT p a3 d ξ dη
(1.202)
−1−1
Végül az elem térfogati és felületi terheléséből származó terhelési vektora összegzéssel megkapható:
f e = f pe + f ρe
(1.203)
Az elemen a teljes potenciális energia a szokásos alakban írható fel 1 Π ep q eT K e q e − q eT f e 2
(1.204)
63
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
2. A MEGOLDÁS ÉS ANNAK HIBÁJA Egy valóságos szerkezet mechanikai modelljének végeselemes közelítő megoldásakor kérdésként merül fel, hogy milyen a közelítés hibája és megoszlása, a megoldás hogyan pontosítható? A hiba származhat a modell hibájából, az anyag és terhelési adatok pontatlanságából, az alkalmazott közelítés jellegéből és az egyenletrendszer megoldása során fellépő kerekítési hibákból is. A végeselemes eljárás numerikus szempontból lineáris algebrai egyenletrendszer megoldására vezeti vissza az eredeti parciális differenciálegyenlet-rendszereket is tartalmazó peremérték feladatot. Az egyenletrendszer szerkezete lényegesen befolyásolja a számítógépi tárolási- és megoldási szükségletét. A végeselemes közelítés alapvetően polinomok alkalmazásával történik. Ha a közelítő függvény lineáris vagy kvadratikus, akkor a végeselemet h-verziójúnak, ha magasabb rendű, akkor p-verziójúnak nevezzük. Az első esetben a megoldás pontosítása az elem jellemző (h) méretének csökkentésével, míg a második esetben az elmozdulást közelítő polinomok (p) rendjének növelésével érhető el. A két módszer ötvözését hp-verzós eljárásnak hívjuk [2].
2.1. Az együttható Mátrix szerkezete Az egyes végeselemekre vonatkozóan, a korábbi fejezetekben bemutattuk az elemek merevségi mátrixának és terhelési vektorának származtatását. Ezen végeselemes mennyiségekkel is felírhatjuk az elemen a teljes potenciális energiát: 1 (2.1) Π ep = q eT K eq e − q eT f e , 2 amely részben alakváltozási energiából és részben a külső erőrendszer munkájából épül fel. Az elem merevségi mátrixa az alakváltozási energiából, a tehervektor a külső erőrendszer munkájából származtatható. Egy végeselemekre felosztott szerkezet teljes potenciális energiája az egyes elemek potenciális energiájának összegéből épül fel. A végeselem-módszer lehetővé teszi nem csak azonos típusú, hanem geometriailag teljesen különböző alakú elemek csatolását is. Az lényeges, hogy az elemek illesztésekor ne sértsük meg az elmozdulás folytonosságát, idegen szóval kompatibilitását. Vagyis a csatlakozó elemek peremein azonos számú és szabadságfokú csomópont csatlakozzon. Ez azt jelenti, hogy a közös csomópontokba befutó elemek elmozdulásai kell, hogy megegyezzenek. A szerkezet csomóponti elmozdulás vektora ncs számú blokkot tartalmaz qT = q1T
qT2
… qTncs .
(2.2)
Hasonlóan a szerkezet f csomóponti terhelési (redukált) vektora is, amely az egyes elemeken értelmezett külsőerők munkájának összegzésével, illetve az elemekhez nem rendelt koncentrált erőhatásokból kapható meg: ncs
∑q
eT
f e + W k = qT f ,
(2.3)
e =1
65
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
ahol f = [fi ] i = 1, … , ncs fi = ∑fie + fik . e∈i
A szerkezet merevségi mátrixa: K a z egyes elemek alakváltozási energiájának összegzéséből származtatható a megfelelő csomópontokhoz-, illetve csomópont párokhoz tartozó merevségi mátrix blokkok összegeként: Ne
∑q
eT
K eq e = qT Kq ,
(2.4)
e =1
ahol K = [K ij ] i, j = 1, … , ncs K ij =
∑K
e ij
.
e∈i , j
Végül a szerkezet teljes potenciális energiája hasonló alakban írható, mint az elem (12.1) teljes potenciális energiája: e 1 1 (2.5) Π p = ∑ q eT Kq e − q eT f e − W k = qT Kq − qT f . 2 e =1 2 A teljes potenciális energia minimuma elvéből nyilvánvalóan következik, hogy az első variációja nulla kell, hogy legyen:
N
δΠ p = δ qT
∂Π p ∂q
=0
(2.6)
Majd a (2.5) képleten végrehajtva a varráció műveletét kapjuk, hogy
δ qT ( Kq − f ) = 0 .
(2.7)
A kinematikai peremfeltételek figyelembe vétele után a szabad csomópontok elmozdulásának a variációja tetszőleges, ebből következik, hogy a (2.7) egyenlet zárójeles kifejezése zérus vektorral egyenlő. Azaz a csomóponti elmozdulások vektorára egy lineáris algebrai egyenletrendszert kapunk: Kq = f .
(2.8)
Ismét megjegyezzük, hogy az együttható mátrixot és tehervektort itt nem részletezett módon módosítottuk a kinematikai peremfeltételek figyelembe vételekor K módosított merevségi mátrix f
módosított terhelési vektor
A csomóponti sorszámok kiosztását a végeselemes szoftverek a hálógenerálás során, esetleges módon, többnyire kedvezőtlen számozással végzik. Az egyenletrendszer megoldása előtt, azonban az átsorszámozó algoritmusok segítségével lényegesen csökkentik a sávszélességet. Ezek az algoritmusok többnyire gráfelméleti alapon találják meg a közel optimális sorszámkiosztást. Kedvezőbb a helyzet, ha a szerkezet alszerkezetekre bontható, mert akkor az egyes alszerkezeteken lényegesen kisebb sávszélesség érhető el külön-külön, mint ha azt a teljes szerkezetre végeznénk el. Kedvező akkor a sávszélesség, ha az egyes elemeken az előforduló sorszámok különbsége a legkisebb. Ezért a sorszámozást úgy célszerű elvégezni, hogy a szerkezet egyik perempontjából kiindulva az oda befutó elem összes további pontját is beszámozzuk, majd ezt követően a korábbi legkisebbtől folytatva ismételjük ezt az eljárást. A 2.1. ábra az előbbiekre mutat egy példát, ez egy optimális sorszámkiosztás.
66
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
2.1. ábra. A K mátrix struktúrája az optimális sorszámozás esetén Az egyenletrendszer megoldásának két nagy csoportja terjedt el a kereskedelmi végeselemes szoftverekben: 1) Direkt eljárások: Gauss elimináció valamilyen alkalmazása, ez a technika kb. ≤ 100 ezer ismeretlenig megbízható, illetve elfogadható sebességű. 2) Iterációs megoldások: mátrixok szorzásával jutunk egyre közelebb a megoldáshoz.
Ez utóbbi eljárás előnye, hogy elegendő csak a zérustól különböző mátrix elemeket tárolni, így valójában érzéketlen a csomóponti sorszám kiosztására. Azonban a direkt megoldóknál az un. feltöltődés miatt a sávon belül a zérus mátrix elemeket is szükséges tárolni. A 2.2. ábra a legelterjedtebb tárolási struktúrákat szemlélteti. A frontális technikánál egyidőben nem készül el a teljes szerkezeti mátrix, csupán egy-egy blokkját találjuk meg a memóriában.
2.2. ábra. Mátrix tárolási formák lehetnek: félsáv blokkolva, aktív oszlop blokkolva, frontális technika (tárolás elemi szinten) A direkt megoldó eljárások előnybe részesülnek az olyan feladatok esetén, amikor a mátrix együtthatóinak nagyságrendje igen nagy különbséget mutat. Az un. locking, azaz rossz mátix kondicionáltság tipikusan ezt a helyzetet eredményezi. De ilyen feladatok lehetnek, pl. az érintkezési feladatok, héjfeladatok és az összenyomhatatlan tulajdonsággal rendelkező gumit tartalmazó problémák is.
2.2. Hibaanalízis A továbbiakban a végeselemes közelítés hibáját elemezzük. A korábbi fejezetekben bemutatott végeselemes elmozdulási módszert a teljes potenciális energia funkcionálra alapoztuk. A megoldás közelítéséből származó hibáját is energia értelemben határozzuk meg. Mint ismeretes a rugalmas feladatok esetén, adott elmozdulási mező alapján az alakváltozásiés feszültségi mező egyértelműen meghatározható: u elmozdulás → ε = ∂ u alakváltozás, σ = D ε feszültség .
67
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
Az alakváltozási energia és ezen keresztül az elmozdulás normája az előbbi mennyiségekkel az alábbi módon fejezhető ki 1
2 1 1 T || u ||= U = ∫εT σ dV = ∫ ( ∂ u ) D ( ∂ u ) dV . 2V 2 V 1 2
(2.9)
Legyen a továbbiakban u = u ex az egzakt megoldás, uVEM a közelítő véges elemes megoldás, akkor az utóbbi hibája, a kettő különbsége e = uVEM − u ex .
(2.10)
Az e hiba pontszerű értelmezése a gyakorlatban ritkán határozható meg, de a (12.9) energia norma alapján értelmezve 1
1 2 || e ||= ∫ ( ∂ e ) D ( ∂ e ) dV 2 V
(2.11)
már léteznek matematikai megalapozottságú becslések. Energia értelemben konvergens a megoldás, ha a hiba normája az ismeretlenek számának növelésével tart a nullához: lim || e ||= 0 ,
(2.12)
N →∞
ahol N az ismeretlenek száma, térbeli izoparametrikus elemeknél ( 3 ⋅ NCS ) . Pontonkénti konvergenciáról beszélünk, ha minden pontban teljesül az alábbi határérték lim e = 0 .
(2.13)
N →∞
Az ismeretlenek számának növelése alapvetően kétféle módon is megvalósítható. Az egyik esetben az elemháló sűrítésével csökkentjük az elemek jellemző h méretét, változatlan lineáris vagy kvadratikus közelítő mező alkalmazása mellett. A másik esetben változatlan felosztás, azaz változatlan elemméret mellett, a közelítő polinomok p rendjét növeljük. A két eljárás ötvözete egyaránt magába foglalja a h méret csökkentését p polinom rendjének növelését. Az alkalmazott módszereket tekintve beszélhetünk h -verziójú-, p -verziójú- és hp -verziójú végelemes eljárásokról. A p -verziónál alkalmazott függvénytér felépíthető Lagrange-féle és Legendre-féle polinomokkal is. Az utóbbi alkalmazása azért előnyösebb, mert az approximációs tér hierarchikusan egymásba ágyazott. Ekkor a magasabb rendű p polinom alapján előállított merevségi mátrix leválaszthatóan tartalmazza az alacsonyabb rendű közelítések mátrixait is, a lineárissal bezárólag. A Legendre-féle polinomok és deriváltjaik ortogonális tulajdonsággal rendelkeznek, ezért a megoldandó egyenletrendszer kondicója is kedvezőbb mint a Lagrangeféle közelítés esetén. A peremérték feladatokat az irodalomban három csoportba szokás sorolni: A. típusról beszélünk, ha megoldás elegendően sima, vagyis a vizsgált tartomány szingularitásokat nem tartalmaz, azaz analitikus jellegű. B. típus esetén szingularitásokat tartalmaz a feladat, de ez a szinguláris hely az elem csomópontjába esik C. típusnál a szingularitások tetszőlegesen helyezkednek el, azaz nem esnek csomópontokba.
68
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A 2.3 ábra példákat mutat be a szinguláris helyekre. Ezek lehetnek, pl. éles saroknál, koncentrált erő támadási helyén, kompozit anyagok peremein. A szinguláris pontok a gyakorlatban a tönkremenetel kiindulási helyeiként igen veszélyesek lehetnek.
2.3. ábra. Példák a szinguláris pontokra A szinguláris pont környezetében az elmozdulási mező lefutását a szingularitás jellege határozza meg. A szingularitás lehet erős és gyenge: ∞
u = ∑ Φi (ϕ ) r λi
r < r0
i =1
ahol
2.4. ábra. A szinguláris pont környezete
r0 az elhalási hossz
< 1 szigorú . ha min λi > 1 nem szigorú
A szigorú szingularitás esetén a feszültség tart a végtelenbe, ha nem szigorú, akkor véges értékű lesz. Az éles sarok geometriája, azaz nyílásszöge döntő befolyással van a szingularitás jellegére. Ha nyílásszög kisebb, mint 120o , akkor a csúcspont veszélyes feszültséggyűjtő hely lehet. A különböző típusú feladatokra vonatkozóan az irodalomban a következő hibabecslő formulák találhatók a h-, p-, hp- verziójú közelítésekre. 2.1. táblázat: Hibabecslő összefüggések A
B p −2
h
|| e ||=≤ k ⋅ N
p
|| e ||≤ exp ( −γ N δ ) ⋅ k
δ> hp
C
|| e ||≤ k ⋅ N
− 12
min ( p , λ )
|| e ||≤ k ⋅ N − λ
|| e ||≤ k ⋅ N
− 12 min ( p , λ )
|| e ||≤ k ⋅ N
− 12 λ
1 2
|| e ||≤ exp ( −γ N δ ) ⋅ k
δ>
1 3
A 2.1. táblázatban szereplő N az ismeretlenek számát, p a közelítő polinom fokszámát, λ a szingularitás mértékét jelenti, k konstans érték.. A 2.1. táblázatból jól látható, hogy a p-verzós közelítés gyorsabb konvergenciával rendelkezik, mint a hagyományos h-verziós. A hp-verziós eljárás exponenciálisan gyors konvergenciájú még B-típusú feladatok, azaz szingularitások csomóponti elhelyezkedése esetén is. Ekkor a felosztást a szingularitás közelében geometria sor szerint szükséges sűríteni. 69
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A konvergencia sebességeket hasonlítja össze a 2.5. ábra is. Az ábra jól mutatja a p- és hpverzió előnyét, mert ugyanolyan hibahatár eléréséhez lényegesen kisebb az ismeretlenek száma a hagyományos h-verziós számításhoz képest [2].
2.5. ábra. A h-, p- és hp-verziós számítások konvergenciája
70
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
3. MODELLEZÉSI KÉRDÉSEK A szerkezetek számításánál számos olyan probléma kerül előtérbe, hogy hogyan lehet a megfogásokból származó peremfeltételeket, az elemek közötti túlfedésből származó hatásokat, a teljes szerkezetnél az ismétlődésből származó periodicitást, a ferde hatásvonalú megtámasztásokat stb. kényelmesen kézben tartani, a számítási időt csökkenteni. A felsorolt problémákat érdemes elem szintjén kezelni. A számítógépes tervezés során általában a teljes szerkezetet nem lehet minden részletre kiterjedően a képernyőn megjeleníteni, ill. gyakran kész szerkezeti elemek, részegységek kerülnek beépítésre, amit a modellezésnél fel kell tudni használni. Ez szintén újabb megfontolást igényel a végeselemes számítás megszervezésére, a modellünk felépítésére [1].
3.1. Alszerkezettechnika Tételezzük fel, hogy a több szerkezeti egységből álló szerkezet egyes részeinek végeselemes felosztását már előállítottuk. A merevségi mátrixot és a csomóponti redukált terhelési vektort kiszámoltuk. Az egyes részeket bizonyos, a tervező által megálmodott felületek mentén össze kell illeszteni [1]. A 3.1. ábrán lévő egyszerű felépítésű szerkezetet két alszerkezetre ( i = 1, 2 ) bontjuk fel. A mechanikai probléma végeselem-módszerrel történő megoldásánál megkövetelt pontosság elérésére megfelelő számú és fokszámú elemeket használunk fel. Az alszerkezetek az Aci felületük mentén csatlakoznak egymáshoz. Az itt található csomóponti elmozdulások vektora q c i , míg a megmaradóké, röviden a belső pontoké qbi . Az alszerkezet összes csomópontjának elmozdulásvektora qi . A teljes potenciális energia minimuma elv szerint fennáll, hogy
=
= A1c fő
1
A2c 2
3.1. ábra. Két alszerkezetre felbontott szerkezet ∂Π ip
K bb = − − = K q f r ∂qi K cb i
i
i
i
K bc K cc
i
i
i
i
q b fb 0 − − = 0, q c f c r
(3.1)
ahol K i az i-edik alszerkezet merevségi mátrixa, f i csomóponti redukált terhelési vektor, r i a csatlakozásnál fellépő belső erőből (hatás-ellenhatás törvénye szerint keletkező) csomóponti vektor. A kapott mátrixegyenlet első blokksora
71
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
∂Π ip ∂q b
amiből
i
= K bbi q bi + K bci q c i − fbi = 0 ,
(3.2)
qb i = ( K bbi ) fbi − ( K bbi ) K bc i qc i . −1
−1
(3.3)
Itt feltételeztük, hogy a csatlakozó pontok száma elegendő ahhoz, hogy a vizsgált alszerkezet i merevségi merevtestszerű elmozdulása le legyen kötve, azaz a belső pontokra vonatkozó Kbb mátrix inverze létezzen. A kapott belső elmozdulások vektorát behelyettesítve a ∂Π ip ∂q ci
= K cb i q bi + K cci q c i − fc i − r i = 0
(3.4)
egyenletbe, nyerjük, hogy
{
K cc i − K cb i ( K bbi ) K bc i −1
} qc i = fc i − K cbi ( K bb i )
−1
fb i + r i ,
(3.5)
ami tömörebben
K ired q c i = f red i + r i
i = 1, 2
(3.6)
i alakban írható fel. Itt K i red , f red az i-edik alszerkezet csatlakozó csomópontokra redukált merevségi mátrixa és redukált csomóponti terhelési vektora. Mivel a hatás-ellenhatás értelmében r1 = −r 2 , továbbá a csatlakozási csomópontokban
q c1 = q c 2 = q c
(3.7)
A csatlakozó pontok egyensúlyát kifejező egyenletek összegzésével a következő végső egyenlethez jutunk a csatlakozó csomóponttokbeli elmozdulás meghatározására:
i i (3.8) ∑ K red qc = ∑ fred i i vagyis megkaptuk a főszerkezet egyensúlyi egyenletét. Az egyenlet megoldásából nyert, a csatlakozási felületen lévő csomópontokra vonatkozó q c = q c1 = q c 2 elmozdulási vektor ismeretében (3.3) alapján az i-dik alszerkezet belső csomópontjainak elmozdulásvektorta ismert lesz. Ezekután a teljes qi vektor ismeretében az i-dik alszerkezet elmozdulás és feszültségállapota számolhatóvá válik, aminek elemzése révén eldönthető annak jósága, avagy további szerkezeti módosítással kell a tervezési feladatot pontosítani.
A módszer előnyei: egyszerűbb az adatelőkészítés, a tipizált alkatrészek, szerkezeti egységek merevségi mátrixait, terhelési vektorait előre ki lehet számolni, azokat el lehet raktározni és újbóli számításnál a teljes rendszerbe könnyen be lehet illeszteni. Az alrészek számításánál a többprocesszorú, párhuzamos számítás technikáját is fel lehet használni jelentős időt megtakarítva. A sávszélesség minimalizálása egyszerűbb alszerkezeti szinten, mint a teljes rendszer vonatkozásában. A számítási eredmények birtokában azon részeken, ahol nem kell i változtatást végrehajtani, az elraktározott K i red , f red mennyiségek újból felhasználhatók, az újraszámítást csak azon részeken kell végrehajtani, ahol a geometriában, anyagban, esetleg a terhelésben álltak be változások. Ezzel gyorsítani lehet a végső tervek elérését. Az alszerkezetekkel kezelt rendszereknél az I/O műveletek száma csökken. Gyakran a
72
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
számítógépi memória korlátja miatt is előnyös használata, mivel nem kell egyszerre a teljes egyenletrendszert tárolni. A gyakorlatban, nagybonyolultságú szerkezeteknél többszintű alszerkezeti struktúra felépítése is javasolt.
3.2. Adott elmozdulások figyelembevétele Az elmozdulásmódszernél, a teljes potenciális energia minimuma elv használatakor a kinematikailag lehetséges elmozdulásmezőnek a kinematikai peremfeltételt ki kell elégítenie. Feltételezéseink értelmében az Au felületre kifutó végeselemek csomópontjainak elmozdulásával a teljes felületen megadott elmozdulás függvényt leírjuk. Így az elem szintjén nagyon egyszerű az adott elmozdulás figyelembevétele. Legyen a teljes potenciális energia Π e p = Π e p ( qe ) =
1 q s e,T q u e,T ( 2
K ess e K us
f se K esu q es 2 − e e e K uu q u fu
)
(3.9)
ahol que a q e csomóponti elmozdulásvektor azon része, amelynél az elmozdulások adottak, míg a q es -el jelöljük, a szabad, ismeretlen elmozdulásokat magában foglaló vektort. A kijelölt műveletek elvégzésével, a minimalizálás szempontjából állandó tagokat elhanyagolva, a minimalizálandó energia Π p = ∑ Π e p ( q es ) = ∑ e
e
1 e,T q s ( K ess q es − 2 ( f se − K esu q u e ) 2
),
(3.10)
vagyis az adott elmozdulás egy kinematikai terhelést jelent, ami arányos az adott elmozdulással −K esu que . Gyakran a kinematikai hatásokat külön terhelésként kezelik, rugalmas szerkezetről lévén szó a szuperpozíció elvének felhasználásával jutunk a teljes - a kölcsönhatásból származó erőhatásokból is származó - terhelés figyelembevételéhez. Ebben az esetben az egyenletrendszert megoldó eljárásnak un. több jobboldalas számításra is alkalmasnak kell lennie. Az alapterhelések megoldásainak lineáris kombinációjával juthatunk el a kívánt terhelések összegzett hatásának az elemzésére. E technikával gépidő takarítható meg. Ugyanis, a szerkezet méretétől függően az alapterhelésekhez tartozó elmozdulások kiszámítása igényel valójában jelentős időt. Azok lineáris kombinációja már gyorsan elvégezhető, melyek variációi a tervezési folyamattól függnek, azok bármikor – az elraktározott futási eredmények birtokában megismételhetők, újjakkal tetszés szerint kiegészíthetők.
3.3. Adott elmozdulásmezőben fennálló szakadás, kezdeti hézag figyelembevétele A gépészmérnöki gyakorlatban gyakran az alkatrészek között túlfedéssel valósítunk meg kötést. Ebben az esetben a kapcsolatot oly módon is tudjuk modellezni, hogy kétoldalú kapcsolatot tételezünk fel, vagyis az alakváltozás után a két test párbaállított pontjai azonos helyet foglalnak el. Ez a modellezésben egy egyszerűsítés, mert az elemek közötti normális érintkezési feszültség előjelére és a súrlódási feltételekre nem vagyunk tekintettel. A geometriai illesztési feltétel kielégítésével a számítás után lehetőség van a feszültségi feltételek ellenőrzésére. Amennyiben az érintkezési tartományon nincs semmiféle adhézió,
73
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
akkor a normál feszültség csak nyomó lehet. Száraz súrlódásnál a COULOMB-féle egyenlőtlenségi feltételnek is fenn kell fennállnia. Nézzük az egyszerűsített modellünket. Tételezzük fel, hogy az A és az F párbaállított csomópontok általánosított csomóponti elmozdulása között a q F = q A + h FA
(3.11)
kapcsolat áll fenn, ahol h FA a szakadásból (kezdeti hézagból) származó vektor. Az F pontot főcsomópontnak, az A pontot alcsomópontnak nevezzük. Az összefüggés értelmében alakváltozás után a két pont a tér egy közös P pontjába kerül (3.2. ábra). Itt is feltételezzük, hogy a csatlakozó Ac tartományon a hézag változását, az elmozdulásmezőt a csomóponti értékek egyértelműen leírják. n z qA
A
P
x hFA qF F
3.2. ábra. Kétoldali kapcsolat x és z irányában Az A csomópontot magába foglaló e jelű elem teljes potenciális energiája az A csomóponthoz q eA és az elem megmaradó pontjaihoz tartozó q em csomóponti elmozdulás vektorokkal írható fel: Π e p = Π e p ( q eA , q em ) =
1 e,T e,T q A q m ( K e 2
q eA f Ae − 2 e e q m f m
),
(3.12)
illetve Πe p =
q − h f e 1 T q F − hTFA , q em,T ( K e F e FA − 2 Ae 2 q m fm
).
(3.13)
A merevségi mátrixot az elmozdulásvektor felbontása szerint négy részre tagoljuk e
K AA K Am Ke = . K mA K mm Ennek felhasználásával nyerjük, hogy Π
e p
=Π
e p
(q
F
,q
e m
) − q
(3.14)
e
T F
,q
e ,T m
K AA h FA , K mA
74
(3.15)
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
vagyis az A pontbeli ismeretlenek F pontba való áthelyezésével az e-edik elem merevségi mátrixa nem módosul, a csomóponti terhelésé azonban igen: f Ae K AA = e+ h FA f m K mA e
f
e mód
(3.16)
Gyakorlatban, számos esetben az al- és főcsomópontok között nem az összes koordináták között van alárendeltség, továbbá az elemnek nem csak egy alcsomópontja van, hanem több. Formálisan, a kapott eredmények ekkor is érvényben vannak: A módosított redukált csomóponti terhelési vektornál a merevségi mátrix megfelelő oszlopait kell megszorozni a h FA vektorral. A kapott eredmények vizsgálatával lehet eldönteni, hogy a kétoldalú kapcsolatú érintkezési feltételek valóban fennállnak-e avagy nem. Ha a testek közötti adhéziótól eltekintünk, akkor az érintkezési felületen keletkező feszültségnek normális összetevője csak nyomó feszültség lehet. Ha a kapott feszültségkép ettől lényegesen eltér, akkor a feladatot az egyoldalú kapcsolatú érintkezési feltételek mellett kell megoldani. Ezzel a kérdés komplexummal a „Dinamikai és nemlineáris feladatok modellezése” c. szakmérnöki jegyet 3. fejezete foglalkozik. Hengeres testek esetén a radiális irányú alárendeléssel lehet a túlfedés hatását figyelembe venni, míg a tengelyirányban az alárendelés hiánya a súrlódásmentességet, alárendelés (elmozdulások azonossága) pedig a végtelen értékű súrlódási tényező felvételét jelenti. A valóságban a súrlódási tényező véges értéke miatt, a két feltételezéssel kapott eredmény között van az „igazság”.
3.4. Excentrikus csatlakozás A valóságos háromdimenziós szerkezet szilárdságtani feladatának megoldásakor gyakran használunk egyváltozós ill. kétváltozós redukált modelleket. Jól ismerjük, hogy egyváltozós feladatként kezeljük a rudakat, kétváltozósként a lemezeket, héjakat. Rúdszerkezeteknél a rudak középvonalai számos esetben kitérők, nem mindig metszik egymást, azt mondhatjuk excentrikusan illeszkednek egymáshoz. A lényeg megértése céljából tekintsük a 3.3. ábrán vázolt esetet. A csomóponti általánosított e
elmozdulásvektor q Ae,T = U , W , ϕ y . Az A és F pont közötti kettős vonallal rajzolt szakaszt A merevnek testnek tekintjük, azaz a merevtestszerű szögelfordulás az A és F pontban azonos. Bevezetve az x = x A − xF , z = z A − zF az A és az F pontok közötti „távolságokat” kapjuk, hogy e
q Ae
U z U 1 0 ˆ EX q . = W = 0 1 − x W = T F ϕy 0 0 1 ϕ y F A
(3.17)
Az e jelű elemen az A csomóponton kívüli helyeken qem a csomóponti elmozdulásvektor. Így e
ˆ EX 0 q F q A T q F EX (3.18) q = = e =T e E q m q m 0 q m és az elem teljes potenciális energiája a főponti és a megmaradó pontokbeli csomóponti elmozdulásvektoron keresztül e
75
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
Π
e p
T 1 = qTF q em,T ( T EX ) ( K e T EX 2
q F f Ae e −2 e q m fm
),
(3.19)
z X~
WF
A Z~
UF
F
ϕ y, F x
3.3. ábra. Excentrikus csatlakozás síkbeli esetben
amiből a transzformált merevségi mátrix K e = ( T EX ) K e T EX T
(3.20)
és a redukált terhelési vektor f e = ( T EX ) f e . T
(3.21)
Térbeli esetben az y = yA − yF elfordulásból 0 ∆u = −z y
z 0 − x
mennyiség további bevezetésével a merevtestszerű
− y ϕx x ϕy = Ωϕ 0 ϕz
(3.22)
elmozdulás keletkezik. Iymódon e
qe A
U V E W = = 0 ϕx ϕy ϕz A
U V Ω W =T ˆ EX q . F E ϕx ϕy ϕz F
(3.23)
A teljes transzformáció a (3.18) alattiak értelemszerű alkalmazásával történik. Térbeli esetekre vonatkozóan további információ az [1] irodalomban olvasható.
3.5. Ferdehatásvonalú támasz figyelembevétele A szerkezetek egy részénél a megtámasztási korlátok nem párhuzamosak a választott koordinátarendszer tengelyeivel. Ilyen megtámasztások, korlátok a ferde hatásvonalú görgős támaszok, különféle csuszkák. Látni fogjuk ezeket az eseteket kényelmes a ferde megtámasztáshoz kötött helyi koordinátarendszerben tárgyalni. A vizsgálatainkat síkbeli esettel kezdjük [1].
76
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
3.5.1. Síkbeli eset A 3.4. ábra alapján a síkbeli ferde hatásvonalú görgős támasz lokális rendszerében értelmezett és az x–z síkban értelmezett i csomópontbeli elmozdulások között az alábbi összefüggések állnak fenn: U cos β q Gi = = W sin β
− sin β u g cos β un
(3.24)
Mivel a görgő elmozdulásának irányára merőlegesen az u n = 0 , úgy az előbbi egyenlet U cos β qGi = = (3.25) u g = TGi qGi W sin β alakot ölti. Az elem összes csomópontjához tartozó görgős támaszok elmozdulásait egybegyűjtve, a görgős megtámasztású pontok qG globális elmozdulásvektora a lokális qG elmozdulásvektorral formálisan kifejezhető q G = TG qG . (3.26) Az elem csomóponti elmozdulási vektorát két részre felbontva kapjuk, hogy q e,T = q Ge,T q em az elem teljes potenciális energiája
Π e p = Π e p ( q Ge , q em ) =
(3.27)
K e 1 qG e,T q m e,T ( GG e 2 K mG
e q Ge fGe K Gm 2 − e K emm q em fm
),
(3.28)
un z
W
ug
n t
β
U
x
3.4. ábra. Síkbeli ferdehatásvonalú görgős támasz
ami a transzformációs összefüggés felhasználásával Π e p = Π e p ( qGe , q em ) = e e (3.29) TGT K GG qGe TGT fGe TG TGT K Gm − 2 , ) e e K emm q em K mG TG fm vagyis a görgős megtámasztáshoz kötött helyi koordináta-rendszerbeli elmozdulásra áttérve, az elem merevségi mátrixának és redukált terhelési vektorának az áttranszformálására van szükség. Látjuk, hogy egy görgőnél a helyi rendszerben az ismeretlenek száma eggyel csökkent, mivel csak az u g szerepel.
1 = qG e,T q m e,T ( 2
77
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
3.5.2. Térbeli eset Ebben az esetben a megtámasztás olyan, hogy a csomópont egy adott, ferdén elhelyezkedő síkban tud elmozdulni, pl. gömbön keresztül csatlakozik a ferde síkkal. A szóban forgó síkban az elmozdulást, tetszőlegesen felvett, két egymásra merőleges e x és e y irányú elmozduláson keresztül fogjuk szemlélni, jelölje ezeket u és v . A ferde sík normálisa e z jobb sodrású koordinátarendszert alkot az e x és e y irányokkal. Ortogonális transzformáció révén a globális x,y,z koordinátarendszerbeli elmozdulás U e x ⋅ e x e x ⋅ e y e x ⋅ e z u e x ⋅ e x e x ⋅ e y u V = e ⋅ e e ⋅ e e ⋅ e v = e ⋅ e e ⋅ e y x y y y z y x y y v W e z ⋅ e x e z ⋅ e y e z ⋅ e z w e z ⋅ e x e z ⋅ e y és így az i jelű csomóponti, csak elmozdulásokat tartalmazó elmozdulásvektor qGi = TGi qGi .
(3.30)
(3.31) ~e z
z ~e y
~e x x y
3.5. ábra Térbeli ferdesík által kijelölt megtámasztás
3.5.3. Csuszka Egy csuszka esetén a csomópont csak egy kitüntetett irányban tud elmozdulni. Legyen ez az e x . Ekkor a csuszka irányú elmozduláson keresztül a globális rendszerbeli elmozdulások a következők módon számolhatók
U e x ⋅ e x qGi = V = e y ⋅ e x W i e z ⋅ e x
,
[u ] = TGi qGi
~e x u~
3.6. ábra. Csuszka
78
(3.32)
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
3.6. Periódikus szerkezet A gépészet, építészet szerkezetei gyakran rendelkeznek szimmetria tulajdonságokkal, vagy ismétlődő részekkel, azonos terhelés és peremfeltételek mellett. A szimmetria és az ismétlődésből származó periodicitás figyelembevétele a számítási igények lényeges csökkenéséhez vezet, mivel a teljes szerkezet viselkedését egy kisebb rész vizsgálatával is tisztázni lehet. Ehhez az adatelőkészítés kevesebb munkája is pozitívan járul hozzá. A szerkezet geometriájából, anyagából, terheléséből és megfogásából származó szimmetria miatt, a szimmetria felületein, vonalain, bizonyos kinematikai mennyiségek zérus értékkel rendelkeznek. Példaként szolgáljon egy olyan téglalap alakú lemez, amely mind a négy oldalán befalazott. A terhelés egyenletes nyomás a lemez teljes felületén. Az egyes oldalak mentén megfelezve a lemezt, a negyedrészének vizsgálatával célhoz érünk, ha a középvonalak mentén a vonalirányú szögelfordulást zérusra állítjuk be, azaz ez lesz a kinematikai peremfeltétel. Egy másik gyakori példa a forgó alkatrészek, szerkezeti elemek-hengerkoordinátarendszerbeli vizsgálata. Ekkor ebben a rendszerben a szerkezet periodicitással rendelkezhet. Pl. egy szivattyú járókereke. A lapátok közötti rész ismétlődik. Egy felvett R sugáron a lapátokat F és A pontban metsszük el. Ehhez a pontokhoz rendre a ϕ F és ϕ A hengerkoordinátarendszerbeli szögek tartoznak. E szögekkel kijelölt helyi koordinátarendszerben a radiális és tangenciális elmozdulások páronként azonosak, azaz u F = u A és vF = vA . Emiatt a periodicitási peremet (felületet) tartalmazó végeselemeknél azon csomópontokban, amelyek ezeken a peremeken helyezkednek el, az x − y rendszerből át kell transzformálni a mennyiségeket a helyi koordinátarendszerekbe, továbbá az F, A pontpár ismeretlenjeit egybe kell ejteni, majd ennek figyelembevételével kell az elemek illesztését elvégezni a végső egyenletrendszer előállítása céljából. Itt is az A pontot tartalmazó elemen kell végrehajtani az alárendelést, hasonlóan, mint az excentrikus csatlakozásnál is tettük. Jelen esetben az A és F pontok sorszámainak nagyobb távolsága miatt a végső egyenletrendszer sávszélessége nő, de az ismeretlenszám lényeges csökkenése e negatívumot kompenzálja.
3.7. Rugalmas ágyazás Amint már a síkbeli tartóknál is említettük, egyszerűsített modellezésnél a vizsgált test megtámasztását biztosító rugalmas testet, Winkler típusú közeggel helyettesítik. A nagymerevségű rúgóállandóval ( c ≈ 106 ), a rugóirányú elmozdulás zérus értékét jól be lehet állítani. A rugalmas megtámasztási felületen elhelyezkedő x, y, z tengelyirányú rugóállandókat jelölje cx , c y , cz , a keletkező felületi terhelést
p = c x , c y , c z [ u v w] = C u , T
(3.33)
e ahol diagonál mátrixot jelöl. Az e jelű elemmel az Arug felület mentén csatlakozó rugalmas támaszban felhalmozott alakváltozási energia e U rug =
1 1 pT u dA = ∫ uT C u dA , ∫ 2 Ae 2 Ae rug
(3.34)
rug
ami az elemenbelüli elmozdulás szokásos közelítése révén ue (x ) = Ne (x ) qe az
79
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
1 1 e U rug = q e,T ∫ NT C N dA q e = q e,T K erug q e 2 2 Ae
(3.35)
rug
végső alakot nyeri. Ezzel a rugalmas közeggel kapcsolódó végeselem teljes potenciális energiája e Πep = Πep (végeselem ) + U rug
(3.36)
alakban számítható.
80
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
4. I-DEAS RENDSZER HASZNÁLATA Az I-DEAS tervező rendszer különböző alkalmazások együttese, melyeket a tervezési folyamat különböző fázisainak megkönnyítésére alkalmazhatunk. Minden egyes gépészeti feladat elvégzése, más-más szoftverrész használatát követeli meg. A program elindítása után megnyíló ablak több részből épül fel, melyek közül a jobb szélsőben találhatjuk meg a különböző alkalmazások kiválasztását engedélyező listaablakot. Ezek például a Design, Simulation, Test, Manufacturing, stb., melyek a különböző feladatok elvégzésére szolgálnak, úgymint a • Design: ebben a szoftverrészben alkatrészek geometriai modelljének létrehozása, módosítása lehetséges. • Simulation: Az I-DEAS végeselemes modulja, olyan feladatok végrehajtása lehetséges, mint a peremfeltételek előírása, végeselemes háló létrehozása, végeselemes számítások eredményeinek elemzése, stb. • Test: Az időkezelés eszközeit tömöríti ez az alkalmazásrész. • Manufacturing: A gyártással kapcsolatos szimulációs programrész. Jelen fejezetben elsősorban az I-DEAS programrendszer Simulation moduljának használatát tekintjük át, de röviden a program működtetéséhez szükséges funkciókat, jellemzőket is megvizsgáljuk [5].
4.1. Általános jellemzők: 1. Parametrikus modellezés. A tervezés során először egy vázlatot kell készíteni, mely nagy vonalakban hasonlít majd az elkészítendő darabhoz, és a méreteket ezután kell pontosan beállítani a kívánalmaknak megfelelően. De természetesen a geometriai elemek pontos koordináták segítéségével is megrajzolhatóak. 2. Tulajdonság alapú modellezés. A bázis alak létrehozása után egyszerűen lehet definiálni kivágást, furatot, beszúrást, stb. 3. Párhuzamos alkatrészfejlesztés. Az alkatrészek közös könyvtárakban helyezhetőek el, melyek a megfelelő tervezők által elérhetők módosíthatók. Az I-DEAS elindítható a parancssorból, menüből vagy ikonnal. Előfordulhat, hogy a program használata speciális account-ot is igényel. A szoftver használatához ki kell választani a rendszerünk által támogatott – lehetőleg a legjobb – grafikai meghajtót, pl. OpenGL, PEX. Az elindítás után egy indító ablak nyílik ki, ahol a következő adatokat kell megadni: • Project neve: mely az adott munkát rendszerezi. Ezt ki is lehet választani a felkínált listából. Vagy behívható egy kiválasztó ablak, az ikonra kattintással. • Model file: a munka során létrehozott objektumhoz tartozó adatok itt tárolódnak el. Ezt segíti egy előhívható lista, mely a file megnyitáshoz, mentéshez hasonló ablakot jelent. • A használni kívánt alkalmazás kiválasztása: Alapértelmezésként felkínálja a program az utójára használtat, illetve a Design csomagot. Ez alatt található az adott alkalmazáson belüli feladat kiválasztására szolgáló legördülő listaablak. Ha az I-DEAS-t parancssorból indítottuk el, akkor lehetőség van megadni opciókat is. Pl: -h az indításhoz használható opciókat jeleníti meg, –d device a grafikus meghajtót lehet vele megadni indításkor. Ha nem adjuk meg a device nevet, akkor egy listát kínál fel, amiből lehet választani, –g a legutóbb végzett munka folytatását teszi lehetővé, –l
81
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
a használni kívánt nyelvet lehet megadni. Ha nem adjuk, meg akkor egy listát kínál fel az elérhető nyelvekkel.
language
4.1.1. Használathoz szükséges alapok A program különböző ablakokat kezel, melyek a következők • Rajzterület: Alapbeállításként egy fekete hátterű rész a monitor legnagyobb bal felső területén. Itt történik a rajzolás. • Ikon ablak: A monitor jobb oldalán húzódó terület, mérete tetszés szerint módosítható. Az itt található parancsikonok három nagy területre oszthatók: a felső 6x3 db ikon, a továbbiakban A mátrixként hivatkozunk rá, a középső 4x3 db ikon, a továbbiakban B mátrix, és az alsó 4x3 db ikon, melyet C mátrixnak nevezünk. • Lista ablak: Az üzenetek, hibák jelzésére használja a program, ha nem használjuk sokszor el lehet rejteni, de hasznos dolog. • Prompt ablak: A kiválasztott parancsnak megfelelően adatok bekérését ezen keresztül végzi a program. Az ikonok használata nem sok magyarázatot igényel. Egy kis gyakorlással könnyen elsajátítatható a kezelésük. Az ikonokról annyit azonban tudni kell, hogy a legtöbb ikon mögött több feladat elvégzése is elérhető. Erre utal, az ikon jobb alsó sarkában egy kis háromszög, jelezve, hogy további funkciók érhetőek el a gyűjtő kinyitásával. Gyors egérkattintással az ikon kiválasztásra kerül, és inverz színben jelenik meg. Ezzel aktiválható a jelzett funkció. Ha egy ikonon lenyomva tartjuk az egérgombot és egy ikon gyűjtőről van szó akkor kinyílik egy lista, melyek közül tetszőleges parancsot lehet indítani. Az egér ikonra pozícionálása megjeleníti az adott elem funkcióját a státuszsorban, mely a grafikus ablak legalsó sorában van. Fontos: Mikor az ikonokat, menüket használjuk érdemes figyelni az üzeneteket, a tájékoztató ablakban, melyet a rajzterület státuszsorában találunk. A program használatához három gombos egér használata az ideális, ahogy ez egy korszerű tervező szoftvertől elvárható. Minden gombnak saját funkciója van. A Bal gomb feladata a parancskiválasztás, a geometriai alakzatok kiválasztása a rajzterületen. A Shift billentyűvel együtt használva csoportos kijelölést tesz lehetővé. (ez pl. törlésnél, méretezésnél hasznos). A Középső gomb ez az Enter vagy a Return billentyűt helyettesíti, azaz a parancsok lezárására szolgál. Tehát a program használatát gyorsítja. A Jobb gomb egy ,,Popup” menüt jelenít meg, ha a rajzterületen használjuk, feladattól függően más-más parancsok aktivizálást gyorsítja. Számtalan billentyűkombináció előre definiált az I-DEAS-ban, melyek felsorolása itt túl hosszú lenne. Most csak az F1 - F12 billentyűket tekintjük. Ezek szerepe átdefiniálható (lásd az ideas.ini állományt), de alapértelmezésben a következő feladatokat gyorsítják: • F1 - F5: eltolás, nagyítás, forgatás, kívánt nézet, reset • F6: az előző 5 funkcióbillentyű szerepét határozza meg a feladatbank kiválasztással, melyet a program jelez a rajzterület jobb alsó sarkánál feltüntetett ikonokkal, háromféle parancsbank definiált. • F7: teljes méretűre állítja a létrehozott rajzot(AU vagy Ctrl-A, ZM: ablakkal nagyít) • F8: egymáshoz közel fekvő rajzelemek közötti választást segíti (Reconsider) • F9: rajzelemek kijelölését megszünteti (Deselect All) • F10: Munkasíkba állítja a rajzot. • F11: Egy dialógus ablakot nyit („Filter”), mely a rajzelemek gyors kiválasztását segíti. • F12: a rajzterület újrarajzolása (Redisplay Ctrl-R) Az F1-F3 funkcióbillentyűk által definiált művelet elvégzését a kívánt billentyű nyomva tartásával és az egér mozgatásával érhetjük el.
82
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A Menü elérése a Ctrl-M kombinációval történik, mely ki/bekapcsolja a menü megjelenítését. Kilépés: a parancssorba írt exit utasítással, vagy a menüből kiválasztva, vagy a Ctrl-E billentyűkombinációval.
4.1.2. Rajzolást könnyítő funkciók Dynamic Navigator – Dinamikus navigátor, mely az alkatrészrajz készítése során nyújt támogatást. A már létező rajzelemekhez viszonyított tulajdonságokat jelzi a program a 4.1. táblázatban felsorolt jelzésekkel. Rajzolás közben az egér jobb billentyűje segítségével további funkciókat aktivizálhatunk, pl. az Align, vagy a Focus parancsokat, melyek egy-egy korábban létrehozott görbéhez való kapcsolást tesznek lehetővé. 4.1. táblázat: Dinamikus rajzolást segítő elemek érintő
végpont
középpont
metszéspont
párhuzamos
függőleges
vízszintes
egybeesés
A rajz készítése közben hasznos a ki/bekapcsolható Grid háló, vagy Snap funkció. A Grid egy általunk definiált diszkrét ponthálót jelenít meg, mellyel a mérethelyes rajzolást könnyíti a program. Bekapcsolt Snap a létrehozandó rajzelemek pontjait, csak a meghatározott diszkrét pontokba engedi elhelyezni. Mindkét parancsot a B mátrix második sorának utolsó oszlopában találjuk a Workplane Appearence parancsot, vagy röviden B(2,1) helyen.
4.2. Rajzolás az I-DEAS-ban program legfontosabb modulja, mely alapvetően a rajzolást támogatja a Master Ennek segítségével tudjuk az alkatrészeket megtervezni, illetve ez ad lehetőséget a két- illetve a háromdimenziós modellek elkészítéséhez, melyeknek például a Simulation modul segítségével a végeselemes analízisét végre tudjuk hajtani. A rajzterületen elindítás után egy kijelölt rajzsíkot látunk, alapértelmezés szerint az x-y sík, de meg lehet változtatni az igényektől függően. A Master Modeler-ben elérhető funkciók a jobb szélső ikon panel felső csoportjában – az A mátrixban, lásd a 4.1. ábrán– találhatók. Ezek tulajdonképpen a rajzelemek létrehozását végzik. Háromszor hat darab gyűjtő ikonból áll, de minden gyűjtőben további funkciók aktivizálhatók. Az
I-DEAS
Modeler.
4.1. ábra. Az A, a B és a C ikongyűjtők a Master Modeler alkalmazásnál Ismertetésük most nem cél, de egy kis gyakorlással, mindenki könnyen megismerheti a különböző parancsok elhelyezkedését. A tanulást segíti az is, hogy az ikon gyűjtők kinyitásakor – az ikont hosszú egérkattintással kell kiválasztani – az ikonok mellett megjelenik a parancs elnevezése is, és a státuszsorban olvasható a parancs rövid leírása. A használat során megfigyelhető, hogy a legutóbb használt parancs ikonja a felső kiválasztott helyzetbe kerül, mely gyorsítja a rajzolást. 83
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
Megjegyzés: Ha valaki nem érti a rövid parancsleírásból a kiválasztandó funkciót, akkor a Help menüben az On Context menüpont segítségével további információt kaphat. Az A mátrixban elérhető funkciók felsorolása: • Rajzsík kijelölő menü: egy tetszőleges sík kiválasztása, munkaasztal kijelölés A(1,1) • Koordináta rendszer menü: referencia sík, pont, vonal A(1,2) • Metszetek menü A(1,3) • Vonalak menü: sokszög, vonal, négyszög, pont létrehozása A(2,1) • Körív menü: különböző körívek rajzolása A(2,2) • 3D-s menü: háromdimenziós rajzelemek létrehozásához A(2,3) • Kör menü: Teljes kör rajzolása, különböző módszerekkel A(3,1) • Görbevonal parancsok: spline-ok, ellipszisek A(3,2) • Leképzés menü: eltolás, leképzések készítése A(3,3) • Méretező parancsok: A(4,1) • Lekerekítés, letörés menü: trimmelés, szétvágás, sarkok készítése A(4,2) • Felület menü: felület kiterjesztés, metszés, stb. A(4,3) • Extrudálás menü: felületek létrehozása, testek extrudálása A(5,1) • 3 D-s lekerekítés, letörés menü: A(5,2) • Halmaz műveletek: metszés, unió, különbség, stb. A(5,3) • Kiosztás menü: négyszög kiosztás, körkiosztás, skálázás, A(6,1) • Szabad felület, él menü A(6,2) • Jellemzők menü: A(6,3) A középső 12 ikon, illetve ikongyűjtő funkciója elsősorban a rajzelemek módosításával, szabályozásával kapcsolatos. A B mátrixban elhelyezkedő parancsokat az alábbiak szerint találjuk, ha a Simulation / Master Modeler alkalmazást használjuk: • Történeti fa B(1,1) • Mozgatás, elforgatás, elrendezés, stb. B(1,2) • Megjelenítés szabályozása, elrejtés, stb. B(1,3) • Módosítás, Undo B(2,1) • Informácó lekérdezése B(2,2) • Megjelenés szabályozása, munkaterület mérete B(2,3) • Újrarajzolás vezérlése B(3,1) • Mérés, jellemzők, anyagok, stb. B(3,2) • Alkatrészek, jellemzők, B(3,3) • Törlés B(4,1) • Alkatrész katalógus kezelése, elnevezés, csoportosítás, stb. B(4,2) • Alkatrész könyvtár kezelés B(4,3) A modell létrehozása során, az I-DEAS minden egyes lépést egy történeti fa segítségével tárol, mely megtekinthető és szerkeszthető. A történeti fa csomópontokat tartalmaz, melyek két gyermekből és egy szülőből épülnek fel. A csomópontok kiválasztásával a grafikus ablakon nyomon lehet követni, hogy melyik elemről van szó. A történeti fa elsődleges szerepe az alapelemek módosíthatósága. A módosítás parancsot a B(2,1) ikonnal indítjuk el, a méretszámok kiválasztása után pedig a felnyíló dialógusablakban módosíthatjuk. Ha nincsenek kint méretek, akkor tetszőleges méretet el tudunk helyezni a méretező segítésével.
84
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A C mátrixban elhelyezett parancsok minden modulnál egységesen megtalálhatóak, ezek a következő megjelenítéssel kapcsolatos feladatokat látják el: • Újrarajzolás parancs C(1,1) • Vonalas ábrázolás ikongyűjtő C(1,2) • Árnyalt megjelenítés C(1,3) • Teljes méret ikongyűjtő C(2,1) • Nagyítás-kicsinyítés parancsok C(2,2) • Rajzfelület menü C(2,3) • Nézet ikongyűjtők C(3,1), C(3,2), C(4,1), C(4,2) • Leállító parancsikon C(3,3) • Nyomtatás parancsok C(4,3)
4.3. Végeselemes analízis Az I-DEAS program végeselemes moduljának a Simulation nevet adták. Ebben a programrészben lehetősége van a felhasználónak tetszőleges peremértékfeladat felállítására, megoldására és a megoldás elemzésére. Itt most csak röviden utalunk rá, hogy melyik feladat elvégzése, melyik programrészben lehetséges. Elsősorban a következő alrészek áttekintését tűzzük ki: • Simulation / Boundary Conditions, a peremfeltételek definiálása, • Simulation / Meshing, a végeselemes háló kialakítása, • Simulation / Model Solution, a peremértékfeladat megoldása, • Simulation / Post Processing, a kapott megoldások vizsgálata.
4.3.1. Peremértékfeladat kitűzése A Simulation modul Boundary Conditions alkalmazás kiválasztása után az A mátrix a 4.2. ábrán jelzett formában jelenik meg. A peremértékfeladat típusának kiválasztása az A(1,1) parancsok segítségével történik. Dinamikai peremfeltételeket az A(2,1) és A(2,2) ikongyűjtőben található utasításokkal írhatunk elő. Az 4.2. ábra. Az A mátrix A(3,2) parancsgyűjtőben találhatók a hőmérsékletmező előírását szolgáló parancsok. A modell szabadságfokainak rögzítését az A(4,2)-ben található ikonok szolgálják. A peremfeltételek nyilvántartását és kezelését szolgáló két utasítás a Sets… és Boundary Conditions… az A(6,2) és A(6,1) helyen van elhelyezve.
4.3.2. Végeselemes háló létrehozása A végeselemes háló létrehozása jelenti a következő fontos lépést a végeselemes analízis során. Ez a lépés azonban megelőzheti a peremfeltételek előírását, mint látni fogjuk a gyakorlatban. A végeselemes hálózat nem csak egy parancskiválasztást jelent, mivel lényeges paraméterek beállítását el kell végezni, úgymint elemtípusok, anyagjellemzők, geometriai jellemzők (pl. héjak esetén falvastagság), csomópontok, elemek konkrét helyen történő definiálása stb. 4.3. ábra. Az A mátrix A Simulation / Meshing modul elindításával az A mátrix a 4.3. ábrán jelzett formában jelenik meg. Az itt kirajzolt ikonok mindegyike további parancsokat takar, melyek elérése a szokásos
85
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
módon történik. Lényeges parancsok az A(5,1) Materials… vagyis anyagjellemzők, A(5,2) Physical Properties… vagyis fizikai tulajdonságok továbbá a hálózást segítő parancsok az A mátrix első sorában.
4.3.3. Feladat megoldása A Simulation / Model Solution alkalmazás választásával jutunk a peremértékfeladat megoldását segítő alkalmazásokhoz. Ez a modul végzi tulajdonképpen a végeselemes számítást, a korábbi lépesek során létrehozott végeselemes modellen. A számítás elindítása előtt néhány beállítást el kell végezni, úgymint az eredmények tárolására alkalmas hely kijelölését, a megoldás módját. Erre szolgál a jelen alkalmazásban a Solution Set… A(1,2) parancs. Itt állíthatjuk be megoldás során szükséges ideiglenes tárolási helyet, a megoldás során alkalmazandó pontossági elvárásokat. A végeselemes feladat kapcsán előállított lineáris algebrai egyenletrendszert az I-DEAS alapvetően kétféle a felhasználó által választható, egyenletrendszer megoldóval, pontosabban egy direkt és egy iteratív technikával képes megoldani. A program dokumentációja szerint a direkt megoldó szinte minden jól felállított feladatot pontosan képes megoldani, bár egy kicsit lassabban, mint az iteratív technikával dolgozó alprogram.
4.3.4. Eredmények megjelenítése Lehetőség van még a Model Solution alkalmazás keretein belül is megtekinteni a kapott megoldást, az A(6,2) alatt található Visualiser segítségével, azonban az eredmények mélyebb elemzése megkívánja egy újabb alkalmazás a Simulation / Post Processing kiválasztását. A megjelenő A mátrix (lásd a 4.4. ábrán) ebben az esetben csak az eredmények kezelésének megfelelő parancsokat tartalmazza. Az A(1,1) Results… ikon a megjeleníteni kívánt eredmény kiválasztását segíti. A felhasználó tetszés szerinti eredményskálát, színezési vagy festési módot állíthat be. Akár elemenkénti eredmény kirajzolást is előírhat. A program képes arra, hogy a kiszámított eredményeket animálja, mely jelentheti a kialakuló végső alak elérésének szemléltetését, vagy a feszültségállapot létrejöttét. A különböző eredmények együttesen különböző skálázással is megtekinthetők, a felhasználói igények szerint. Természetesen a program az eredmények kirajzolása során támogatja a különböző 4.4. ábra. Az A mátrix transzformációkat is. Figyelembe kell azonban venni, hogy a viszonylag sűrű végeselem hálózat, illetve a túl sok csomópont a kirajzolást korlátozhatja.
86
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
5. C ÁLLVÁNY VIZSGÁLATA Szerszámgépállványok egy lehetséges kialakítási módját vizsgáljuk meg végeselem módszer segítségével. Jelen esetben a kétdimenziós végeselemek egy típusának használatát kívánjuk áttekinteni, annak érdekében, hogy meghatározzuk a felvett peremértékfeladatnál kialakuló alakváltozást és feszültségmezőt. A numerikus vizsgálat elvégzéséhez az I-DEAS programrendszert használjuk a 2D-s síkfeszültségi modell elemzéséhez. A feladathoz tartozó geometriai kialakítást a 5.1. ábra szemlélteti.
5.1. ábra. C–állvány geometriai adatai A szerkezet kialakítása hegesztéssel történik, melynek a hatásait jelen modellben elhanyagoljuk. A szerszámgépállvány talajhoz való rögzítése csavarozással történik. A vizsgált fém keret vastagsága 25 mm . Az anyagát jellemző paraméterek a következők; rugalmassági modulus: E = 2 ⋅105 MPa , Poisson tényező: ν = 0.3 . A végeselemes számítások célja, hogy meghatározzuk a legnagyobb terheléskor kialakuló elmozdulásmezőt, melyből számítható a szerkezetben ébredő feszültség eloszlása is. Ennek minden valóságos esetben egy meghatározott korlát alatt kell lennie. Így tervezéskor előre lehet tudni azt, hogy mely helyeken és hogyan kell esetlegesen módosítani a szerkezetet a kívánt feszültségi korlátok betartása miatt. Az I-DEAS használata során a következő főbb lépéseket kell megtenni a feladat végrehajtásához, melyeket a Simulation modulrészben találunk meg: • Geometriai modell létrehozása: a Master Modeler segítségével történik. • A peremértékfeladat definiálása: a Boundary Conditions által. • Végeselemek létrehozása: a Meshing programrészben. • Feladat numerikus megoldása: a Model Solution alkalmazásban. • Az eredmények értékelése: a Post Processing programrésszel. A lépések végrehajtási sorrendje kötött, azonban a végeselemes háló kialakítása és a peremfeltételek előírása tetszőleges sorrendben történhet. Természetesen, ha a geometriai modellen helyezzük el a peremfeltételeket tehát az ott definiált élekhez, pontokhoz rendelünk kinematikai vagy dinamikai peremfeltételeket akkor az I-DEAS rendszer ezt átalakítja – a végeselem hálózat létrehozásakor – oly módon, hogy a peremfeltételek egyértelműen a csomópontokhoz legyenek társítva.
87
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
5.1. Geometria létrehozása A 5.1. ábrán vázolt szerkezet létrehozása a célunk az Irendszerben, melyhez a Simulation modul Master Model alkalmazását válasszuk ki, ahogy a 5.2. ábra szemlélteti. DEAS
5.2. ábra. Geometria definiálása
A program alapbeállításait módosítani kell. Egyrészt az Options menüben a Units parancsot válasszuk ki, és a kinyíló Popup menüben a mm[newton] lehetőséget jelöljük ki. Ezáltal a program a hosszúsági mérteteket mm-ben, az erőket N-ban várja el, illetve az eredményeket is ennek megfelelően mutatja. Másrészt át kell állítani a munkaterület méretét is. Ehhez válasszuk a Workplane Appearance parancsot, melyet a B(2,3) ikongyűjtőben találunk. A kinyíló dialógusablakban beállíthatóak a rajzterület jellemző X és Y méretek, ahogy a 5.3. ábra mutatja.
5.3. ábra. Rajzterület méretbeállítása A C–állvány külső kontúrjának létrehozásához a Polylines parancs az A(2,1) ikongyűjtőben található. A méretek pontos beállításához a Modify Entity B(2,1) utasítás szükséges, melyet a méret kiválasztása után kell elindítani. A méretek beállítása után mentsünk: Ctrl+S. Következő lépésként a külső kontúr alapján létre kell hozni egy felületetet, melyből a kitöréseket kell majd kivágni. A Surface by Boundary A(5,1) parancs a külső kontúr kijelölése után létrehozza a kívánt felületet. A kivágandó négyzet és kör kontúrok létrehozásához szükséges a Center Edge A(3,1) körrajzoló és Polyline A(2,1) parancs. A méretek beállítása után válasszuk a Trim at Curve A(4,1) utasítást, mely ezen utóbb létrehozott kontúrokat vágja ki a korábban definiált felületből. Végül nevezzük el a létrehozott modellt, a Name Parts B(4,2) utasítással, pl. callvany-ra. Ezzel készen van a geometriai modell, melyet a továbbiakban használni fogunk a végeselemes modell létrehozásához. Mielőtt tovább haladunk mentsünk: Ctrl+S.
5.2. Végeselemes modell Térjünk át a Simulation modul Meshing programrészbe (lásd. a 5.4. ábrát). Létrehozunk egy véges-elem modellt a Create FE Model B(4,2) utasítással. Az anyagjellemzők beállítása a Materials A(5,1) paranccsal lehetséges, ahol az új anyagtábla létrehozását kell 5.4. ábra. Végeselem hálózás elindítani. Ekkor jutunk a 5.5. ábra által bemutatott dialógus ablakokhoz, ahol beállítjuk a kívánt anyagjellemzőket.
88
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
Ezt követően megadjuk a 2D-s elemekhez tartozó fizikai jellemzőket a Physical Property A(5,2) parancs segítségével. Itt 2D-s Plane Stress táblázatot kell kitölteni, szem előtt tartva, hogy a C–állvány vastagsága 25 mm (lásd a 5.6. ábrát).
5.5. ábra. Anyagtábla létrehozása Ezután létrehozzuk a végeselem hálót, azaz definiáljuk a csomópontokat és a végeselemeket. Ehhez a Create Shell Mesh A(1,1) parancs elindítása szükséges, mely hatására megjelenik a 5.6. ábra által mutatott dialógus ablak.
5.6. ábra. Fizikai jellemzők beállítása A végeselem hálózat létrehozása ebben a feladatban a következő paraméterek beállításával történt: Element Length: 100, Element Family: Plane Stress, Element Type 6 csomópontú háromszögelem. Továbbá kiválasztásra került a már korában létrehozott anyag és fizikai jellemzők táblája.
89
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
5.7. ábra. Végeselem háló létrehozása
5.8. ábra.
A létrehozott végeselem hálózatot a 5.8. ábra szemlélteti. A végeselem modellről ezek után már csak a peremfeltélelek előírása hiányzik. Ehhez át kell lépnünk a Simulation modul Boundary Conditions alkalmazásába. Mielőtt azonban továbbhaladunk végezzünk egy mentést: Ctrl+S. A peremfeltételek közül először a kinematikai peremfeltételeket adjuk meg, ehhez válasszuk a Displacement Restraint A(4,1) utasítást. A rögzíteni kívánt csomópontokat kijelölve és Done
Végeselem hálózat paranccsal elfogadva jutunk a 5.9. ábrán jelzett dialógusablakhoz. Itt kiválasztható, hogy milyen típusú megfogásokat kívánunk alkalmazni az adott helyeken. Jelen feladat kapcsán a C–állvány bal oldalán Clamp míg jobb oldalán Slider(X) típusú rögzítést adtunk meg. A dialógus ablak bezárásával a végeselemes hálón megjelennek a megfogást jelző szimbólumok, ahogy ezt az 5.8. ábra mutatja. A dinamikai peremfeltételek megadása a Force A(2,1) paranccsal történik. A csomópontok kijelölése majd elfogadása után megadjuk az adott helyen működő terhelés nagyságát. Jelen feladatban a lefelé, illetve felfelé működő erők eredője is 250 kN. A továbblépés előtt mentsünk: Ctrl+S.
5.9. ábra. Kinematikai peremfeltételek előírása
5.3. A feladat megoldása
5.10. ábra. Feladat megoldás
A felállított peremértékfeladat megoldása a Simulation modul Model Solution alkalmazás segítségével végezhető el. Itt válasszuk a Solution Set A(2,1) parancsot, mellyel létrehozzuk a megoldás tárolásához szükséges táblázatot. Ezzel rendeljük hozzá a feladathoz a megoldási módszert, kiválasztható a tárolni kívánteredmény. Ha többféle peremfeltételt is definiáltunk akkor ezek közül is ez a dialógusablak – lásd a 5.11. ábra-át – engedi meg a választást. Ezeket a
90
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
táblázatokat a program mentés után eltárolja, és ha már került bele megoldás, akkor a peremfeltételek nem módosíthatók egy esetleges újabb futtatás miatt addig, amíg a táblázatokat nem töröljük.
5.11. ábra. Megoldáshoz szükséges paraméterek beállítása A feladat megoldását a Solve A(2,1) parancs indítja el. A számítási folyamatról a program az I-DEAS List ablakba üzen. Sikeres befejezéskor a program a No warnings no errors encountered in last run üzenetet adja.
5.4. Számítási eredmények Az eredmények megtekintéséhez térjünk át a Simulation modul Post Processing alkalmazásra. Ezen programrész számtalan segédeszközt ad a számítási eredmények elemzésére. A Results A(1,1) szolgál a megjeleníteni kívánt eredmény kiválasztására, míg a Display Template A(1,2) és Color Bar A(2,2) parancsok a kirajzolási paraméterek beállítását segítik. A Display A(2,1) utasítással indítjuk el a kirajzolást.
5.12. ábra. Kiszámított elmozdulásmező abszolút értéke
91
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
5.13. ábra. Redukált feszültségmező A program által kiszámított elmozdulás maximális értéke 8.7 mm, míg maximális redukált feszültség 293 MPa, ahogy ezt a 5.12. és 5.13. ábrák szemléltetik.
92
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
6. TENGELYSZIMMETRIKUS FELADAT ELEMZÉSE 6.1. A peremértékfeladat Vizsgálatunk tárgya egy belsőnyomással terhelt, függőlegesen álló alumínium gázpalack. A szerkezet geometriája és terhelése egyaránt tengelyszimmetrikus. Ebből következik, hogy a fellépő mechanikai jellemző mennyiségek függetlenek a forgásszögtől. Azaz az eredetileg háromdimenziós feladat valójában kétdimenziósra redukálódik. Ilyen esetben elegendő a feladat egyetlen meridián metszetét vizsgálni. Ez igen nagypontosságú elemzést tesz lehetővé, hiszen egyetlen metszeten igen sűrű felosztás vehető fel relatíve kis ismeretlenszám mellett. A tengelyszimmetrikus feladatokat az I-DEAS programrendszer a z forgástengelyű x-z koordinátarendszer első és negyedik síknegyedében elhelyezett meridián metszetén értelmezi. A szerkezet geometriáját a jellemző méretek feltüntetésével a 6.1. ábra szemlélteti. A tartály egész metszetét egyetlen zárt felületként definiáljuk. A tartományt 6 csomópontú, kvadratikus approximációt biztosító háromszög alakú elemekre osztjuk fel. A látszólag síkbeli elemek mechanikai szempontból - a metszetnek megfelelő gyűrű elemeket írnak le. A tartályt 3 mm falvastagságú alumínium lemezből készült. Az alumínium (AL_2014) ötvözet lineárisan rugalmas anyagjellemzői adottak: E=77000 MPa, v=0.33. A tartály belső falán 7 bar túlnyomást írunk elő. Az elvégzendő számítás után az eredmények közül megtekinthetjük többek között az elmozdulást és a feszültségeket.
6.1. ábra. A gázpalack geometriai méretei
6.2. A geometria létrehozása Az I-DEAS Open GL ikonjára kattintva elindítjuk a programrendszert, majd beírjuk a feladat nevét és kiválasztjuk a megfelelő alkalmazást és feladatot: Model File name: csocsat Application: Simulation Task: Master Modeler 93
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A kívánt dimenzió beállítása: Options/units/mm[newton]
Beállítjuk a munkaterület méreteit és a raszter sűrűséget:
6.2. ábra. Munkaterület paraméterei A munkaterületet át kell, hogy helyezzük az x-z síkba. Ehhez először alkalmazzuk az A(1,2) Coordinate System, majd az A(1,1) alatti Sketch in Place utasításokat az x-z sík kiválasztásával. Továbbá célszerű a C(3,1) menü pont alatt a Bottom View utasítást választani. Ekkor a z tengely függőlegesen felfelé fog állni. Először célszerű a z-tengely mentén egy 500 mm hosszú segédegyenest megrajzolni az A(2,1) Line, A(4,1) Dimension és a B(2,1)Modify parancsok segítségével. Ezen egyenes végpontjai kijelölik a palackot lezáró gömbök, azaz a körök középpontjait. A A(3,1) paranccsal két-két centrikus kört hozunk létre a kijelölt középpontokban. Az egyik legyen 300 mm átmérőjű a második 306 mm. Majd az A(2,1)Line paranccsal rajzolunk egy vízszintes és egy függőleges egyenest, amelyek átmennek a körök középpontjain. Ezután képezzük a körök és egyenesek metszését az A(4,2) alatti Divide At utasítással. A felesleges köríveket a B(4,1)-el töröljük. A további egyeneseket és lekerekítést a 6.1. ábra szerint az A(2,1) Line, A(4,1) Dimension, valamint az A(4,2) Fillet és a B(2,1)Modify parancsok segítségével értelemszerűen megrajzoljuk. A görbék metszése után megmaradt felesleges szakaszokat a B(4,1)Delete paranccsal töröljük. A peremgörbék egy egyszeresen összefüggő zárt felületetet határoznak meg, amelyet az A(5,1) Surface by Boundary utasítással hozunk létre. Az utasítás végrehajtása során, a határoló görbéket a választott körbejárás iránynak megfelelően egymásután kiválasztjuk.
6.3. Végeselemes hálózás A feladatok közül kiválasztjuk a Meshing menüpontot, és létrehozzuk a végeselemes feladatot a B(4,2)Create FE Modela utasítással.
6.3. ábra. Végeselemes felosztás definiálása Az A(1,1)Define Shell Mesh utasítással a tartományt tengely-szimmetrikus háromszögelemekre osztjuk (lásd a 6.3 ábrát). Az átlagos elem méretét 1 mm -re, az elem családot Axisymmetric Solid –ra és az elem típusát 6 csomópontú háromszögelemre állítjuk. A megtekintő ikonra kattintva megszemlélhetjük és megerősíthetjük a végeselem-hálót.
94
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A palack anyagát A(5,1)Material ikonjának meghívásával definiáljuk. Ekkor a 6.4. ábrán megjelenő Materials ablakban a felkínált acél anyagát átnevezzük AL2014-re, majd a mellette lévő módosító ikonra kattintunk és a megfelelő anyag-jellemzőket beállítjuk.
6.4. ábra. Az anyag tulajdonságok definiálása
6.4. peremfeltételek előírása Áttérünk a
menüpontra és a lineáris szilárdsági peremértékfeladatot az utasítással választjuk ki. A könnyebb kezelhetőség érdekében célszerű letiltani a végeselemes háló megjelenítését. Ezt a B(1,3)Display Filter utasítás alatt hajtjuk végre a Fem Display előtti kapcsoló üresre állításával. A szerkezetet a palackhoz kapcsolódó szoknya alján támasztjuk meg. Az A(4,2) Displacement Restraint gombra kattintva és a megfelelő vonalszakaszt kiválasztása után a 6.5. ábrán látható ablak jelenik meg. A Set All Free utasítással először minden szabadságfokot felszabadítunk, majd z irányban korlátozást írunk elő. A tartály belső falán 7 bar nyomáskülönbség hat. Az A(2,1) Force parancs aktivizálása után a megfelelő vonalszakaszokat megfogjuk és a 6.6. ábrán látható ablakban felkapcsoljuk az Intensity menüpontot, és az In Plane Force adatmezőbe MPa-ban kifejezett 0.7 nyomást írunk. Ezzel a peremérték feladat kitűzése befejeződött, a 6.7. ábrán láthatjuk a terhelést és megfogást is tartalmazó szerkezetet. Ezután a megoldás aktivizálása következik. Boundary Condition
A(1,1)Linear Statics
6.5. ábra. A szerkezet megtámasztásának előírása
95
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
6.6. ábra. Belsőnyomás megadása
6.7. ábra. A peremfeltételek megjelenítése
6.5. A feladat megoldása A feladatok közül kiválasztjuk a Model módszerrel oldjuk meg:
Solution
menüpontot. A feladatunkat lineáris
A(1,1) Linear
Hivatkozunk a kijelölt peremfeltételekre: A(1,2) Solution Set [Create]
A megoldó elindítása A(2,1) Solve
6.6. Az eredmények kiértékelése A feladatok közül kiválasztjuk a Post Processing menüpontot. Az eredmények megtekintése során célszerű elrejteni nemcsak a végeselemhálót hanem a szerkezeti rajzokat is. A rajzok elrejtése: B(1,2) Display Filter [ ]Wireframe [ ]Parts [ ]Assembly [ ]Work Plane [ ]Model views Border [ ] FE Model
az elrejtést az objektumok előtti megjelölés megszüntetésével érhetjük el. A szemléltetendő eredmények kiválasztása: Elmozdulás aszolút értéke
96
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A(1,1) Results [1-B.C. 1,DISPLACEMENT…] [Magnitude] [1-B.C. 1,DISPLACEMENT…] A(2,1) Display Közép Billentyű
Az elmozdulást 6.8. ábra mutatja. Redukált feszültség A(1,1) Results [3-B.C. 1,STRESS…] [von Mises] [1-B.C. 1,DISPLACEMENT…] A(2,1) Display Közép Billentyű
A redukált feszültséget a kinagyított szerkezeti részeken a 6.9. ábra szemlélteti. A többi feszültségi komponens és a főfeszültségek hasonló módon megjeleníthetők és ezzel a feladat teljes feszültségi analízise elvégezhető.
6.8. ábra. Az elmozdulás abszolút értéke [mm]
97
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
6.9. ábra. Redukált feszültség eloszlása [MPa]
98
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
7. TENGELYSZIMMETRIKUS FELADAT ELEMZÉSE 7.1. A peremértékfeladat ismertetése Egy-egy speciális esettől eltekintve a csövek csatlakozásait gyakran héjfeladatként modellezhetjük. Itt egy belső nyomással terhelt két különböző átmérőjű csőből álló rendszer szilárdsági analízisét mutatjuk be. A 7.1. ábrán egy 200 mm ámérőjű hajlított csőbe egy 100 mm átmérőjű szintén hajlított cső csatlakozik. A kétféle cső geometriája az alábbi módon van kialakítva. A nagyobbik 200 mm átmérőjű cső egyenes szakaszának hossza 1000 mm, amely egy 90 fokos könyökben folytatódik. A könyök középvonalának sugara 500 mm. A kisebb 100 mm átmérőjű cső merőlegesen csatlakozik a nagyobb átmérőjű csőbe, úgy, hogy a tengelyeik is merőlegesen metszik egymást. Az egyenes szakasz a csövek áthatásánál kezdődik, a végpontja a középvonalak metszésétől mérve 250 mm távolságra van. Ez a cső is 90 fokos könyökben folytatódik, amelynek középvonala egy 250 mm sugarú körívet ír le. A könyök után a cső 500 mm hosszú egyenes szakasszal rendelkezik. A 200 mm átmérőjű cső falvastagsága 20 mm, a 100 mm átmérőjűé 14 mm. A csövek anyaga acél. Az acél rugalmassági modulusa E= 200000 MPa, a Poisson tényező v=0.33. A csövek végei be vannak falazva, azaz az elmozdulási és a szögelfordulási koordináták egyaránt kötöttek. A csöveket 50 bar belsőnyomás terheli. A szerkezet megrajzolása során két alkatrészt (Part-ot) hozunk létre egyiket a 200 mm átmérőjű cső alkotja, másikat a 100 mm átmérőjű. A két Part-ból egy összeállítást (Assemblyt) építünk fel és az áthatást, ez utóbbin hajtjuk végre. A végelemes feladatot is az Assembly– hez rendeljük. A megfelelő héjvastagságot előírjuk, és az elemhálózást létrehozzuk az összeállításon. A szerkezet szabad peremein az elmozdulást és a szögelfordulást egyaránt megakadályozzuk, a felületeken pedig a belsőnyomást írjuk elő. A kitűzött peremértékfeladat megoldása után az eredményeket szemléltetjük.
7.1. ábra. Egy 200mm és egy 100mm átmérőjű cső csatlakozása
99
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
7.2. A geometria felépítése Az I-DEAS Open GL ikonjára kattintva elindítjuk a programrendszert, majd beírjuk a feladat nevét és kiválasztjuk a megfelelő alkalmazást és feladatot: Model File name: csocsat Application: Simulation Task: Master Modeler A kívánt dimenzió beállítása: Options/units/mm[newton]
A munkasík kijelölése: B(2,3) Workplane Appearence [-500] [-500] [ 500] [ 500]
200 mm átmérőjű kör rajzolása és átmérőjének pontosítása: A(3,1) Center Edge B(2,1) Modify Entity [200]
1000 mm hosszú henger extrudálása: A(5,1) Extrude [1000]
Függőleges segédegyenes rajzolása az eredeti munkasíkon a kör középpontjától 500 mm távolságra lefelé, valamint méretének pontosítása: A(2,1) Lines/Egér jobb gomb/Focus (kör középpontját kiválasztjuk)/Locate start/ Locate end/ A(4,1) Dimension [500]
Egy vízszintes segédegyenes (forgástengely) rajzolása A(2,1) Lines
A csőkönyök létrehozása forgatással
A 200 mm átmérőjű
A(5,1) Revolve Angle [90] *Protrude cső Part elnevezése: B(4,2) Name Parts [D200]
Függőleges segédegyenes rajzolása az eredeti munkasíkon a kör középpontjától 250 mm távolságra felfelé, méretének pontosítása: A(2,1) Lines/Egér jobb gomb/Focus (kör középpontját kiválasztjuk)/Locate start/ Locate end A(4,1) Dimension [250]
Vízszintes segédegyenes rajzolása A(2,1) Lines
Segédsík extrudálása a megrajzolt egyenes segítségével: A(5,1) Extrude [500]
A munkasík áthelyezése a segédsíkra: A(1,1) Sketch in place
A síkra merőleges nézet kiválasztása: C(3,1)
100mm átmérőjű kör rajzolása az új munkasíkon A(3,1) Center Edge B(2,1) Modify Entity [100]
250 mm hosszú henger extrudálása lefelé: A(5,1) Extrude *Newpart [250]
100
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A 100 mm átmérőjű cső Part elnevezése: B(4,2) Name Parts [D100]
A munkasík áthelyezése az extrudált hengert lezáró felületén: A(1,1) Sketch in place
Merőleges segédegyenes rajzolása az eredeti munkasíkon a kör középpontjától 250 mm távolságra kifelé, azaz a 200 mm átmérőjű egyenes csőre merőlegesen, méretének pontosítása: A(2,1) Lines/Egér jobb gomb/Focus (kör középpontját kiválasztjuk)/Locate start/ Locate end A(4,1) Dimension [250]
Vízszintes segédegyenes (forgástengely) rajzolása: A(2,1) Lines
A csőkönyök létrehozása forgatással: A(5,1) Revolve Angle [90] *Protrude
500 m hosszú egyenes csőszakasz rajzolása: A(5,1) Extrude *Join [500]
Az összeszerelés feladat kiválasztása Master Assembly
Összeszerelés: A(1,1) Hierarchy a(1,1) Name/[Assembly1] a(2,1) Add to Assemnbly [D200, D100]
Az áthatás végrehajtása: A(5,3) Join A(5,3) Cut
7.3. Végeselemes hálózás A feladatok közül kiválasztjuk a Meshing menüpontot. A végeselemes feladat létrehozása: B(4,2) Create FE Model [Assembly1] [Part3]
Hálózás átlagosan 15 mm méretű 6 csomópontú háromszög héjelemekkel: A(1,1) Define Shell Mesh(felületek kiválasztása) [15] [Thin Shell] [6 csomópontú háromszögelem] [@@]megtekintés/Keep Mesh
Falvastagságok megadása: A(7,2) Surface Thickness [20] A(7,2) Surface Thickness [14]
7.4. Peremfeltételek előírása A feladatok közül kiválasztjuk a Boundary A hálózat elrejtése:
Condition
B(1,2) Display Filter [ ] FE Model
101
menüpontot.
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
Befalazás előírása: A(4,2) Displacement Restraint [OK]
50 bar belső nyomás előírása a csövekben: A(2,2) Pressure (felületek kiválasztásával) [5]
7.2. ábra. A végeselemes háló 15 mm átlagos méretű 6 csomópontú háromszögelemeket tartalmaz
7.3. ábra. A csövekben 5 MPa belsőnyomást, a csővégeken befalazást írunk elő
7.5. A feladat megoldása A feladatok közül kiválasztjuk a Model Solution menüpontot. A feladatunkat lineáris módszerrel oldjuk meg: A(1,1) Linear
Hivatkozunk a kijelölt peremfeltételekre: A(1,2) Solution Set [Create]
A megoldó elindítása
102
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A(2,1) Solve
7.6. Az eredmények szemléltetése A feladatok közül kiválasztjuk a Post Processing menüpontot. Az eredmények megtekintése során célszerű elrejteni nem csak a végeselemhálót hanem a szerkezeti rajzokat is. A rajzok elrejtése: B(1,2) Display Filter [ ]Wireframe [ ]Parts [ ]Assembly [ ]Work Plane [ ]Model views Border [ ] FE Model
az elrejtést az objektumok előtti megjelölés megszüntetésével érhetjük el. A szemléltetendő eredmények kiválasztása: Elmozdulás abszolút értéke A(1,1) Results [1-B.C. 1,DISPLACEMENT…] [Magnitude] [1-B.C. 1,DISPLACEMENT…] A(2,1) Display Közép Billentyű
Redukált feszültség A(1,1) Results [3-B.C. 1,STRESS…] [von Mises] [1-B.C. 1,DISPLACEMENT…] A(2,1) Display Közép Billentyű
A számítást célszerű megismételni sűrűbb felosztással is, azaz átlagosan 7.5 mm, valamint 3.75 mm méretű háromszög alakú elemekkel is. A vizsgálat célja, hogy meggyőződjünk arról, hogy a kapott feszültségcsúcsok vajon tartanak-e egy képzeletbeli egzakt határértékhez?
7.4. ábra. Az elmozdulás abszolút értéke mm-ben
103
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
7.5. ábra. A redukált feszültség eloszlása MPa-ban
104
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
8. TÉRBELI FELADAT ELEMZÉSE Ebben a fejezetben térbeli feladatok modellezésére mutatunk be példákat. Megismerkedhetünk az I-DEAS programrendszerben alkalmazható térbeli elemtípusokkal és a térbeli hálógenerálással. Egy fali konzolt fogunk modellezni, amelyet két csavarral rögzítünk a falhoz. Feltételezzük, hogy két konzol felső felületére fektetünk egy polcot, amit 500 N erővel terhelünk. A konzolok és a terhelés szimmetrikus elrendezése esetén mindkét konzol felső felületére 250 N eredőjű egyenletesen megoszló terhelés jut. A konzolok és a polc helyett csak egyetlen konzolt modellezünk: előállítjuk a geometriát, megadjuk a peremfeltételeket, előállítjuk a térbeli végeselemes hálót, lefuttatjuk a számítást, majd kiértékeljük az eredményeket.
8.1. A geometria előállítása A fali konzol előállításához az I-DEAS programrendszer Simulation moduljának Master Modeler alkalmazását válasszuk ki. Az Options/Units paranccsal állítsuk be az alkalmazott mértékegységeket milliméterre és Newtonra. 8.1. ábra A konzol szelvényének kontúrját a Polylines paranccsal hozzuk Master Modeler létre (A(2, 1) ikongyűjtő). Töröljük a felesleges méreteket a B(4, 1) Delete parancs segítségével, majd állítsuk elő a mérethálót (A(4, 1) Dimension). Állítsuk be a kívánt méreteket a B(2, 1) Modify Entity paranccsal. A modellt a CTRL-A billentyűzetkombinációval nagyíthatjuk látható méretűre. A lekerekítéseket az A(4, 2) Fillet paranccsal (Trim/Extend opció bekapcsolva) hozhatjuk létre. Extrudáljuk a kontúrt 30 mm-es vastagságúra (A(5, 1) Extrude).
8.2. ábra. Fali konzol szelvénye A következőkben létrehozzuk a furatokat a konzol falhoz csatlakozó felületén. Ehhez jelöljük ki a kérdéses felületet szerkesztési síknak (A(1, 1) Sketch in place). Rajzoljunk két kört a furatoknak (A(3, 1) Center Edge), majd méretezzük be (Dimension, Modify Entity). Ezután „kifúrjuk” a furatokat (Extrude (Cut opció)).
105
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
8.3. ábra. A furatok kivágása
8.2. Peremfeltételek megadása A peremfeltételek megadásához válasszuk a Boundary Conditions alkalmazást.
8.4. ábra. Boundary Conditions Kinematikai peremfeltételként kössük meg a furatok belső felületének szabadságfokait (A(4, 2) Displacement Restraint). A dinamikai peremfeltétel a konzol felső felületén előírt normál irányú megoszló terhelés (A(2, 2) Pressure, Total Force: 250 N opció).
8.5. ábra. Kinematikai peremfeltétel
106
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
8.6. ábra. Dinamikai peremfeltétel
8.7. ábra. Peremfeltételek
8.3. Végeselemes háló generálása A hálógenerálást a Meshing alkalmazásban találjuk. Térbeli feladatokhoz négycsomópontú (lineáris) vagy tízcsomópontú (kvadratikus) tetraéder elemet használhatunk. Az utóbbi alkalmazása kedvezőbb, mert megbízhatóbb megoldást kapunk.
8.8. ábra. Meshing Generáljunk hálót először 15 mm-es átlagos elemméretű négycsomópontú tetraéder elemekkel (A(1, 1) gyűjtő, Define Solid Mesh). Anyagként általános izotrop acélt válasszunk (ez az alapbeállítás).
8.9. ábra Térbeli elemtípusok 8.10. ábra. Hálógenerálás
107
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
8.11. ábra. Végeselemes háló
8.4. Feladat megoldása, eredmények kiértékelése A Model Solution alkalmazásban tudjuk a kitűzött feladatot megoldani és a kapott eredményeket kiértékelni.
8.12. ábra. Model Solution Hozzunk létre egy új megoldáshalmazt az eredményeknek (A(1, 2) Solution Set, Create). Itt állíthatjuk be a megoldás módjára vonatkozó paramétereket is (direkt vagy iteratív megoldás, stb.). A feladat megoldását az A(2, 1) Manage Solve paranccsal indíthatjuk. A generált végeselemes háló a nagyméretű lineáris elemek miatt rosszul közelíti a görbült felületeket, és a megoldás pontosságára is hatással van. Ezen elemek használata ezért is kerülendő. Erre utal a figyelmeztető üzenet.
8.13. ábra. A négycsomópontú térbeli elemek használata kerülendő Az eredényeket az A(6, 2) Visualizer eszközzel tekinthetjük meg. A feszültség egyenetlen képet mutat, ami a durva felosztásnak és a lineáris elemek használatának tudható be.
108
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
8.14. ábra. Lineáris elemek használatából adódó egyenetlen feszültség
8.5. A háló finomítása A továbbiakban új, sűrűbb hálót generálunk tízcsomópontú kvadratikus elemekkel. Először tehát töröljük a legutóbbi megoldást (Solution Set, Solution Set1 kijelölése, Delete). A Meshing alkalmazás B(1, 3) Display Filter parancsával a megjelenítés tulajdonságait úgy állítjuk be, hogy csak az elemek és a csomópontok látszódjanak (Parts elől a pipa törlése, FE Models alatt csak a Node és Element legyen kijelölve).
8.15. ábra. Az entitások megjelenítésének tulajdonságai Ezután jelöljünk ki egy elemet, majd a jobbgombos menüből válasszuk az All menüpontot, ezzel kijelölve az összes elemet. Ezután a B(4, 1) Delete FE Entities parancssal töröljük az elemeket (a műveletet meg kell erősíteni). Hasonlóan, csak csomópontot kijelölve töröljük az összes csomópontot is. Ezzel töröltük a végeselemes hálót, de megtartottuk a geometriához kötött (de nem a csomópontokhoz és elemekhez kötött) peremfeltételeket. A Display Filterben állítsuk vissza a megjelenítés tulajdonságait. Most a Define Solid Mesh paranccsal 10 mm-es átlagos elemmérettel, kvadratikus elemeket használva hozzuk létre az új hálót. Az előzőekben leírt módon, a Model Solution alkalmazásban oldjuk meg a feladatot, és indítsuk el a Visualizert. A kapott feszültségeloszlásból (Mises-féle redukált feszültség) látható, hogy a felső furat pereme a feszültséggyűjtő hely (kb. 50 MPa a maximális feszültség). Érdekes még a lekerekített hajlat,
109
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
ahol megközelítőleg 30 MPa a feszültség. A maximális elmozdulás 0.569 mm (a konzol végén).
8.16. ábra. Kvadratikus elemek használatával kapott feszültségeloszlás
8.17. ábra. Elmozdulásmező A következő fejezetben megváltoztatjuk a konzol geometriáját, hogy kedvezőbb feszültségeloszlást kapjunk. A felosztást is sűrítjük a pontosabb megoldás reményében.
110
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
9. TÉRBELI FELADAT ELEMZÉSE Ezúttal az előző fejezetben modellezett fali konzolban ébredő feszültségeloszlást próbáljuk megváltoztatni egy nagyobb sugarú lekerekítés alkalmazásával.
9.1. A geometria módosítása Az előző fejezetben létrehozott modellből először töröljük a végeselemes modellt. Ehhez a Master Modeler alkalmazásban válasszuk a B(4, 2) Manage Bins parancsot. A felnyíló párbeszédablak View opcióját állítsuk FE Modelre. Jelöljük ki a korábban létrehozott végeselemes modellt, majd töröljük a Delete parancssal.
9.1. ábra. Manage Bins párbeszédablak A lekerekítés módosításához használjuk a B(1, 1) History Access parancsot, és jelöljük ki az alkatrészt. Ennek hatására egy fastruktúrában láthatjuk az alkatrész „történetét”, vagyis a geometria kialakításához felhasznált műveleteket.
9.2. ábra. Az alkatrész „története”
111
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
A konzol szelvényének extrudálását a 9.2. ábrán látható módon kijelölve nyomjuk meg a Wireframe gombot. Ennek hatására szerkeszthetővé válik a szelvény. A B(2, 1) Modify Entity parancssal módosítsuk a konzol „hajlatának” lekerekítési sugarát 50 mm-re.
9.3. ábra. A módosított lekerekítési sugarú szelvény A módosítások érvényre juttatása a B(3, 2) Complete Update parancssal, és az alkatrész kijelölésével történik. A furatok peremeinél kialakuló feszültséggyűjtő helyeknél alkalmazzunk 1,5 mm-es élletörést, a peremek kijelölésével és az A(5, 2) Chamfer (élletörés) parancs kiadásával. A módosított geometria 9.4. ábrán látható.
9.4. ábra. Az alkatrész a módosított lekerekítéssel és az élletörésekkel
9.2. Peremfeltételek A kinematikai és dinamikai peremfeltételek megegyeznek az előző fejezetben alkalmazottakkal, vagyis a furatok belső felületén nem engedünk meg elmozdulást, a konzol felső felületén pedig 250 N eredőjű megoszló terhelést írunk elő. Az I-DEAS programrendszerben a peremfeltételek megadási módját az előző fejezet példájánál ismertettük (Boundary Conditions alkalmazás).
112
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
9.5. ábra. A fali konzol peremfeltételei
9.3. Végeselemes háló Az előző fejezetben ismertetett módon, 10 mm átlagos elemméretű végeselemes hálót generálunk tízcsomópontú térbeli kvadratikus elemekkel (Meshing alkalmazás).
9.6. ábra. Végeselemes háló generálása kvadratikus térbeli elemekkel
9.4. Eredmények A számítás elvégzése után (Solution Set, Manage Solve (Model Solution alkalmazás)) a 9.7. ábrán látható feszültségeloszlást kapjuk (Visualizer).
113
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
9.7. ábra. A módosított alkatrész feszültségeloszlása A Mises-féle redukált feszültség eloszlása a „hajlatban” egyenletesebb lett, valamint a feszültség értéke is alacsonyabb (kb. 10 MPa). A felső furat továbbra is feszültséggyűjtő hely, de az eloszlás egyenletesebb. A maximális elmozdulás 0,378 mm-re csökkent.
9.8. ábra. Elmozdulás a módosított modellen A számítást 5 mm-es átlagos elemmérettel megismételve a következő feszültségeloszlást és elmozdulást kapjuk. A feszültséggyűjtő helyen a maximális redukált feszültség tovább nő, de a többi eredmény megbízhatóbbnak vehető.
114
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
9.9. ábra. Feszültségeloszlás 5 mm-es átlagos elemméretű háló esetén
9.10. ábra. A sűrített hálóval kapott elmozdulások
115
A VÉGESELEMES MODELLEZÉS KONTINUUMMECHANIKAI ALAPJAI
10. IRODALOMJEGYZÉK [1] PÁCZELT I.: Végeselem-módszer a mérnöki gyakorlatban, I. kötet, Miskolci Egyetemi Kiadó, Miskolc, 1999. [2] SZABÓ, B. BABUSKA,I.: Finite Element Analysis, John Wiley & Sons Inc., New York, 1991. [3] BATHE, K.J.: Finite Element Procedures, Prentice-Hall, Inc., New Jersey, 1996. [4] HUGHES, T.J.R.: The Finite Element Method Linear Static and Dynamic Finite Element Analysis, Dover Publications, Inc., Mineola, New York, 2000. [5] I-DEAS Online Program dokumentáció, http://www.i-deas.hu
117