Content extract
http://www.doksihu Eötvös Loránd Tudományegyetem Természettudományi Kar Peremfeltételek hatása a végeselem módszer pontosságára Szakdolgozat Írta: Szabó András Matematika BSc szak Témavezet®k: Simon Péter, egyetemi docens Alkalmazott Analízis és Számításmatematikai Tanszék Eötvös Loránd Tudományegyetem, Természettudományi Kar Horváth Tamás, egyetemi tanársegéd Matematika és Számítástudomány Tanszék Széchenyi István Egyetem, M¶szaki Tudományi Kar Budapest, 2009 http://www.doksihu Tartalomjegyzék 1. Bevezetés 1 2. Egydimenziós végeselem 2 2.1 Homogén Dirichlet peremfeltétel 2.11 Diszkrét probléma 3 . 6 2.12 Megvalósítás 8 2.2 További peremfeltételek 9 2.21 Inhomogén Dirichlet peremfeltétel 9 2.22 Neumann peremfeltétel 11 2.23 Harmadfajú peremfeltétel 12
3. Magasabb dimenziós végeselem 13 3.1 Homogén Dirichlet peremfeltétel 14 3.11 Megvalósítás 16 3.2 További peremfeltételek 18 3.21 Inhomogén Dirichlet peremfeltétel 18 3.22 Neumann peremfeltétel 19 3.23 Harmadfajú peremfeltétel 19 4. Peremfeltételek vizsgálata 21 4.1 Hibaszámítás 22 4.2 L-alakú tartomány 22 4.3 Vizsgálat 23 4.31 Végtelen szingularitású függvény vizsgálata 24 4.32 Oszcilláló függvény vizsgálata 25 ii http://www.doksihu 4.33 Hirtelen ugrású függvény vizsgálata 27 4.4 Konklúzió 30 Köszönetnyilvánítás 31 Irodalomjegyzék 32 iii http://www.doksihu 1. fejezet Bevezetés Jelen dolgozat célja az, hogy megismertesse az
Olvasót a végeselem módszerrel, amely az egyik leggyakrabban használt eljárás parciális dierenciálegyenletek numerikus megoldására. Ezentúl különböz® peremfeltételek esetén megvizsgáljuk a módszer pontosságát a hiba tesztelés alapjait bemutatva. A tartományt illet®en Solin, Cerveny és Dolezel cikke ([6]) motiválja, hogy a számításokat a (−1, 1)2 (−1, 0)2 úgynevezett L-alakú tartományon végezzük. A vizsgálat során három különböz® típusú megoldásra fogjuk számolni a hibákat, nevezetesen amikor a megoldásnak a sarokpont körül végtelenbeli határértéke van, a megoldás oszcilláló függvény és végül amikor hirtelen ugrása van. Az els® fejezetben tárgyaljuk az egydimenziós esetet egy példafeladaton keresztül, azaz amikor egy másodrend¶ közönséges dierenciálegyenletnek keressük a közelít® megoldását. Ezek után a második fejezetben bemutatjuk a módszert elliptikus parciális dierenciálegyenletek esetén Az
elliptikus egyenleteket analitikusan csak néhány speciális esetben tudjuk megoldani, ezért kézenfekv® a végeselem módszer alkalmazása. Az els® két fejezetben el®ször homogén Dirichlet peremfeltételt csatolunk az egyenletünkhöz, majd a többi peremfeltételnél fellép® módosításokat külön-külön tárgyaljuk. A harmadik fejezetben a már említett L-alakú tartományon vizsgáljuk az egyenletet oly módon, hogy a pontos megoldás ismeretében viszonyítjuk a számolt megoldást különböz® peremfeltételek esetén. 1 http://www.doksihu 2. fejezet Egydimenziós végeselem Jelen fejezetben tárgyaljuk a végeselem módszert másodrend¶ közönséges dierenciálegyenletek esetén. El®ször homogén Dirichlet peremfeltétel mellett mutatjuk be az elméleti alapokat, majd a módszer gyakorlati megvalósítása után részletezzük a többi peremfeltételt. Legyen p ≥ p0 > 0, q ≥ q0 > 0, p, q ∈ L∞ (I) és tekintsük az I = (0, 1) intervallumon a
következ® másodrend¶ közönséges dierenciálegyenletet: −(pu0 )0 + qu = f Az egyenletnek kell adnunk valamiféle peremfeltételt, amely három f® csoportból kerülhet ki: Dirichlet peremfeltétel: u(0) = ξ , u(1) = η . Neumann peremfeltétel: −p(0)u0 (0) = ξ , p(1)u0 (1) = η . Harmadfajú peremfeltétel: −u0 (0) + σ0 u(0) = c0 , u0 (1) + σ1 u(1) = c1 , ahol σ0 -ról és σ1 -r®l kés®bb teszünk megkötéseket. Ha ξ = η = 0 akkor homogén, különben inhomogén peremfeltételr®l beszélünk. 2 http://www.doksihu 2.1 Homogén Dirichlet peremfeltétel Tekintsük a következ® egyenletet és a hozzá tartozó peremfeltételt: −(pu0 )0 + qu = f (2.1) u(0) = u(1) = 0 A fenti feladat klasszikus megoldásán egy u ∈ C 2 (I) függvényt értünk, de ez sokszor túl er®s feltétel az u-ra nézve, ezért szeretnénk egy b®vebb függvénytérben keresni a megoldást. Legyen v ∈ C 2 (I) olyan függvény, melyre v(0) = v(1) = 0 és szorozzuk meg v -vel a
(2.1) egyenletet, valamint integráljunk 0-tól 1-ig Ekkor Z 1 Z 1 0 0 f v ∀v ∈ C 2 (I) pu v + quv = {z } | 0{z } |0 a(u,v) (2.2) ϕ(v) hiszen, ha −(pu0 )0 v -t parciálisan integráljuk, akkor Z 1 Z 1 h i1 Z 1 0 0 0 0 0 (−pu )v = pu0 v 0 −(pu ) v = − pu v − 0 0 0 0 | {z } 0 A (2.2) egyenletben azonban már csak u els® deriváltja szerepel, vagyis kereshetjük a megoldást egy b®vebb térben. Felvet®dik a kérdés, hogy a továbbiakban milyen függvénytérben keressük a feladat megoldását, ezért deniáljuk a H 1 (I) és H01 (I) tereket. Legyen H 1 (I) := {u : u ∈ L2 (I), u0 ∈ L2 (I)}, ahol a deriválás disztribúciós értelemben értend® ([3]). 2.11 Állítás A H 1 (I) tér elemei folytonosak Ezek után deniáljuk a H01 (I) teret, legyen H01 (I) := {u : u ∈ H 1 (I), u(0) = u(1) = 0} amely értelmes, hiszen H 1 (I) elemei folytonosak. A peremfeltételt gyelembe véve az u megoldást a C 2 (I) tér helyett mostantól a H01 (I) térben keressük,
ez a gondolatmenet motiválja a következ® deníciót: 2.11 Deníció Egy u ∈ H01 (I) függvényt a (21) peremértékfeladat gyenge megoldásának nevezzük, ha Z Z 1 1 0 0 (pu v + quv) = 0 fv 0 teljesül minden v ∈ H01 (I) esetén. 3 http://www.doksihu 2.11 Megjegyzés Ha u egy klasszikus megoldás, akkor egyben gyenge megoldás is, fordítva is majdnem igaz, hiszen ha u egy C 2 -beli gyenge megoldás, akkor klasszikus megoldás is. A továbbiakban szeretnénk megfogalmazni egy tételt a megoldás létezését és unicitását illet®en, ehhez azonban ki kell mondjunk néhány segédállítást ([3]). H 1 (I) egy Hilbert tér, azaz teljes a h·, ·iH 1 skalárszorzatból indukált k·kH 1 normával, ahol Z 1 hu, viH 1 (I) = (u0 v 0 + uv) 0 és kukH 1 (I) s Z q = hu, uiH 1 (I) = 1 (|u0 |2 + |u|2 ) 0 2.12 Állítás A H01 (I) tér teljes a h·, ·iH 1 (I) skalárszorzattal Deniáljuk a H01 (I) téren a h·, ·iH01 (I) skalárszorzatot és a k · kH01 (I)
normát a következ®képpen: Z hu, viH01 (I) := valamint kukH01 (I) 1 (u0 v 0 ) 0 p := hu, ui = s Z 1 |u0 |2 0 1. Tétel (Sztyeklov egyenl®tlenség) Létezik c > 0, hogy minden u ∈ H01 ([a, b]) függvényre Z Z b b 2 |u| ≤ c a |u0 |2 . a 2.12 Megjegyzés Visszafelé nem igaz, azaz k·kL2 norma nem ekvivalens k·kH01 (I) normával. 2.13 Állítás k · kH 1 norma ekvivalens k · kH01 (I) normával a H01 (I) térben Bizonyítás Z kukH01 (I) = 1 Z 1 0 2 |u | ≤ 0 Z 0 2 1 2 (|u | + |u| ) = kukH 1 ≤ (1 + c) 0 0 4 |u0 |2 = (1 + c)kukH01 (I) http://www.doksihu A továbbiakban általánosítjuk a problémát, azaz nem konkrétan a H01 (I) térben keressük a megoldást, hanem egy tetsz®leges V térben, amir®l csak annyit követelünk meg, hogy Hilbert tér legyen. 2.12 Deníció (Folytonosság) Legyen V egy normált vektortér az a(·, ·) bilineáris forma folytonos, ha a(u, v) ≤ CkukV kvkV teljesül minden u, v ∈ V esetén, ahol C
konstans. 2.13 Deníció (Elliptikusság) Legyen V egy normált vektortér az a(·, ·) bilineáris forma V -elliptikus, ha a(u, u) ≥ ckuk2V teljesül minden u ∈ V esetén, ahol c konstans. Legyen adott tehát egy V Hilbert tér, ami teljes a k · kV normával, valamint legyen a(·, ·) : V × V R szimmetrikus, bilineáris függvény, mely folytonos és V elliptikus. Ekkor h·, ·ia := a(·, ·) skalárszorzat V -ben és a bel®le származtott k · ka norma ekvivalens az eredeti k · kV normával. Továbbá legyen adott egy ϕ : V R folytonos lineáris funkcionál. 2.14 Deníció (Absztrakt variációs feladat) Keresend® u ∈ V , melyre a(u, v) = ϕ(v) ∀v ∈ V 2. Tétel (Egzisztencia és unicitás) Az absztrakt variációs feladatnak minden ϕ(·) folytonos lineáris funkcionál esetén létezik egyértelm¶ megoldása. Bizonyítás Mivel V Hilbert tér és ϕ(·) mind k · ka mind k · kV normában folytonos lineáris funkcionál ezért a Riesz reprezentációs
tétel szerint egyértelm¶en létezik egy u ∈ V függvény melyre a(v, u) = ϕ(v), mivel a(·, ·) szimmetrikus ezért a(u, v) = ϕ(v) 3. Tétel A (21) peremértékfeladatnak mindig létezik egyértelm¶ megoldása 5 http://www.doksihu R1 Bizonyítás Legyen most V = H01 (I), belátjuk, hogy ϕ(v) = funkcionál: Z s Z 1 |ϕ(v)| = | |v|2 = kf kL2 kvkL2 0 H ölder f v folytonos lineáris 1 |f |2 f v| ≤ |{z} 0 s Z 1 0 0 ∀v ∈ H01 (I) ≤ (kkf kL2 )kvkH01 (I) = M kvkH01 (I) ugyanis kvkL2 ≤ kvkH 1 ≤ kkvkH01 (I) továbbá ϕ(·) linearitása triviális. R1 Nyilván a(u, v) = 0 pu0 v 0 + quv lineáris és szimmetrikus, belátjuk, hogy folytonos és elliptikus: Folytonosság: Z Z 1 pu v + quv| ≤ |a(u, v)| = | ³Z ³ Z 1 ≤M s ³ Z s s 1 0 Z 1 |u0 |2 ≤M( 0 0 |v|2 ≤ |{z} Sztyeklov Z 1 0 ´ 1 |u0 |2 D C + 0 ´ 1 0 Z 1 |v 0 |2 0 Z + 0 |u0 |2 ≤M Z 1 s Z 1 |uv| ≤ 0 |u|2 0 ´ 1 0 0 s Z |v 0 |2 0 Z 1
|u v | + 1 |u0 |2 0 0 M Z |quv| ≤ 0 ≤ max (kpkL∞ , kqkL∞ ) {z } | 1 0 0 |pu v | + 0 s Z 1 0 0 |v 0 |2 ≤ 0 |v 0 |2 ) = M 0 kukH01 (I) kvkH01 (I) ∀u, v ∈ H01 (I) Elliptikusság: Z Z 1 |a(u, u)| = p|u | + 0 Z 1 0 2 1 2 q|u| ≥ 0 0 p0 |u0 |2 = p0 kuk2H 1 (I) 0 ∀u ∈ H01 (I) Teljesülnek a 2. tétel feltételei, tehát a (21) peremértékfeladatnak mindig létezik egyértelm¶ megoldása. 2.11 Diszkrét probléma A továbbiakban (2.1) peremértékfeladat egy közelít® megoldását keressük, ehhez deniálunk egy véges dimenziós teret és ebben keressük az absztrakt variációs feladat 6 http://www.doksihu megoldását. Legyen Vh ⊂ V véges dimenziós altér. Keresünk egy uh ∈ Vh függvényt, melyre (2.3) a(uh , vh ) = ϕ(vh ) ∀vh ∈ Vh . Nyilván a Vh térre megszorított a(·, ·) és ϕ(·) továbbra is eleget tesz a 2. tételnek, vagyis a diszkrét problémának mindig létezik egyértelm¶ megoldása. 2.14
Állítás u − uh ⊥ Vh az a(·, ·) skalárszorzat szerint, azaz: a(u − uh , vh ) = 0 ∀vh ∈ Vh Bizonyítás Tudjuk, hogy a(u, vh ) = ϕ(vh ) ∀vh ∈ Vh ⊂ V , ebb®l ha kivonjuk a (2.3) egyenletet, akkor azt kapjuk, hogy a(u − uh , vh ) = 0 ∀vh ∈ Vh 2.15 Állítás uh az u legjobb Vh -beli közelítése az k · ka norma szerint Bizonyítás Legyen vh ∈ Vh , ekkor ku − vh k2a = ku − uh + uh − vh k2a = k u − uh k2a + k uh − vh k2a ≥ ku − uh k2a | {z } | {z } ∈Vh⊥ ∈Vh 2.11 Lemma (Céa lemma) uh kvázioptimális közelítése u-nak k·kV norma szerint, azaz ku − uh kV ≤ C inf ku − vh kV vh ∈Vh ahol C = M . m Bizonyítás mku − uh k2V ≤ |a(u − uh , u − uh )| = |a(u − uh , u − vh + vh − uh )| = |a(u − uh , u − vh ) + a(u − uh , vh − uh ) | = |a(u − uh , u − vh )| = M ku − uh kV ku − vh kV | {z } | {z } ∈V ∈V ⊥ {z } | 0 7 http://www.doksihu 2.12 Megvalósítás A továbbiakban szeretnénk
valamilyen számítógépen programozható eljárást adni, amely segítségével könnyen megtalálhatjuk az uh közelít® megoldást. Tekintsük a [0, 1] intervallum egy egyenletes felosztását: 0 = x0 < x1 < . < xn−1 < xn = 1, h = xi+1 − xi Legyenek a Vh teret generáló bázisfüggvények: Φ1 (x), Φ2 (x), . , Φn−1 (x), ahol x−xi−1 x ∈ [xi−1 , xi ] xi −xi−1 Φi (x) := xi+1 −x x ∈ [xi , xi+1 ] xi+1 −xi 0 különben. A vizsgálataink során mi csak els®fokú bázisfüggvényekkel fogunk dolgozni, de megjegyezzük, hogy a fokszám növelésével sok esetben pontosabb megoldást kaphatunk, ezt hívják p-végeselemnek, ahol a p a bázisfüggvények fokszámára utal. A fent deniált Φi (x) függvényeket szokták kalapfüggvényeknek is nevezni. Tehát Vh = span {Φ1 (x), Φ2 (x), . , Φn−1 (x)}, továbbá keressük a közelít® megoldást 1 0.8 0.6 0.4 0.2 0 0 0.2 0.4 0.6 0.8
1 2.1 ábra Φi bázisfüggvénynek uh (x) = n−1 X αi Φi (x) i=1 alakban, vagyis a feladatunk meghatározni az αi együtthatókat. Pn−1 P Ekkor ϕ(Φj ) = a(uh , Φj ) = a( n−1 i=1 αi a(Φi , Φj ), j = 1, . , n − 1 i=1 αi Φi , Φj ) = vagyis az αi együtthatókra egy Aα = F lineáris egyenletrendszer adódik, ahol Ai,j = a(Φi , Φj ) (i, j = 1, . , n − 1), α = (α1 , , αn−1 )T , F = (ϕ(Φ1 ), , ϕ(Φn−1 ))T 8 http://www.doksihu Könny¶ meggondolni, hogy a(Φi , Φj ) = 0, ha |i − j| > 1 így A tridiagonális mátrix, ezért az egyenletrendszer O(n) id®ben megoldható ([1]). A módszer konvergenciájára vonatkozóan igazolható az alábbi tétel ([2]): 4. Tétel Legyen u a variációs feladat megoldása, uh pedig a diszkrét feladat megoldása Ha u00 létezik, akkor fennáll a következ® egyenl®tlenség: ku−uh kL2 (I) ≤ Ch2 ku00 kL2 (I) . 2.2 További peremfeltételek Más peremfeltétel esetén a módszer rendkívül
hasonló az el®z® fejezetben tárgyaltakhoz, mind az elméletet, mind a megvalósítást illet®en. Azonban jelen dolgozat f® célja különböz® peremfeltételek esetén bemutatni a módszer pontosságát, ezért röviden bemutatjuk a f® különbségeket a homogén Dirichlet peremfeltételhez képest. 2.21 Inhomogén Dirichlet peremfeltétel Tekintsük a következ® egyenletet és a hozzá tartozó peremfeltételt az I = (0, 1) intervallumon: −(pu0 )0 + qu = f (2.4) u(0) = ξ, u(1) = η Legyen u0 egy olyan H 1 (I) beli függvény, melyre u0 (0) = ξ és u0 (1) = η és legyen u = w + u0 ahol w ∈ H01 (I). Szorozzuk meg egy v ∈ H01 (I) függvénnyel a (2.4) egyenletet és integráljunk 0-tól 1-ig Ekkor Z 1 Z (−(pu ) v + quv) = 0 Z 1 0 0 1 0 0 (−(pw ) v + qwv) + 0 0 Z (−(pu00 )0 v 1 + qu0 v) = fv 0 Kihasználva, hogy v ∈ H01 (I) adódik Z 1 h i1 Z 1 0 0 0 pw0 v 0 −(pw ) v = − pw v + 0 0 0 | {z } 0 Hasonlóan w helyett u0 ra is igaz a fenti
egyenl®ség, azaz egyenletünk a következ®: Z 1 Z 1 Z 1 0 0 0 0 f v ∀v ∈ H01 (I) (−pw v + qwv) + (−pu0 v + qu0 v) = 0 0 0 9 http://www.doksihu a korábbi jelöléseket használva: a(w, v) + a(u0 , v) = ϕ(v) ∀v ∈ H01 (I) átrendezve a(w, v) = ϕ(v) − a(u0 , v) ∀v ∈ H01 (I) mivel u0 ismert, ezért a feladatunk megtalálni a w megoldást. Könny¶ meggondolni, g = ϕ(v)−a(u0 , v) továbbra is folytonos lineáris funkcionál, ezért a feladatot hogy ϕ(v) visszavezettük a homogén Dirichlet peremfeltétel esetére, mivel ha megtaláltuk a w megoldást, akkor az eredeti megoldást u = w + u0 összeg adja meg, ahol u0 ismert. 2.21 Állítás A megoldás független u0 megválasztásától Most megvizsgáljuk a megvalósítást, azaz, hogy miként tudjuk algoritmikusan megközelíteni a feladatot. A 21 szakaszban tárgyalt els®fokú bázisfüggvényeket itt is jól tudjuk használni, mindössze annyi a feladatunk, hogy deniáljuk a 0 és az 1 ponthoz
tartozó Φ0 (x) illetve Φn (x) bázisokat a következ®képpen: x1 −x x ∈ [x0 , x1 ] x1 Φ0 (x) := 0 különben x−xn−1 Φn (x) := 1−xn−1 x ∈ [xn−1 , xn ] 0 különben 1 0.8 0.6 0.4 0.2 0 0 0.2 0.4 0.6 0.8 1 2.2 ábra Φ0 és Φn bázisfüggvény Legyen u0 (x) = ξΦ0 (x) + ηΦn (x), valamint Vh = span {Φ0 , Φ1 , . , Φn } és keP ressük a w közelít® megoldását wh (x) = n−1 i=1 αi Φi (x) alakban. 10 http://www.doksihu P Pn−1 Ekkor ϕ(Φj ) = a(wh , Φj ) = a( n−1 i=1 αi Φi , Φj ) = i=1 αi a(Φi , Φj ), (j = 0, . , n), tehát az αi együtthatókra ismét egy Aα = F lineáris egyenletrendszert kaptunk, ahol A továbbra is tridiagonális mátrix, ezért az algoritmus O(n) id®ben megoldható. Ebb®l u közelít® megoldását úgy kapjuk meg, hogy hozzádajuk az u0 függvényt, azaz uh = wh + u0 . 2.22 Neumann peremfeltétel Legyen most egyenletünk és a hozzátartozó peremfeltétel a
következ®: −(pu0 )0 + qu = f (2.5) p(0)u0 (0) = ξ, p(1)u0 (1) = 0 ahol q > 0. Továbbra is az I = (0, 1) intervallumon fogjuk vizsgálni a feladatot Tekintsük a következ® függvényteret: V = {v ∈ H 1 (I) : v(1) = 0}. Szorozzuk egy v ∈ V függvénnyel a (2.5) egyenletet, majd integráljunk 0-tól 1-ig Ekkor Z 1 h i Z 0 0 0 −(pu ) v = − pu v + 0 1 Z 0 0 pu v = p(0)u (0)v(0) + 0 Z 1 0 1 0 0 pu v = −ξv(0) + 0 pu0 v 0 0 Továbbá −ξv(0) = −ξδ0 (v), ahol δ0 a 0-hoz tartozó Dirac disztribúció. Ezt felhasználva és az egyenletet átrendezve: Z 1 Z 1 0 0 f v + ξδ0 (v) pu v + quv = {z } |0 {z } |0 a(u,v) ϕ(v) Belátható, hogy ϕ(·) ez esetben is folytonos lineáris funkcionál, továbbá a(·, ·) folytonos és V-elliptikus, vagyis a (2.5) feladatnak is létezik egyértelm¶ megoldása 11 http://www.doksihu 2.23 Harmadfajú peremfeltétel −(pu0 )0 + qu = f (2.6) −p(0)u(0)0 + σ0 u(0) = c0 p(1)u(1)0 + σ1 u(1) = c1 ahol σ0 ,
σ1 > 0. Szorozzunk egy v ∈ H 1 (I) függvénnyel és integráljunk az I intervallumon, ekkor: Z 1 Z 1 h i1 Z 1 0 0 0 0 0 0 0 −(pu ) v = − pu v + pu v = −p(1)u(1) v(1) + p(0)u(0) v(0) + pu0 v 0 = 0 0 0 0 Z 1 = (σ1 u(1) − c1 )v(1) + (σ0 u(0) − c0 )v(0) + pu0 v 0 ∀v ∈ H 1 (I) 0 Ezt felhasználva és átrendezve adódik, hogy Z 1 Z 1 0 0 f v − (σ1 u(1) − c1 )v(1) − (σ0 u(0) − c0 )v(0) pu v + quv = 0 0 {z } | {z } | a(u,v) ∀v ∈ H 1 (I) ϕ(v) Ismét belátható, hogy a(·, ·) folytonos és V -elliptikus, valamint ϕ(·) folytonos lineáris funkcionál. 12 http://www.doksihu 3. fejezet Magasabb dimenziós végeselem Ebben a fejezetben tárgyaljuk a végeselem módszer megvalósítását elliptikus parciális dierenciálegyenletek esetén. Az el®z® fejezet felépítését követve el®ször a homogén Dirichlet esetre mutatjuk be az elméleti alapokat, majd a megvalósítást követ®en részletezzük a többi esetet. Célszer¶
el®ször egy általánosabb függvényegyütthatós lineáris egyenletet vizsgálni, melynek speciális esetét véve kapjuk az elliptikus egyenletek kanonikus alakját. Legyen p ≥ p0 > 0, q ≥ 0, p, q ∈ L∞ (Ω) és tekintsük a következ® parciális dierenciál egyenletet az Ω ⊆ Rn tartományon: −div (p∇u) + qu = f Az egyenlethez csatolnunk kell valamilyen peremfeltételt, melyek három f® csoportból kerülhetnek ki: Dirichlet peremfeltétel: u|∂Ω = u0 . Neumann peremfeltétel: p ∂u | = a. ∂η ∂Ω Harmadfajú peremfeltétel: p ∂u | + σu|∂Ω = c. σ ≥ σ0 ≥ 0 ∂η ∂Ω ahol η az Ω tartomány kifelé mutató normálisa és ∂u ∂η a normális irányú iránymenti derivált. Megjegyezzük, hogy ha p ≡ 1 és q ≡ c, akkor a következ® alakra egyszer¶södik az egyenletünk: −4u + cu = f 13 http://www.doksihu vagyis a vizsgálni kívánt kanonikus alakot kapjuk. 3.1 Homogén Dirichlet peremfeltétel Tekintsük a
következ® peremérték feladatot −div (p∇u) + qu = f (3.1) u|∂Ω = 0 A feladat klasszikus megoldásán, egy u ∈ C 2 (Ω)∩C(Ω) függvényt értünk. Az egy dimenziós esethez hasonlóan gyengíteni szeretnénk a megoldásra vonatkozó feltételt, azaz egy b®vebb függvénytérben szeretnénk keresni u-t. Ehhez szükségünk lesz a Green tétel alábbi változatára: 5. Tétel (1 Green tétel) Tegyük fel, hogy u ∈ C 2 (Ω) és v ∈ C 1 (Ω) Ekkor Z Z Ω [div (p∇u) + qu]v = − Z Z p∇u∇v + Ω quv + Ω p ∂Ω ∂u vdσ. ∂η Legyen v ∈ C 1 (Ω)-beli függvény melyre v|∂Ω = 0, szorozzuk meg v -vel a (3.1) egyenletet és integráljunk Ω-n. Ekkor a Green tétel felhasználásával adódik, hogy Z Z Z Z ∂u fv p∇u∇v + quv − p vdσ = Ω Ω Ω ∂Ω ∂η {z } | 0 Megjegyezzük, hogy a Green tétel ezen alakja tekinthet® a parciális integrálás magasabb dimenziós változatának, vagyis a fenti egyenl®ség az egy
dimenziós esetben tárgyalt átalakítás. Mivel v az Ω határán elt¶nik, ezért a bal oldalon szerepl® felületi integrál értéke nulla, azaz az egyenletünk: Z Z Z p∇u∇v + {z | Ω Ω a(u,v) fv quv = Ω } | {z } ϕ(v) Ebben az egyenletben azonban már csak u els® rend¶ parciális deriváltjai szerepelnek, tehát gyengíthetjük a megoldásra tett feltételünket. 14 http://www.doksihu Az egy dimenziós esethez hasonlóan deniáljuk a H 1 (Ω) Hilbert teret, azaz legyen H 1 (Ω) = {u ∈ L2 (Ω) : ∂1 u ∈ L2 (Ω), . , ∂n u ∈ L2 (Ω)}, ahol a deriválást ismét disztribúciós értelemben értjük ([3]) A skalárszorzat hu, viH 1 (Ω) = hu, viL2 (Ω) +h∇u, ∇viL2 (Ω) és az indukált norma kuk2H 1 (Ω) = kuk2L2 (Ω) + k∇uk2L2 (Ω) . Továbbá az egy dimenziós esethez hasonlóan legyen H01 (Ω) = {u : u ∈ H 1 (Ω), u|∂Ω = 0}, ahol u|∂Ω az úgynevezett nyom operátor ([3]). 3.11 Állítás H 1 (Ω) elemei nem
feltétlenül folytonosak 3.11 Deníció Egy u ∈ H01 (Ω) függvényt a (31) peremértékfeladat gyenge megoldásának nevezzük, ha Z Z Z p∇u∇v + Ω quv = Ω fv Ω teljesül minden v ∈ H01 (Ω) esetén. A továbbiakban szükségünk lesz az egy dimenzióban már használt Sztyeklov egyenl®tlenség magasabb dimenziós változatára. 6. Tétel (Magasabb dimenziós Sztyeklov egyenl®tlenség) Létezik c > 0, hogy minden u ∈ H01 (Ω) függvényre Z Z 2 |∇u|2 . |u| ≤ c Ω Ω Ezek után szeretnénk kimondani az egy dimenziós esethez hasonlóan egy egzisztencia illetve unicitás tételt. A 2 fejezetben tárgyaltak miatt elég belátnunk a(·, ·) elliptikusságát illetve folytonosságát Megjegyezzük, hogy ϕ(·) változatlanul folytonos lineáris funkcionál hiszen csak az integrációs tartomány változott meg. a(·, ·) elliptikussága: Z Z 2 |a(u, u)| = 2 (p|∇u| + q|u| ) ≥ Ω Ω p0 |∇u|2 = p0 |u|2H 1 (Ω) 0 a(·, ·)
folytonossága: Z Z |a(u, v)| = Z (p∇u∇v + quv) ≤ c Ω (∇u∇v + uv) ≤ c Ω sZ sZ |∇u|2 |∇v|2 = c|u|H01 (Ω) |v|H01 (Ω) ≤c Ω Ω 15 2∇u∇v Ω http://www.doksihu 3.11 Következmény A (31) peremértékfeladatnak mindig létezik egyértelm¶ megoldása. Ezek után diszkretizáljuk a problémát a 2. fejezetben elmondottak alapján, azaz ismét meg kell adnunk egy Vh véges dimenziós teret és ebben kell az uh közelít® megoldást keresnünk. 3.11 Megvalósítás Egy dimenzióban az I = (0, 1) invervallum egy felosztását tekintettük és ezen deniáltunk bázisfüggvényeket. Több dimenzióban is hasonlóan szeretnénk eljárni, ehhez azonban célszer¶ konkrétan megválasztani a dimenziót és az Ω tartományt. Legyen n = 2 és Ω = (0, 1) × (0, 1) vagyis az egységnégyzet. Vegyük a tartomány egy felosztását, mégpedig az alábbi rácsot (3.1 ábra): 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 0.2 0.4 0.6 0.8 1 3.1 ábra Ω
rácsozása Meg kell adnunk a Vh teret illetve egy bázisát. Az el®z® fejezet mintájára most gúlákat kell deniálnunk a bels® pontokra (hiszen a perem pontokon a megoldás értéke nulla). A feladat összetettebb az egy dimenziós esetnél, ezért célszer¶ megszámozni külön a pontokat és külön az Ωi háromszögeket Legyen i a tartomány egy bels® pontja. Az i pontot 6 darab háromszög veszi körül, tehát az i-hez tartozó Φi 16 http://www.doksihu bázisfüggvény 6 darab síkból fog állni. A síkokat úgy kell meghatároznunk, hogy az i pontban 1 legyen az értékük, a maradék 6 pontban viszont nulla.(32 ábra) 1 0.8 0.6 0.4 0.2 0 0 2 1.5 0.5 1 1 0.5 1.5 2 0 3.2 ábra Φi bázisfüggvény Meg kell határoznunk a Φi bázisfüggvényt alkotó síkok egyenletét. Ezt lokálisan fogjuk végezni, azaz csak az i pont körüli 6 darab háromszöget tekintjük és ezekre adjuk meg a megfelel® sík egyenleteket. Érdemes deniálnunk egy referencia
háromszöget, Ω0 -t, melynek minden csúcsához hozzárendelünk egy függvényt, mégpedig úgy, hogy az adott pontban legyen 1, a maradék két pontban pedig nulla. A sík egyenleteket úgy kapjuk meg, hogy egy lineáris transzformációval az adott háromszögre transzformáljuk Ω0 -t és az i ponttal megegyez® csúcsához rendelt függvény lesz a keresett egyenlet. Ezzel a módszerrel tehát megkaptuk a Φi bázisfüggvényeket Mostantól legyen N a rácspontok száma és keressük a közelít® megoldást a következ® alakban: uh (x, y) = N X αi Φi (x, y) i=1 A feladatunk ismét az, hogy az αi együtthatókat meghatározzuk. PN P ϕ(Φj ) = a(uh , Φj ) = a( N i=1 αi a(Φi , Φj ), j = 1, . , N vagyis az i=1 αi Φi , Φj ) = αi együtthatókra egy Aα = F lineáris egyenletrendszer adódik, ahol Ai,j = a(Φi , Φj ) (i, j = 1, . , N ), α = (α1 , , αN )T , F = (ϕ(Φ1 ), , ϕ(ΦN ))T 17 http://www.doksihu 3.2 További peremfeltételek A homogén
Dirichlet esetben láttuk a módszert, a 3. fejezetben azonban az inhomogén Dirichlet és Neumann peremfeltételeket is fogjuk használni Ebben a részben ezen esetekre röviden bemutatjuk a módszert, els®sorban a különbségeket hangsúlyozva. 3.21 Inhomogén Dirichlet peremfeltétel Tekintsük az alábbi peremérték feladatot: −div (p∇u) + qu = f (3.2) u|∂Ω = c ahol c ismert függvény. Keressük a megoldást u = w + u0 alakban, ahol u0 |∂Ω = c és w|∂Ω = 0. Tegyük fel, hogy v ∈ C 1 (Ω)-beli függvény, melyre v|∂Ω = 0, szorozzuk meg v -vel a (3.2) egyenletet és integráljunk Ω-n A Green tétel felhasználásával adódik, hogy Z Z Z Z ∂u p∇u∇v + quv − p vdσ = fv Ω Ω ∂Ω ∂η Ω azaz Z Z Z ∂(w + u0 ) p∇(w + u0 )∇v + q(w + u0 )v − p vdσ = ∂η Ω Ω ∂Ω Z fv Ω felhasználva, hogy w az Ω határán elt¶nik adódik, hogy Z Z Z Z Z ³Z ´ ∂u0 p∇w∇v + qwv − p p∇u0 ∇v + qu0 v vdσ = fv − ∂η Ω
Ω ∂Ω Ω Ω Ω vagyis az eddigi jelölést felhasználva a(w, v) = ϕ(v) − a(u0 , v) ∀w ∈ H01 (Ω) g = ϕ(v) − a(u0 , v) továbbra is folytonos lineáris ahol u0 ismert függvény, és ϕ(v) funkcionál. 3.21 Állítás A megoldás független u0 választásától A megvalósítás során meg kell adnunk a peremen lév® pontokhoz tartozó bázisfüggvényeket, ez megtehet® a már ismert referencia háromszöges módszerrel. Szemléletesen ezek a bázisok a bels® pontok gúláinak a fele, illetve a negyede 18 http://www.doksihu 3.22 Neumann peremfeltétel Legyen most a vizsgált feladat a következ®: −div (p∇u) + qu = f p (3.3) ∂u |∂Ω = c ∂η ahol c ismert függvény, továbbá q > 0. Ismét szorozzuk meg a (3.3) egyenletet egy v ∈ H 1 (Ω) függvénnyel, melyre v|∂Ω = 0 és integráljunk az Ω tartományon. Felhasználva a Green tételt, azt kapjuk, hogy Z Z Z Z ∂u p∇u∇v + quv − p vdσ = fv Ω Ω ∂Ω ∂η Ω azaz Z Z
Z p∇u∇v + Ω Z quv − cvdσ = Ω ∂Ω fv Ω Legyen ψ(·) folytonos lineáris funkcionál, melyre Z ψ(c) = cvdσ ∂Ω Ekkor az egyenletünk a(u, v) = ϕ(v) + ψ(c) ∀v ∈ H01 (Ω) g = ϕ(v) + ψ(c) nyilván folytonos lineáris funkcionál. ahol c ismert és ϕ(v) 3.23 Harmadfajú peremfeltétel Most a következ® peremérték feladatot vizsgáljuk −div (p∇u) + qu = f p (3.4) ∂u |∂Ω + σu|∂Ω = c ∂η ahol σ ≥ σ0 > 0. Az egyenletet megint szorozzuk egy v ∈ H 1 (Ω) függvénnyel, majd integrálunk és a Green tételt alkalmazva adódik, hogy Z Z Z p∇u∇v + quv − Ω Ω Z ∂Ω 19 (c − σu) vdσ = | {z } p ∂u ∂η fv Ω http://www.doksihu azaz Z Z p∇u∇v + | Ω Z Z quv + {z Ω ∂Ω Z fv + cvdσ σuvdσ = Ω ∂Ω } | {z } a(u,v) ϕ(v) Igazolható, hogy a(·, ·) folytonos és elliptikus, valamint ϕ(·) folytonos lineáris funkcionál. 20 http://www.doksihu 4. fejezet Peremfeltételek
vizsgálata Az el®z® fejezetekben bemutattuk a végeselem módszert különböz® peremfeltételek esetén. Kérdés, hogy mit mondhatunk a pontosságról, ha ugyanazon egyenlethez más-más peremfeltételeket csatolunk és az esetleges hibákat hogyan lehet kiküszöbölni. Tekintsük az alábbi egyenletet −4u + u = f (4.1) Amely az el®z® fejezetben szerepl® −div (p∇u) + qu = f egyenlet speciális esete p ≡ q ≡ 1 esetén. Ezentúl mindig a fenti (41) egyenlet közelít® megoldását fogjuk vizsgálni egy speciális L-alakú Ω = (−1, 1)2 (−1, 0)2 tartományon. Megmutatjuk, hogy az origó körüli úgynevezett sarokpont környezetében a közelítés pontossága rosszabb mint a tartomány más részein. A számításokhoz el®ször megadjuk az u pontos megoldást, majd az egyenletbe behelyettesítve meghatározzuk az f függvényt, végül az így kapott egyenlet uh közelít® megoldását u-hoz viszonyítva következtetünk a hibára. Kérdés, hogy mi
legyen a mér®szám u és uh között, az alábbiakban bevezetünk két deníciót. 4.01 Deníció (Abszolút hiba) Legyen u az egyenlet pontos megoldása, uh pedig a közelít® megoldása. Ekkor a végeselem módszer abszolút hibája Eabs (Ω) = ku − uh kH 1 (Ω) Az abszolút hiba képlet az el®z® fejezetben tárgyaltak miatt természetesen adódik. A szakirodalomban azonban sokszor a relatív hibával számolnak és a vizsgálataink során mi is ezt fogjuk használni. 21 http://www.doksihu 4.02 Deníció (Relatív hiba) Legyen u az egyenlet pontos megoldása, uh pedig a közelít® megoldása. Ekkor a végeselem módszer relatív hibája Erel (Ω) = ku − uh kH 1 (Ω) kukH 1 (Ω) 4.1 Hibaszámítás A numerikus számításokhoz a COMSOL Multiphysics nev¶ szoftvert fogjuk használni, röviden ismertetjük ennek menetét. A program indításakor ki kell választanunk a megoldani kívánt egyenlet típusát, majd meg kell adnunk az együtthatókat és el kell
végeznünk perem beállításokat, végül a COMSOL el®állítja az uh közelít® megoldást. Esetünkben az alábbi Helmholtz egyenletet választjuk −∇(c∇u) + au = f az együtthatók legyenek a ≡ c ≡ 1 ekkor a vizsgálni kívánt (4.1) egyenletet kapjuk A COMSOL-ban lehet®ség van a közelít® megoldást tartalmazó kifejezés integrálásra, ezt fogjuk használni az abszolút illetve relatív hiba kiszámításához. A rácsépítéshez az alapértelmezett beállításokat használjuk, ez azt jelenti, hogy a tartományt hegyesszög¶ háromszögekkel rácsozzuk. Lehet®ségünk van a rácsozás nomítására ezt a COMSOL automatikusan elvégzi a kívánt területen Bázisunk az el®z® fejezetben tárgyalt lineáris bázis lesz. 4.2 L-alakú tartomány A vizsgálatok során mindvégig az Ω = (−1, 1)2 (−1, 0)2 tartományon fogjuk tekinteni az egyenletünket. Felvet®dik a kérdés, hogy miért választunk ilyen speciális tartományt a megszokott (0, 1)2
egységnégyzet helyett. Ennek az az oka, hogy az u megoldás sarokpont (vagyis az origó) körüli viselkedése (pl. szingularitása) nagyban befolyásolja a megoldás pontosságát Kérdés, hogy a peremfeltétel ismeretében hogyan válasszuk meg a rácsfelosztás nomságát. Legyen Sr := {(x, y) ∈ R2 : x2 + y 2 ≤ r} vagyis az origó közep¶ r sugarú körlap és Kr := Sr ∩ Ω. 22 http://www.doksihu A hibákat két tartományon fogjuk számítani az egész Ω tartományon, valamint a K 1 körlapon. 2 Arra számítunk, hogy ha az egész tartományon nagy a relatív hiba akkor a K 1 2 körben még nagyobb lesz, vagyis az origó körül "romlik el" a közelít® megoldás. 4.3 Vizsgálat Bontsuk fel az Ω tartomány határát diszjunkt módon hat darab részre az alábbi módon: Γ1 = {1}×[−1, 1), Γ2 = [0, 1)×{−1}, Γ3 = {0}×(−1, 0], Γ4 = [−1, 0)×{0}, Γ5 = {−1} × (0, 1], Γ6 = (−1, 1] × {1} ekkor nyilván ∂Ω = Γ1 ∪ Γ2 ∪ Γ3
∪ Γ4 ∪ Γ5 ∪ Γ6 . A vizsgálatok során csak az origóval érintkez® Γ3 és Γ4 szakaszokon fogjuk változtatni a peremfeltételt, a többi szakaszon végig Dirichlet peremfeltételt fogunk megadni. Mindig az alábbi peremérték feladatot fogjuk megoldani: −4u + u = f u|∂ΩΓ = g ∂u |Γ = f ∂η 23 (4.2) http://www.doksihu ahol az f , g , h függvényeket úgy állítjuk be, hogy az el®re adott u0 legyen a megoldás, azaz f = −4u0 + u0 , g = u0 |∂ΩΓ , f = ∂u0 | . ∂η Γ A Γ halmaz tartalmazza a perem azon részeit, ahol Neumann peremfeltételt adunk meg. 4.31 Végtelen szingularitású függvény vizsgálata Legyen a pontos megoldásunk u0 (x, y) = ln((x − x0 )2 + (y − y0 )2 ). Ennek a függvénynek az (x0 , y0 ) pontban szingularitása van, pontosabban a határértéke −∞. Válasszuk az (x0 , y0 ) pontot úgy, hogy az L-alakú tartomány sarkán kívül legyen, 1 1 , − 1000 ). de közel az origóhoz, például legyen (x0 , y0 ) =
(− 1000 Tehát a pontos megoldás u0 (x, y) = ln((x + 1 )2 1000 + (y + 1 )2 ). 1000 2 0 −2 −4 −6 −8 −10 −12 −14 1 0.5 0 −0.5 −1 −1 −0.5 0 0.5 1 4.1 ábra Végtelen szingularitású függvény El®ször legyen Γ az üres halmaz, vagyis ekkor az egész peremen Dirichlet peremfeltételt adunk meg. Az alábbi táblázat tartalmazza a relatív illetve abszolút hibákat Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 0.71 8.00 0.75 7.94 2952 0.58 6.56 0.61 6.55 47232 0.38 4.32 0.41 4.31 24 http://www.doksihu Azt látjuk, hogy a relatív hiba nagyobb a K 1 körben mint az egész Ω tar2 tományon, ez pont azt jeletni, hogy az szingularitás körül rosszabbul tudjuk közelíteni a megoldást. Nézzük mi történik ha Γ = Γ3 azaz az egyik origóval érintkez® szakaszon Neumann, a másikon pedig Dirichlet peremfeltételt adunk meg Megjegyezzük, hogy az u pontos megoldásunk forgás szimmetrikus,
ezért a Γ = Γ4 választással ugyanezeket az adatokat kapnánk. Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 0.68 7.67 0.71 7.52 2952 0.55 6.27 0.59 6.21 47232 0.37 4.23 0.40 4.22 A hibák ugyanazt mutatják mint az el®z® esetben, vagyis az origó körül még mindig nagyobb az eltérés, bár valamelyest kisebb hibákat kaptunk. Nézzük mi történik ha Γ = Γ3 ∪ Γ4 tehát amikor az origóval érintkez® szakaszokon Neumann peremfeltételt adunk meg a ∂Ω többi részén pedig Dirichlet peremfeltételt. Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 0.84 9.51 0.85 9.02 2952 0.85 9.57 0.86 9.13 47232 0.53 6.02 0.56 5.88 Látható, hogy a hibák jelent®sen megn®ttek és K 1 -ben még mindig nagyobb a 2 relatív hiba mint az egész tartományon. Milyen következtetést vonhatunk le ebb®l? Összehasonlítva a három táblázat adatait, megállapíthatjuk, hogy Neumann
peremfeltétel esetén érdemes mindig nomabb felosztást alkalmazni mint, amit Dirichlet feltétel esetén alkalmaznánk. 4.32 Oszcilláló függvény vizsgálata Legyen most a pontos megoldásunk u0 (x, y) = sin(nπx)sin(mπy) ahol n és m választásával az oszcilláció nagyságát tudjuk változtatni. Például n = 10 és m = 20, ekkor u0 (x, y) = sin(10πx)sin(20πy) adódik. 25 http://www.doksihu 1 0.5 0 −0.5 −1 1 0.5 1 0.5 0 0 −0.5 −0.5 −1 −1 4.2 ábra Oszcilláló függvény Legyen Γ most is el®ször az üres halmaz, a hibákat az alábbi táblázat mutatja Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 1.23 75.29 1.33 35.74 2952 0.79 47.82 0.82 22.14 47232 0.23 13.92 0.24 6.45 Azt kapjuk, hogy az oszcillációt csak nagy felosztásra sikerül közelíteni, viszont a relatív hiba jelent®sen csökken a háromszögek számának növelésével. A tartományokat illet®en továbbra is azt látjuk,
hogy az origó körül nagyobb a hiba. Tekintsük a Γ = Γ3 esetet Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 1.36 83.00 1.50 40.50 2952 0.79 47.82 0.82 22.15 47232 0.23 13.92 0.24 6.45 Látszik, hogy kevés háromszögre nagyobb hibák adódnak, de a nomítás után ugyanazokat a hibákat kapjuk. Hasonlóan Γ = Γ4 esetén is ezt látjuk az alábbi táblázat adataiból: 26 http://www.doksihu Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 1.26 76.80 1.38 37.15 2952 0.79 47.79 0.82 22.11 47232 0.23 13.92 0.24 6.45 Legyen Γ = Γ3 ∪ Γ4 , ekkor a hibák Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 1.40 85.25 1.49 40.23 2952 0.79 47.79 0.82 22.11 47232 0.23 13.92 0.24 6.45 Azt kaptuk, hogy ha mindkét szakaszon Neumann peremfeltételt adunk meg, akkor a hibák nem változnak jelent®sen, vagyis u0 oszcillációjának közelítése
szempontjából csak a felosztás nomsága a lényeges és a peremfeltétel nem befolyásolja a pontosságot. Ennek az a magyarázata, hogy jelen esetben a Neumann peremfeltételben konstans függvények szerepelnek vagyis nem jelent nehézséget a közelítésben 4.33 Hirtelen ugrású függvény vizsgálata Legyen u0 (x, y) = arctan(m( p (x − x0 )2 + (y − y0 )2 − r), ahol (x0 , y0 ) pont je- lenti, a falat alkotó henger alsókörének középpontját, r pedig a sugarát. Válasszuk ezt p 7 , m = 60, ekkor u0 (x, y) = arctan(60( x2 + y 2 − a pontot az origónak és legyen r = 10 7 )). 10 Legyen Γ = ∅ ekkor az alábbi adatokat kapjuk Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 1.09 19.33 4.43 5.20 2952 0.46 8.14 0.84 0.99 47232 0.13 2.37 ≈0 ≈0 Jelent®s különbség van a K 1 tartomány és az egész Ω tartományon kapott 2 hibák között, azaz a sarokpont körül lényegesen rosszabb a közelítés, még 2952 27
http://www.doksihu 2 1 0 −1 −2 1 1 0.5 0.5 0 0 −0.5 −0.5 −1 −1 4.3 ábra Hirtelen ugrású függvény 1 háromszögre is. Azonban ha jelent®sen megnöveljük a nomítást, akkor éppen a sarokpont körül lesz jobb a közelítés, igaz az egész tartományon is nagy javulást tapasztalunk. Ugyanezt tapasztaljuk ha Γ = Γ3 , ami megegyezik a Γ = Γ4 esettel. Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 1.11 19.68 4.56 5.36 2952 0.46 8.18 0.84 0.98 47232 0.13 2.37 ≈0 ≈0 Az oszcilláló függvényhez hasonlóan ha Γ = Γ3 ∪Γ4 , akkor is lényegében a peremfeltételt®l függetlenül ugyanezeket az adatokat kapjuk. Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 1.11 19.78 4.55 5.34 2952 0.46 8.22 0.79 0.93 47232 0.13 2.37 ≈0 ≈0 Jól látszik, hogy a hibák a felosztás növelésével gyorsabban csökkenek mint az el®z® két függvény esetében, vagyis a
hirtelen ugrás jól kezelhet®. Azonban meglep®, 28 http://www.doksihu hogy a csökkenés sebessége az origó körül sokkal gyorsabb mint az egész tartományon. Felvet®dik a kérdés, hogy ha a falat alkotó henger alapköre a K 1 tartományba 2 esik, akkor is ilyen jó közelítést kapunk-e vagy éppen fordítva, azért javult ilyen gyorsan a közelítés mert a fal a K 1 tartományon kívül esett. A kérdés megválas2 zolásához újabb vizsgálatokat kell végezünk ezúttal r < 1 2 sugárral. 2 1 0 −1 −2 1 1 0.5 0.5 0 0 −0.5 −0.5 −1 −1 4.4 ábra Hirtelen ugrású függvény 2 Legyen u0 (x, y) = arctan(60( Háromszögek száma p x2 + y 2 − 41 )) és tekintsük el®ször a Γ = ∅ esetet: Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 0.95 10.32 0.95 10.09 2952 0.46 5.03 0.47 5.01 47232 0.13 1.42 0.13 1.42 Látható, hogy megsz¶nt az el®z®ekben tapasztalt tendencia, miszerint a bels® körön kapott hibák
gyorsabban csökkennek az egész tartományon számolt hibákhoz képest, s®t csak nagyon kis eltérést tapasztalunk. Ennek oka abban rejlik, hogy az u pontos megoldás az ugrás után közelít®leg konstansként viselkedik, amit alacsony felosztás mellett is jól tudunk közelíteni. Tehát az Ω K 1 tartományon már nagyon 2 jó a közelítés, ezért a tényleges hibát csak a bels® körben fellép® rossz közelítés adja. 29 http://www.doksihu Ugyanezt tapasztaljuk ha Γ = Γ3 , ami megegyezik a Γ = Γ4 esettel. Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 0.95 10.34 0.95 10.10 2952 0.47 5.10 0.48 5.06 47232 0.13 1.42 0.13 1.41 Végül a Γ = Γ3 ∪ Γ4 eset. Háromszögek száma Erel (Ω) Eabs (Ω) Erel (K 1 ) Eabs (K 1 ) 2 2 738 1.20 13.02 1.11 11.74 2952 0.48 5.16 0.48 5.11 47232 0.13 1.42 0.13 1.42 Azt kapjuk, hogy a hibák valamelyest n®ttek de fennáll a tendencia, miszerint a két vizsgált
tartományon azonos hibák adódnak. 4.4 Konklúzió Láttuk, hogy mindhárom függvénytípus esetén a sarokpont körül rosszabb a közelítés mint a tartomány más részein, ugyanakkor ez a felosztás nomításával orvosolható. A különböz® peremfeltételekre azonban csak a végtelen szingularitású függvény esetében kaptunk lényegesen nagyobb hibákat. Tehát ha a tartományunk a vizsgált L-alakhoz hasonló, mindenképpen érdemes a sarokpont környezetében nomabb felosztás választanunk. 30 http://www.doksihu Köszönetnyilvánítás Szeretném megköszönni témavezet®imnek, Simon Péternek és Horváth Tamásnak, hogy felhívták gyelmemet erre az érdekes témára. Valamint, hogy mindig szakítottak rám id®t és rengeteg hasznos tanáccsal láttak el. 31 http://www.doksihu Irodalomjegyzék [1] Gergó Lajos: Numerikus módszerek, egyetemi jegyzet, ELTE, Eötvös Kiadó, (2000) [2] Horváth Tamás: Végeselem módszerek alapjai, Budapest (2009)
(elektronikus jegyzet) [3] Simon László, E. A Baderko: Másodrend¶ lineáris parciális dierenciálegyen- letek [4] Stoyan Gisbert, Takó Galina: Numerikus módszerek 3., ELTE-Typotex, Budapest (1997) [5] P. Solin: Partial Dierential Equations and the Finite Element Method, John Wiley and Sons, (2005) [6] P. Solin, J Cerveny, I Dolezel: Arbitrary-level hanging nodes and automatic adaptivity in the hp-FEM, Math. Comput Simul 77 117 - 132, (2008) 32