Femlab – nový nástroj pro simulaci

Femlab – nový nástroj pro simulaci Úvod

Výpočetní metody umožňují řešit množství náročných a složitých úloh v různých oborech. V současné době však nestačí „jen“ úlohy řešit, ale je třeba je také modelovat, simulovat a měnit různé parametry, abychom tak získali představu o tom, jak by se při jejich úpravě změnilo výsledné řešení. Často hovoříme o přístupu what if.

Femlab všechny tyto možnosti nabízí. V tisku již bylo o tomto programu uveřejněno několik informací, my se však v následujícím článku budeme Femlabem zabývat podrobněji.

Podstatou Femlabu je řešení parciálních diferenciálních rovnic (PDR), které popisují mnoho úloh ve fyzice a v různých vědních i technických oborech. Femlab tyto úlohy diskretizuje a k řešení používá metodu konečných prvků. PDR popisují úlohy v oblasti obecné fyziky, chemie, pružnosti a pevnosti, akustiky, elektromagnetismu, prostupu tepla, proudění tekutin, optoelektroniky a dalších. Při řešení parciálních diferenciálních rovnic je nezbytná znalost okrajových podmínek a vlastností prostředí uvnitř takto ohraničených oblastí. Z pohledu praxe tak modelujeme vnější vlivy, které působí na vytvořený model systému, a zároveň určujeme vlastnosti prostředí uvnitř modelu. Tímto postupem můžeme modelovat vnější síly, které působí na určitou součást nebo konstrukci, teplotu okolí, která ovlivňuje dané prostředí, nebo magnetické vlnění působící v definované oblasti.

Výhodou Femlabu je možnost zahrnout do jednoho řešení několik fyzikálních vlivů současně, např. zatěžování silou při současném působení okolní teploty nebo proudění s prostupem tepla.

Při práci s Femlabem máme k dispozici srozumitelné a intuitivní grafické dvojrozměrné nebo trojrozměrné prostředí, a nejsme tedy odkázáni na nenázorné textové údaje doplňující složité rovnice. Přesto, že je Femlab plný matematiky, může svoji úlohu s úspěchem a rychle vyřešit i člověk, který sice není v matematice příliš zběhlý, ale svému oboru dobře rozumí. Grafické výstupy jsou názorným vodítkem při hledání optimálního řešení zadané úlohy.

Femlab pracuje společně s Matlabem, využívá jeho pracovní prostor (workspace) a množství matematických a hlavně grafických funkcí. Pro chod Femlabu je instalace Matlabu nutná.

Práci ve Femlabu lze rozdělit do tří hlavních etap:

1. Volba typu úlohy, vytváření geometrického modelu řešené úlohy a zadání okrajových podmínek.
2. Generování sítě a řešení úlohy.
3. Zpracování výsledků řešení, data, grafy, animace.

Vstupní bránou do Femlabu je tzv. Model Navigator (obr. 1). Je to dialogové okno, ve kterém si zvolíme pracovní oblast – např. prostup tepla, proudění, pružnost a pevnost, elektromagnetismus atd. Touto volbou určujeme tvar PDR, které bude Femlab řešit, a zároveň tím dáváme pokyn, aby se v grafickém editoru zobrazovaly dialogy spojené s touto problematikou. Současně také volíme dimenzi geometrického prostoru, ve kterém budeme úlohu řešit. Může být jednorozměrný, dvojrozměrný nebo trojrozměrný.

Jestliže charakter naší úlohy nesplňuje žádná z již definovaných aplikací, máme možnost použít obecný tvar PDR a vhodným zadáváním koeficientů do těchto rovnic si vytváříme vlastní aplikaci. Femlab nabízí čtyři skupiny aplikací:

• fyzikální aplikace (vodiče stejnosměrného proudu, elektromagnetismus střídavého proudu, elektrostatika, magnetostatika, vodivost, difuze, nestlačitelné proudění, prostup tepla, pružnost a pevnost),

• režimy PDR (ve tvaru koeficientů, v obecném tvaru),

• slabé režimy (pro subdomény, hraniční body a omezené hranice),

• klasické PDR (Laplaceova, Helmholzova, Poissonova a Schrödingerova rovnice, rovnice vedení tepla a vlnění).

Po otevření grafického editoru, který lze v této fázi přirovnat k malému systému CAD, máme k dispozici všechny prostředky pro vytvoření geometrického modelu (obr. 2). Základem jsou geometrická primitiva, jako je kružnice, obdélník, bod, úsečka, oblouk, Bézierovy křivky, koule, kvádr, válec, kužel, krychle a ovoid. Složitější geometrické útvary pomáhají vytvářet tři booleovské operace, union, intersection a difference, a dále rotace nebo protlačování geometrických profilů. Geometrický model však můžeme do Femlabu také importovat: v ploše přes formát DXF, v trojrozměrném prostoru přes formát IGES.

Okrajové podmínky
Dalším krokem řešení je zadání okrajových podmínek. Nejedná se o nic jiného než o přiřazení určitých hodnot nebo funkčních závislostí obrysům tělesa. Pracujeme-li ve dvojrozměrném prostoru, definujeme okrajové podmínky pro body nebo hrany nakresleného geometrického útvaru. V trojrozměrném prostoru k těmto dvěma podmínkám přibude okrajová podmínka pro plochy. Každé geometrické těleso je Femlabem automaticky rozděleno na několik částí, což umožňuje zadávat různým částem modelu různé okrajové podmínky. Skládáme-li model z více částí, můžeme podmínky zadávat i uvnitř geometrického tělesa.

Vedle okrajových podmínek zadáváme také vlastnosti tzv. subdomény, která např. charakterizuje vlastnosti použitého materiálu uvnitř modelu. I zde platí možnost definovat pro různé vnitřní oblasti (plocha nebo objem) různé mechanické, tepelné nebo chemické vlastnosti. Můžeme tak řešit model složený z více druhů materiálů. Při práci s dialogovými okny (obr. 3) je zajištěna vzájemná vazba částí geometrického útvaru s údaji v dialogu. Současně se v dialogovém okně zobrazuje rovnice odpovídající zvolenému aplikačnímu režimu a symboly v editačních polích odpovídají zobrazené rovnici. Vždy tak máme kontrolu nad tím, kterou položku právě editujeme.

Generování sítě
V prvním kroku našeho řešení můžeme síť generovat automaticky. Vyžaduje-li to přesnost výpočtu, máme možnost síť dále zjemňovat nebo nastavením vhodných parametrů dále upravovat. Femlab také umožňuje generovanou síť statisticky vyhodnocovat, abychom mohli posoudit její kvalitu s dopadem na přesnost výpočtu. Při zobrazení sítě modelu v trojrozměrném prostoru máme možnost kontrolovat síť uvnitř tělesa nastavením vhodných řezů tělesem nebo nastavením jeho průhlednosti.

Řešení úlohy
Podle typu úlohy si Femlab zvolí řešič, jehož parametry můžeme dále upravovat. Nastavený řešič však můžeme kdykoliv změnit. Používají se tyto řešiče:

• lineární řešič pro řešení lineárních nebo linearizovaných PDR přímou metodou,

• nelineární řešič pro řešení nelineárních PDR,

• řešič pro řešení úloh v čase, lineární a nelineární PDR,

• řešič pro řešení úloh vedoucích na vlastní čísla, lineární a nelineární PDR,

• adaptivní řešič pro lineární a nelineární úlohy PDR, které jsou určeny malou plochou, s požadavkem na velkou přesnost (adaptivní generování sítě se pokouší identifikovat oblasti, ve kterých je třeba velká přesnost),

• iterační řešič lineární a nelineární PDR, a pro výpočet se používá iterační algoritmus zabudovaný uvnitř Femlabu,

• parametrický řešič (není to samostatný řešič, ale při jeho použití se využívá vlastnost řešičů pro lineární a nelineární PDR),

• vícesíťový řešič (multigrid) pro lineární a nelineární PDR ve tvaru koeficientů nebo v obecném tvaru (k řešení úlohy v jednoduchém geometrickém tělese používá lineární Lagrangeovy elementy; při řešení mohou být použity Lagrangeovy elementy od lineárního až do 5. řádu, Hermitovy 3. až 5. řádu a Argyrisovy 5. řádu).

S proměnnými pracuje Femlab místně, tzn. že jejich hodnoty se počítají v tom bodě, ve kterém jsou všechny potřebné informace. Můžeme však definovat i tzv. svázané proměnné (coupling variables), kterým lze pro výpočet přiřadit informace z jiného místa geometrického tělesa, nebo dokonce z jiného geometrického tělesa. Takové proměnné mohou být použity např. v koeficientech PDR nebo v průběhu vizualizace výsledků, tj. v následném zpracování. Při použití v koeficientech PDR jsou výsledkem nemístní závislosti (non-local dependencies), což nazýváme rozšířená multifyzikální úloha (extended multiphysics) – jako protiklad k obyčejné multifyzikální úloze s místními závislostmi v místních bodech.

Následné zpracování
Femlab nabízí zobrazení výsledků v různých řezech, plochách, obrysech, vektorech, proudnicích nebo izoplochách. V této části jsou plně využity grafické možnosti Matlabu. Vizualizovaná data lze zobrazovat ve dvojrozměrném nebo trojrozměrném prostoru, animace můžeme zobrazovat v klasickém okně Matlabu nebo ukládat do externího souboru. Femlab umožňuje zjišťovat grafické závislosti vypočtených hodnot v definovaných řezech nebo hodnoty v definovaných bodech. Používá k tomu grafické okno Matlabu.

V post-režimu můžeme využívat integraci nad oblastí nebo nad okrajovou podmínkou, abychom tak mohli např. získat průměrné hodnoty požadovaných veličin.

Kromě importu grafických souborů DXF a IGES můžeme do Femlabu načítat geometrická data zapsaná v M-souborech Matlabu. Celou modelovanou úlohu standardně zapisujeme do binárního MAT-souboru nebo v podobě znaků ASCII do M-souboru. Obsah tohoto souboru odpovídá datové struktuře FEM, která tvoří databázi informací o řešené úloze ve Femlabu. Její export do pracovního prostoru Matlabu nám umožní pracovat s jakýmikoliv daty: geometrickými tělesy, generovanou sítí nebo s výsledky řešení. Další možností je exportovat do pracovního prostoru Matlabu jenom vybrané části modelu, jako jsou části geometrických těles, okrajové podmínky, koeficienty PDR, síť, výsledky řešení nebo animace úlohy. Geometrická tělesa lze exportovat do souboru DXF nebo IGES, popř. do M-souboru.

Femlab může pracovat se Simulinkem, což je nadstavba Matlabu určená pro simulaci dynamických dějů. Při exportu se data ukládají do datové struktury, která je součástí bloku FEM v Simulinku. Model Femlabu se tak stává součástí schématu dynamického systému.

Femlab ke své činnosti potřebuje Matlab, protože využívá některé jeho funkce a je propojen s jeho pracovním prostorem. Femlab má k dispozici řadu svých vlastních funkcí, které lze rozdělit do tří hlavních skupin.

1. Funkce pro vytváření geometrických tvarů zkoumaného modelu. Jsou to funkce pro kreslení základních geometrických útvarů ve dvojrozměrném a trojrozměrném prostoru, funkce pro orientaci v prostoru a booleovské operace, které umožňují vytvářet složitější modely. Do této oblasti lze zařadit také funkce pro import a export modelů těles ve formátech DXF a IGES nebo využití binárních MAT-souborů a skriptových M-souborů Matlabu.

2. Funkce související s řešením úlohy. Do této oblasti patří funkce související s přípravou modelu pro jeho řešení, jako je zadávání okrajových podmínek, generování geometrické sítě potřebné k diskretizaci úlohy a nakonec funkce – řešiče, které provádějí vlastní numerický výpočet.

3. Funkce určené ke zpracovaní výsledků. Jsou zaměřeny na grafické zobrazování vypočtených veličin např. v řezech nebo v izoplochách, je možná animace úlohy jako movie v Matlabu nebo zápis výsledků do externího souboru.

Všechny uvedené funkce jsou přístupné dvěma způsoby: zaprvé z grafického rozhraní, které se podobá běžnému prostředí CAD, kdy uživatel pracuje s myší a používá k práci ikony nebo roletová menu, a zadruhé z příkazového řádku Matlabu, do kterého zkušenější uživatel zadává model, vstupní hodnoty, okrajové podmínky, geometrickou síť, parametry řešiče nebo zobrazování výsledků.

Uvedli jsme, že volbou aplikačního režimu a vstupem do Femlabu zároveň zvolíme tvar PDR, které řeší úlohu jako aplikaci. Avšak jádro Femlabu řeší jenom některé úlohy z oblasti pružnosti a pevnosti, fyziky, chemie a elektromagnetismu. Pro specialisty ve vybraných oborech jsou určeny přídavné moduly Femlabu, ve kterých jsou tvary PDR přizpůsobeny některým speciálním úlohám. V současné době jsou přídavné moduly zaměřeny na profese v oborech pružnosti a pevnosti – Structural Mechanics Module, průmyslové chemie – Chemical Engineering Module a elektromagnetismu – Electromagnetics Module. Všechny moduly obsahují aplikace ve dvojrozměrném a trojrozměrném prostoru.

Structural Mechanic Module – SMM
Femlab sám o sobě z oboru pružnosti a pevnosti řeší rovinnou deformaci, rovinnou napjatost a běžné zatížení těles. Modul SMM tyto oblasti rozšiřuje v ploše o osově symetrickou deformaci, prostup tepla, Eulerovy nosníky a zatížení desek. V trojrozměrném prostoru je navíc možnost výpočtu skořepin.

Chemical Engineering Module – Chem
Přestože název modulu může přitahovat zájemce hlavně z oboru chemie nebo průmyslové chemie, určitě uspokojí i uživatele v jiných profesích. Modul obsahuje aplikace použitelné například při řešení úloh v oblasti stlačitelného a nestlačitelného proudění. Modul Chem ve svém principu vychází z momentové, energetické a hmotnostní rovnováhy. Tomu odpovídá proudění tekutin (Navierova-Stokesova rovnice, Eulerova rovnice), nenewtonovské proudění, Brinkmanova rovnice, Darcyho pravidlo, prostup a vedení tepla, difuze, difuze a konvekce, Maxwellova-Stefanova difuze a konvekce a definice Nerstova-Planckova elektrokinetického modelu.

Electromagnetics Module – EM
Možnosti výpočtu elektrostatiky a magnetostatiky v samotném Femlabu rozšiřuje EM o řadu aplikací ve dvojrozměrném a trojrozměrném prostoru, jako jsou např. statická a elektrická pole, proudy v rovině, kolmé proudy v rovině a ve vodičích, šíření vln v rovině a ve vodičích a elektomagnetické vlnění.

V tomto příkladu se řeší proudění ohřívané vody ve výměníku tepla. Voda proudí ze spodní části výměníku kolem topných těles směrem nahoru. Trojrozměrnou úlohu zjednodušíme na úlohu v ploše podle obr. 4.

Na vytvářený model se zaměříme jako na úlohu, v níž je třeba zároveň řešit prostup tepla a proudění vody popsané Navierovou-Stokesovou rovnicí pro proudění nestlačitelných tekutin. Řešíme tedy multifyzikální úlohu se čtyřmi neznámými proměnnými: složky rychlosti proudění ve směru x a y, tlaku a teploty.

K řešení použijeme Navierovu-Stokesovu rovnici ve tvaru:

kde F je objemová síla, r hustota vody, h dynamická viskozita, u vektor rychlosti proudění.

Tepelná rovnice je dána vztahem

kde c_p je měrná tepelná kapacita za stálého tlaku, T teplota, Q teplo zdroje.

V Model navigatoru zvolíme multifyzikální aplikace Heat transfer a Incompressible Navier-Stokes (obr. 5). Zvolením obou režimů zároveň vybereme uvedené PDR, které bude Femlab řešit. U obou režimů jsou zobrazeny závisle proměnné veličiny (T, u a p)¹⁾. Implicitní konečné prvky pro řešení proudění jsou smíšené trojúhelníkové Lagrangeovy elementy prvního a druhého řádu, pro aplikaci prostupu tepla jsou to Lagrangeovy elementy druhého řádu. V dalším kroku vytvoříme jednoduché geometrické těleso znázorňující část prostředí, ve kterém proudí ohřívaná voda. Použijeme obdélník o rozměru 0,005 × 0,04 m. Kruhový výřez znázorňující polovinu trubky výměníku vytvoříme odečtením kružnice o poloměru 0,0025 m použitím booleovské operace difference.

V následujícím kroku zadáme okrajové podmínky pro každý z aplikačních režimů zvlášť. Volbou režimu se zároveň mění obsah dialogového okna pro zadání těchto podmínek. V režimu prostupu tepla zadáváme teplotu vody na vstupu (ve spodní části tělesa) a teplotu trubky, která vodu vyhřívá (půlkruhový výřez v obdélníku). Ostatní hrany tělesa označíme jako symetrické. Dále zadáme okrajové podmínky pro režim proudění vody. Na vstupu vody do ohřívače zadáme rychlost proudící vody a na výstupu položíme tlak roven nule. V okolí kruhového průřezu geometrie zadáme nulové složky rychlosti proudící vody (v_x = 0, v_y = 0), ostatní části průřezu označíme opět jako symetrické.

Volbou subdomény zadáváme fyzikální vlastnosti proudícího média. Z hlediska prostupu tepla je třeba znát hustotu vody r, její měrnou teplenou kapacitu c_p a tepelnou vodivost k_c. Jako tepelný zdroj uvedeme výraz pro šíření tepelného toku, protože v modelu není rozložený tepelný zdroj:

Hodnoty derivací T podle x a y jsou ve Femlabu automaticky k dispozici jako vnitřní proměnná. Pro řešení úlohy v čase zadáme hodnotu závisle proměnné T v čase t₀.

Zvolením režimu proudění podle Navierovy-Stokesovy rovnice získáme dialog vyžadující zadat hustotu proudící vody r a její dynamickou viskozitu h.

Objemová síla se projeví v souvislosti se změnou hustoty vody vlivem teploty. Objemovou sílu ve směru osy y vyjádříme výrazem: F = agr (T – T₀)

kde a je teplotní roztažnost vody, g tíhové zrychlení, T₀ teplota vody na vstupu do ohřívače.

Pro řešení úlohy v čase zadáme hodnotu vstupní rychlosti u_y(0) v čase t = 0 (obr. 6).

Předchozími úkony jsme připravili úlohu k řešení – viz první etapa pracovního postupu. V následujícím kroku úlohu vyřešíme. Vygenerujeme síť bez jakýchkoliv dalších úprav (obr. 7) a spustíme řešení. Využijeme přitom implicitně nastavený nelineární řešič. Výsledkem je barevné znázornění teploty v proudící vodě.

Chceme však také znát, jakým způsobem voda ve výměníku proudí. Nastavením parametrů v dialogu Post (jako postprocesor) zvolíme vykreslení složek rychlosti ve směru os x a y pomocí vektorů – položka Arrow. Výsledkem je současně znázornění proudění vody a rozložení teploty (obr. 8). Kromě teploty máme možnost znázorňovat teplotní gradienty a tepelné toky, složky vektoru rychlosti u, tlak p, rychlostní pole a vířivost. Pomocí řezů můžeme zobrazit např. graf profilu teploty T v libovolném řezu modelem. Abychom zjistili průměrnou teplotu vody na výstupu z výměníku, můžeme využít další funkci Femlabu: integraci hranice nebo subdomény. V tomto příkladu použijeme integraci hranice na výstupu z výměníku. Průměrná teplota T_p je dána podílem dvou integrálů:

Využitím dialogu pro integraci hranic postupně zjistíme hodnoty jednotlivých integrálů a zápisem výsledků do pracovního prostoru Matlabu vypočteme průměrnou teplotu ohřívané vody.

V poslední části uvedeného příkladu vyřešíme úlohu v závislosti na čase, abychom zjistili, za jak dlouho po zapnutí ohřívač dosáhne stabilního stavu. Předpokladem tedy bude zapnutí ohřívače, ve kterém již proudí voda. Nejprve samostatně vypočteme proudění vody (Navierovou-Stokesovou rovnicí), což bude sloužit jako počáteční podmínka pro další výpočet závislý na čase. V dialogu pro nastavení parametrů zvolíme pouze režim Incompressible Navier-Stokes a spustíme řešení.

Po skončení výpočtu nastavíme časově závislý řešič (Time dependent solver) a v jeho parametrech zadáme rozsah časů od 0 do 30 sekund s krokem jedna sekunda. Dále je nutné do řešení zahrnout oba fyzikální režimy, tj. prostup tepla a proudění. Femlab umožňuje použít pro další výpočet předchozí řešení jako počáteční podmínku v tzv. režimu Restart. V postprocesoru si potom můžeme zobrazit řešení v každé sekundě výpočtu. Na obr. 9 je viditelné měnící se rozložení teploty po šířce tělesa v jednom průřezu těsně nad ohřívající trubkou v prvních čtrnácti sekundách ohřevu, kdy se teplotní rozložení stabilizuje.

Na druhou otázku, za jak dlouho ohřívač dosáhne nového stabilního stavu, si můžeme odpovědět nakreslením grafu, ve kterém zobrazujeme teploty v jednotlivých bodech vedeného řezu. Grafem je potom závislost teploty v každém bodě na čase. Z výsledků je zřejmé, že ve vybraném řezu se ohřev stabilizuje asi po deseti sekundách.

Hlavní výhodou Femlabu je jeho všestranná použitelnost v různých profesních oborech a oblastech fyziky, přičemž je možné kombinovat několik fyzikálních pohledů v jediné úloze. Jednoduchý příklad uvedený v tomto článku může posloužit jako první krok k vytvoření představy o možnostech Femlabu. Nejsou-li čtenáři proudění a prostup tepla právě blízké, CD s Femlabem obsahuje knihovnu podrobných příkladů s více než osmdesáti dalšími řešenými úlohami z oblasti akustiky, elektromagnetismu, pružnosti a pevnosti, polovodičů, geofyziky, šíření elektromagnetických vln atd.

Co je další vlastností Femlabu? Existuje mnoho softwarových balíků pro řešení úloh metodou konečných prvků, které jsou úzce profesně zaměřeny a jejich ceny mohou být několikanásobně vyšší. Femlab je určen pro obecné použití s možností vytvořit vlastní uživatelskou aplikaci „šitou na míru“.

Program pracuje společně s Matlabem, proto jsou jeho požadavky na operační systém stejné. Velikost operační paměti se doporučuje od 512 MB do 2 GB podle typu řešené úlohy.

¹⁾ Pozn. red.: V Model navigatoru – viz obr. 5 – jsou složky rychlosti u označeny takto: u_x = u, u_y = v. V článku se budeme držet obvyklého označení složek vektoru indexy.

Karel Bittner,
Humusoft s. r. o.