Softwarová pomůcka pro 2D generaci sítě konečných prvků Ing. Filip Hejnic*, doc.Ing. Petr Štemberk, Ph.D.** *České Vysoké Učení Technické v Praze, Thákurova 7, 166 29 Prague 6, E-mail:
[email protected] ** České Vysoké Učení Technické v Praze, Thákurova 7, 166 29 Prague 6
Abstrakt V inženýrské praxi je v dnešní době nejvíce používanou metodou pro výpočet vnitřních sil Metoda konečných prvků (MKP). Pro využití této metody je nezbytné rozdělení potřebné oblasti na konečné prvky. Všechny komerční programy si chrání svoje vlastnictví autorskými právy a je těžké využít jejich kód pro generování sítě zdarma. Na internetu lze najít různé knihovny napsané v programovacích jazycích, ale pro jejich použití je nutné umět je použít. Tento článek se snaží vyhnout těmto úskalím a poskytnout nástroj pro generování sítě zdarma. Uložení souřadnic sítě do textového souboru je samozřejmostí. V programu je zahrnuto zjemňování a vyhlazování sítě Klíčová slova: generace sítě, metoda konečných prvků, trojúhelníky, síť 1. Úvod Jak již bylo řečeno v abstraktu článku tak dnes ve všech technických odvětví při numerickém modelování je nejvíce rozšířenou metodou metoda konečných prvků. Obzvláště ve stavebnictví je tato metoda široce využívána, jak v praxi, tak i při vědeckých výzkumech. Základem této metody je rozdělení zkoumané oblast na konečný počet částí tak, aby „dobře“ pokrývaly uvažovanou oblast. Zkoumaný prvek může být například táhlo (1D prvek), železobetonová deska (2D prvek) nebo ostatní prvky zkoumané ve třech rozměrech. Článek se zabývá druhou oblastí těles 2D prvky. Do této skupiny se zahrnují většinou prvky, které lze výpočetně zjednodušit z trojrozměrného prostoru, jako jsou desky, skořepiny a další podobné typy. Rozdělení dvojrozměrné oblasti na konečné prvky je závislé na druhu prvků, které při výpočtu budeme chtít využít. Mezi dnes používané 2D prvky patří trojúhelníky, čtyřúhelníky (Obr. 1) a jejich kombinace.
Obrázek 1. Typy používaných 2D prvků [3]
Jednotlivé druhy prvků jsou spojeny s typem metody, kterou lze použít pro automatické vytvoření celé sítě konečných prvků. Pokud chceme rozdělit danou oblast na prvky, je snadnější využít trojúhelníku, jakožto základního prvku. Trojúhelníkové sítě se snáze vytvářejí z důvodu jednoduššího matematického popisu trojúhelníku. Při použití trojúhelníkových prvků je nutné k pokrytí téže oblasti větší množství prvků než při použití čtyřúhelníkových, jak je vidět z obrázku 1. Při použití čtyřúhelníkových prvků dostáváme lepší řešení než u trojúhelníkových z důvodu většího počtu vrcholů. Jelikož při potřebě vygenerování sítě, pro universitní projekt, nebyl k dispozici žádný nástroj, který by bylo možno použít zdarma, vznikl nápad naprogramovat aplikaci, která umožní vygenerování sítě a bude dostupná na internetu. Hlavním cílem bylo vytvořit jednoduchý program, ve které bude možno zadat polygon a vygenerovat síť s možností uložení vrcholů sítě do souboru.
2. Základní metody generování sítí Základním dělením sítí konečných prvků jsou sítě strukturované a nestrukturované. Strukturované jsou takové, jejichž prvky jsou pravidelně uspořádány (pravidelná obdélníková síť). Nestrukturované sítě jsou složeny převážně z trojúhelníkových, ve vyšších dimenzích pak ze čtyřstěnných prvků, nepravidelného vzoru. Tento článek se dále zabývá vytvářením druhého typu sítí. Pro vytváření nestrukturovaných sítí existují tři základní metody. Metoda Quad/Octree je založena na postupném dělení dané oblasti a uchovávání informací o těchto částech ve struktuře stromu (Obr. 2a). Každá část, která je rozdělena na menší části z důvodu zjemnění sítě je uložena, jako podoblast své nadřazené oblasti. Oblast je dělena pravidelně na obdélníkové nebo, v případě prostorového útvaru, na hranolové části. Zjemňování probíhá, dokud není dosaženo požadované jemnosti v místech, které jsou pro to určené. Po takto připravené síti je možné do každého pravidelného obrazce vložit šablonu trojúhelníků, popř. čtyřstěnů, které nám danou oblast vyplní sítí. Při 2D rozměru je počet šablon, které je možné použít 16 (Obr. 2b). Na obrázku jsou zobrazeny pouze šablony, které jsou jiného typu. Pro obecný 3D prostor narůstá počet šablon na počet 4096. Výhodou této metody je rychlost, lze ji použít pro jakoukoliv oblast a vytváří dobrou kvalitu sítě. Nevýhodou je převážně opakující se velikost prvků. Sít lze rozdělit pouze na násobky jednotlivých prvků.
Obrázek 2.(a) Struktura dat metody Quad/Octree; (b) Šablony pro 2D síť
Další metodou pro tvorbu sítí je metoda postupující fronty. Metoda je založena, jak sám název napovídá, na postupném doplňování trojúhelníkových obrazců od zvoleného počátečního bodu. Základem je rovnoměrné rozdělení hranice oblasti na úseky. Tyto části poté slouží jako strany rovnostranných trojúhelníků, které jsou v první fázi generovány. Problém vzniká převážně v rozích, kde se potkávají dva obrazce. Toto lze řešit tak, že pokud v rohu již jeden obrazec je, tak zde nevznikne nový, ale využije se stávajícího již vloženého vrcholu. Zbývající část, která ještě není pokryta elementy, lze doplnit stejným způsobem pomocí vytváření rovnostranných obrazců. Velkou výhodou metody postupující fronty je velká kvalita a spolehlivost 2D sítě (nejkvalitnější prvky kolem hranice). K nevýhodám patří zejména rychlost algoritmu a možnost tvorby prvků malých objemů.
Obrázek 3. Postup tvorby sítě metodou postupující fronty
Poslední ze základních typů metod je Delaunayho triangulace. Metoda byla vynalezena ruským matematikem a fyzikem B. N. Delaunay. Principem postupu je vytvořit síť (Obr. 4a), která bude složena pouze z takových trojúhelníku, jejichž kružnice opsané nebudou obsahovat žádný jiný bod sítě (pouze vrcholy daného prvku). Pokud toto pravidlo nebude splněno, nejedná se o Delaunayho triangulaci (Obr. 4b). Prvním krokem tvorby sítě je vytvoření ohraničující obálky (tzv. bounding box) kolem dané oblasti a rovnoměrné rozdělení hranice. Následuje vytvoření počáteční sítě, která splňuje víše popsané podmínky. Dalším krokem je vkládání vnitřních bodů. Pro tento postup existuje několik
možností jak tento problém vyřešit. Zaprvé je možné vkládat uzly do těžišť již vzniklých trojúhelníků. Tento způsob není moc vhodný, protože vytváří obrazce s malým vnitřním úhlem, ale je snadno použitelný. Za druhé je možné strany obrazců rozdělit na poloviny a získat tak čtyři trojúhelníky v jednom obrazci. Nejrozšířenějším postupem je použití iterativního Bowyer-Watsonova algoritmu. Základním principem je vložení bodu do sítě (např. do těžiště největšího obrazce) a zjištění, jestli tento bod neporušuje zásadu triangulace. Pokud ano, všechny obrazce, jejichž opsaná kružnice obsahuje daný bod, jsou smazány a vytvoří se nové, které začlení tento bod jako vrchol sítě. Proces se opakuje, dokud velikost největšího útvaru není přijatelná (např. plocha všech trojúhelníku je v rozmezí zadaných hodnot). Hlavními výhodami této triangulace je malá výpočetní náročnost algoritmu, dobrá kvalita sítě (dle vkládání bodů) a velmi silný matematický popis, který dovoluje metodu rozšířit jak ve 2D, tak i do 3D. Nevýhodou je především obtížné vyjmutí hranice, z ohraničující obálky, pokud je složitého tvaru. Delaunayho triangulace byla zvolena jako vhodná metoda pro implementaci do vyvíjeného programu.
Obrázek 4. (a) Delaunayho triangulace, (b) Rozdíl mezi triangulacemi
3. Postup vytvoření sítě zadaného obrazce Základním vstupním požadavkem před samotným procesem je nutné určit polygon nebo čtvercovou oblast, kterou budeme chtít pokrýt sítí. Dochází pouze k určení x-ových a y-ových souřadných bodů zadávané oblasti. Tyto souřadnice jsou uloženy a dále vstupují do výpočtu. V prvním kroku výpočtu je určena obálka, která obklopuje celý obrazec. Obálka je tvaru obdélníka a jeho vrcholy jsou vypočteny z již dostupných souřadnic a přidány do seznamu bodů. S tímto vstupním seznamem začíná první triangulace. Obálka je následně rozdělena na dva stejné trojúhelníky, tím začíná proces generování sítě. Algoritmus procesu triangulace je uveden v (1.1) Výstupními hodnotami je seznam, kde každý řádek odpovídá vytvořenému trojúhelníku. Každý řádek obsahuje tři čísla, která označují pozici souřadnice v seznamu souřadnic. Vytvoření prázdného seznamu trojúhelníků a zařazení dvou počátečních trojúhelníků Pro každý bod v seznamu Vytvoření vyrovnávacího seznamu Pro každý trojúhelník v seznamu Spočti střed a poloměr kružnice opsané Když bod leží uvnitř kružnice Přidej hrany trojúhelníku do vyrovnávacího seznamu Odstraň trojúhelník ze seznamu trojúhelníků End End Odstraň všechny zdvojené okraje z vyrovnávacího seznamu Přidej do seznamu všechny trojúhelníky vytvořené mezi bodem a okrajem polygonu End Odstraň trojúhelníky a vrcholy, které jsou mimo obrazec End
(1.1)
Po vytvoření základní sítě dochází k postupnému zjemňování. V závislosti na velikosti obrazce dojde k rovnoměrnému rozdělení hranice a tím přidání bodů. Následným přidáváním bodů do těžišť trojúhelníků dojde ke zjemnění. Po každém novém vložení bodů dojde k novému přetransformování seznamu trojúhelníku podle víše uvedeného postupu (1.1). 3.1. Vyhlazování sítě Pro lepší kvalitu sítě byl použit algoritmus změny tvaru jednotlivých trojúhelníku neboli „vyhlazovací“ algoritmus známí jako Lapalceovo vyhlazování. Před uvedením postupu vyhlazovacího algoritmu je vhodné se zmínit, co znamená kvalita trojúhelníku. V tomto kontextu je kvalita považována jako změna polohy vrcholů trojúhelníka, tak aby došlo ke změně vnitřních úhlů, pro zlepšení poměru mezi těmito úhly. Pro hodnocení kvality byla použita rovnice (1.2). 4A√3 (1.2) q= 2 ℎ1 + ℎ22 + ℎ32 kde A je plocha trojúhelníka, h 1 , h 2 , h 3 jsou délky stran trojúhelníka a q je výsledná hodnota kvality. Pokud q > 0,6 má daný trojúhelník vyhovující tvar. Nejlepší tvar je pokud q = 1, tj. budou všechny stany stejně dlouhé, jedná se tedy o rovnostranný trojúhelník, a tento tvar bude považován za optimální. Základním principem vyhlazování je vzít uzel uvnitř sítě. Vybrat všechny trojúhelníky, které mají tento uzel jako svůj vrchol a zhotovit polygon, který ohraničuje původní vnitřní uzel. Nová souřadnice vybraného uzlu sítě se poté spočítá jako souřadnice těžiště polygonu vytvořeného sousedními body (1.2; 1.3;1.4). 𝑛−1
1 𝐶𝑥 = �(𝑥𝑖 + 𝑥𝑖+1 )(𝑥𝑖 𝑦𝑖+1 − 𝑥𝑖+1 𝑦𝑖 ) 6𝐴 𝐶𝑦 = 𝐴=
(1.3)
𝑖=0
𝑛−1
1 �(𝑦𝑖 + 𝑦𝑖+1 )(𝑥𝑖 𝑦𝑖+1 − 𝑥𝑖+1 𝑦𝑖 ) 6𝐴
(1.4)
𝑖=0
𝑛−1
1 �(𝑥𝑖 𝑦𝑖+1 − 𝑥𝑖+1 𝑦𝑖 ) 2
(1.5)
𝑖=0
Kde C x , C y jsou nové souřadnice uzlu; A je plocha polygonu; n je počet uzlů tvořících polygon; x,y jsou jednotlivé souřadnice vymezující polygon. Souřadnice polygonu musejí být seřazeny podle postupu vzniku polygonu. Nesmějí být uloženy v přeházeném pořadí. Pro zajištění této skutečnosti byl implementován algoritmus Konvexní obálky. Konvexní obálka je jediný konvexní polygon, jehož vrcholy obsahují všechny body, které jsou v dané množině zadány. Tento problém byl nejvíce studován v počítačové geometrii. Nejvíce se používá v různých aplikacích, jako je detekce kolizí ve hrách, geografických informačních systémech a v robotice. Algoritmus výpočtu je ukázán v rovnici (1.5) Nahrání všech čar mezi danými body for i= 0 to i < celkový počet čar for j=0 to j< celkový počet bodů Získá body, které nejsou v aktuální čáře end if Zkontroluje, zda je čára na okraji konvexní obálky tak, že některý z bodů se nachází na jeho levé straně. Pokud ne odebere čáru ze seznamu, nastaví i a pokračuje else Zvýší i a pokračuje end end
(1.6)
Pokud se výše popsaný proces provede pro všechny uzly konečně prvkové sítě, tak dojde k „vyhlazení“. Rozdíl mezi vyhlazenou a nevyhlazenou sítí je patrný z obrázku 5. Všechny trojúhelníky se posunou blíže k optimálnímu výsledku. Tento proces lze opakovat, dokud není dosaženo
přijatelného výsledku. Bylo vypozorováno, že po 15 iteračních cyklech je posun vrcholů již zanedbatelně malý a další vyhlazování již není efektivní.
Obrázek 5. Rozdíl mezi vyhlazenou (tmavá) a nevyhlazenou (světlá) sítí
4. Použité prostředí Pro provedení výše zmíněných myšlenek byla zvolena platforma .NET Framework (Obr. 6) od společnosti Microsoft. Pod platformou .NET se skrývá mnoho možností jaký programovací jazyk využít. Pro tento projekt byl zvolen programovací jazyk C#, který se objevil počátkem roku 2006. Jazyk C# byl vytvořen zcela nově jako moderní, jednoduchý, objektově orientovaný a typově bezpečný programovací jazyk odvozen od C a C++.Relativní jednoduchost vytvoření uživatelského prostředí dělá z tohoto prostředí velmi silný nástroj. Platforma .NET je bohužel zatím spustitelná pouze pod operačními systémy Windows. Jelikož tento systém je využíván 95% počítačů neměl by být problém tento projekt spustit na jakékoliv dnešním počítači či notebooku s operačním systémem Windows XP či novějším, kde .NET Framework je nainstalován.
Obrázek 6. Schéma platformy .NET
4.1. Modul CLR Jádrem platformy je její běhové prostředí, které se označuje jako modul CLR (Common Language Runtime). Kód spuštěný pod kontrolou modulu CLR se označuje jako řízený kód (Manage code). Než tento modul může spustit námi napsaný program je nutné zdrojový kód přeložit do jazyka IL (Intermediate Language) a poté z jazyka IL do kódu pro cílovou platformu. Jazyk IL se překládá metodou Just-In-Time (JIT). Tento postup poskytuje dvě základní výhody. Zvýšení výkonu a budoucí přenositelnost aplikace na různé platformy. 4.2. Jazyk C# Pro vytvoření uživatelského prostředí je možné použít dvě možnosti. Starší, ale o to rozšířenější grafické rozhraní Windows Forms (Application Programming Interface) nebo novější rozhraní Windows Presentation Foundation (Rich User Interface) založeného na XAML (Extensible Application Markup Language). Pro tento projekt byla zvolena druhá z možností. Předpokládá se, že tento jazyk se bude velice rychle rozvíjet a vzniknout další obdobné projekty, které by bylo možné jednoduše spojit. Další důvodem je relativní jednoduchost vytvoření uživatelsky přátelského prostředí.
5. Program Uživatelské prostředí programu se skládá z nabídky, kde lze najít základní příkazy pro vytvoření sítě, jimiž jsou možnosti vytvořit, vyhladit a smazat síť. Dále je zde možnost zavřít program a uložit souřadnice jednotlivých obrazců sítě do textového souboru. Pod kontextovou nabídkou je umístěn panel, na kterém jsou umístěny příkazy pro vkládání 2D geometrie (bod, vícenásobná čára, obdélník). Geometrii lze zadávat pomocí myši, ale i pomocí příkazového řádku, kde lze zadat potřebnou hranici v souřadném systému. Jako počátek souřadnicového systému byl zvolen pravý horní roh kreslící plochy. Největšími přednostmi zmíněného programu je jednoduchost a dostupnost zdarma bez jakýchkoliv omezení. Na obrázku 6 lze vidět základní rozhraní programu a názorný příklad při postupu generování sítě. V pravém horním rohu je zobrazen základní obrazec pro vytvoření sítě. V levé dolní části je zobrazena počáteční síť, která je vytvořena prvotně a poté je z ní vytvořena finální síť (vpravo dole).
Obrázek 7. Postup generování sítě v prostředí FEMSoftware
Program je možné zdarma stáhnout z osobních stránek autora na Katedře betonových a zděných konstrukcí Fakulty stavební ČVUT v Praze (http://concrete.fsv.cvut.cz/~hejnic). Program není vázán žádnou licencí pro použití. 6. Závěr Článek ukazuje na možnost využití, dnes zdarma dostupných prostředků, pro vytvoření užitečných programů pro zjednodušení řešení inženýrských problémů. Byly zde popsány základní metody pro generování sítí konečných prvků a jejich základní vlastnosti. Je zde ukázáno zpracované řešení daného problému s využitím moderního programovacího jazyka C#, které je dostupné široké veřejnosti. Podrobněji jsou rozebrány detaily realizace vytvoření sítě konečných prvků při zadání libovolného polygonu. Dalším krokem, pro zlepšení funkčnosti programu, bude vylepšená implementace pro vkládání vnitřních bodů do sítě. Poděkování Článek vznikl za finanční podpory grantu Studentské grantové soutěže číslo SGS11/107/OHK1/2T/11 na fakultě stavební ČVUT v Praze. Literatura [1] P.J. Frey, P-L. George, Mesh Generation, HERMES Science Europe Ltd, 2000 [2] P-L. George, H. Borouchaku, Delaunay Triangulation and Meshing, Editions HERMES, Paris, 1998 [3] GAMBIT 2.3 Documentation,{http://my.fit.edu/itresources/manuals/gambit2.3/help/index.htm } [4] MATLAB Online Help Library, {http://mathworks.com } [5] Microsoft Online Help Library, {http://msdn.microsoft.com }
Seznam obrázků: Obrázek 8. Typy používaných 2D prvků [3] Obrázek 9.(a) Struktura dat metody Quad/Octree; (b) Šablony pro 2D síť Obrázek 10. Postup tvorby sítě metodou postupující fronty Obrázek 11. (a) Delaunayho triangulace, (b) Rozdíl mezi triangulacemi Obrázek 12. Rozdíl mezi vyhlazenou (tmavá) a nevyhlazenou (světlá) sítí Obrázek 13. Schéma platformy .NET Obrázek 14. Postup generování sítě v prostředí FEMSoftware