Generator mreže za modeliranje strjevanja teles zapletenih oblik z dvojno recipročno metodo robnih elementov
Mesh Generator for Modelling Complex Bodies Solidification by the Dual Reciprocity Boundary Element Method
Košir A\ B. Šarler, LFDT, FS, Univerza v Ljubljani
l/ članku Izvirno predlagamo, kako z optlmlzacijsko metodo konstruirati koordinate nestrukturirane mreže točk, primernih za diskretizaoijo zapletenih območij. Uspešnost diskretizaoije preverimo na modelu strjevanja neskončnega vogala z numerično metodo robnih elementov z dvojno recipročnostjo.
Ključne besede: generator mreže, optimizacija, ulivanje, strjevanje, metoda robnih elementov, metoda dvojne recipročnosti, Štefanova naloga, prenos toplote in snovi
An original method for optimal unstructured mesh generation, suitable for complex shape discretization, is proposed in the paper. Efficiency of discretization is checked numerically by the boundary element method vvith dual reciprocity on the infinite corner solidification model.
Key vvords: mesh generation, optimization, casting, boundary element method, dual reciprocity method, Štefan problem, heat and mass transfer
1. Uvod
Z uporabo dvojno recipročne metode (DRM) robnih elementov (BEM) smo doslej uspešno numerično izračunali testne primere strjevanja teles preprostih oblik1 in se prepričali, daje to učinkovita in natančna numerična metoda, primerna tudi za tovrstno družino močno nelinearnih problemov.
Strjevanje kovinskih ulitkov fizikalno najpreprosteje opišemo z ohranitveno enačbo energije v idealni kontinuumski mešanici trdne in tekoče faze v izbrani snovi. Izbrani fizikalni model omogoča upoštevati temperaturno odvisno toplotno prevodnost in specifično toplotno kapaciteto snovi, ki se strjuje, vse tri vrste običajno zanimivih robnih pogojev in poljubne procesne parametre ulivanja. Ker so v industrijski praksi običajno bolj zanimiva telesa zapletenih oblik, je ključna ustrezna učinkovita diskretizacija območja in njegovega roba.
Metoda robnih elementov omogoča v primerjavi z ostalimi klasičnimi numeričnimi metodami večjo svobodo pri diskretizaciji območja. Lege mrežnih točk niso nujno razporejene ortogonalno in enakomerno, pač pa povsem svobodno, običajno bolj zgoščeno na območjih, kjer pričakujemo večjo dinamiko procesov. Čeprav je mreža točk lahko poljubna, nekatere diskretizacije zmanjšajo velikost napake pri enakem številu mrežnih točk.
Namen te raziskave je konstruirati z ustrezno minimizacijo izbranih ciljnih funkcij generator za nestrukturirane mreže.
1 Ale5 KOŠIR. dipl. inž. liz.
Laboratorij 7.a dinam. fluidov in termodin.. FS. Univerza v Ljubljani Aškerčeva 6. 61000 Ljubljana
2. Ciljne funkci je
Kartezične koordinate ni točk v n razsežnem prostoru označimo z .v/, j = 1.....m,i = 1.....n.Z izrazom
r
-xk i i
7=1
min, V
j* k
xi.-xk)2
(D
definiramo indeksni vektor jk, k = 1.....m. ki za vsako točko k
pove indeks njene najbližje sosede.
Informacijo o obliki območja, na katerem generiramo mrežo, zajema funkcija
Award( a ) =
|l. s'eS!ur.
(2)
kjer množica Q predstavlja odprto povezano območje in T njegov rob. Množica ii naj bo zaprtje območja, fi = £2 u f.
Definirajmo nekaj ciljnih funkcij in si oglejmo njihove lastnosti. Najpreprostejša je prva ciljna funkcija
ObjFunc,(x) =
1
A = i Award( xk )
(3)
Funkcija je neomejena, če katera od mrežnih točk leži zunaj območja, sicer pa doseže svoj minimum m pri katerikoli porazdelitvi točk po območju in njegovem robu. Ker ciljna funkcija upošteva le obliko območja in ne medsebojnih leg mrežnih točk, je zanimiva le primerjalno z ostalimi ciljnimi funkcijami.
Definicija druge ciljne funkcije je
m 111
ObjFunc;(x)= XX =i/=i j* k
A = l,=l\ X" ,(A/-xA)2Award(A'A)
\ ^^ l—l I	I
(4)
Ta funkcija je do multiplikativne konstante enaka potencialni energiji v območje zaprtih enako nabitih točkastih nabojev. Ker se naboji medsebojno odbijajo, je njihova porazdelitev gostejša na robu območja.
Kot prejšnja ciljna funkcija tudi tretja
odvodih in z uporabo Fletcher-Reeves-Polak-Ribierove metode konjugiranega gradienta1" učinkovito poiščemo njen globalni minimum.
Čeprav tretja ciljna funkcija znotraj območja ni povsod gladka, smo jo poskusili tako kot drugo funkcijo optimizirati z omenjeno metodo konjugiranega gradienta. Navedimo njen parcialni odvod po /-ti koordinati p-le točke znotraj območja
50bjFunc,(A) _
ObjFunc,(x ) = - f] min \ ( x ' -\k )2 Award(xA ) (5)
j*k
pri minimizaciji sili mrežne točke vsaksebi, a ima to zanimivo in uporabno lastnost, da doseže na n razsežni kocki svoj minimum, če je m točk enakomerno porazdeljenih in če velja m'" e N.
3.	Vsebnostni test
Opišimo, kako smo presodili, ali je testirana točka vsebovana v zaprtju območja, .v* e 11. Ta postopek imenujemo vsebnostni test.
Območje v dvorazsežnem prostoru opišemo s kartezičnimi koordinatami oglišč poligona, ki se območju najbolj prilega. Najprej preskusimo, če leži točka v* znotraj minimalnega pra-vokotnika. ki vsebuje vse območje. Nato z vsebnostnim testom presodimo, ali točka ne leži zunaj poligona.
Stabilni numerični vsebnostni testi so v ravnini /.aradi končne aritmetike še danes trd oreh numerične analize". Primerjali smo
•	paniostni vsebnostni test, pri katerem tvorimo daljico med testirano točko in znano točko zunaj poligona ter izračunamo pamost števila presečišč daljice s stranicami poligona'1,
•	vsebnostni test : barvanjem rasterizirane notranjosti poligona,
•	vsebnostni test r določanjem orientacije trikotnikov, ki jih tvorimo s testno točko in sosednjima ogliščema poligona4.
•	vsebnostni test : ovijalnim številom, ki meri. v kateri smeri poligon ovije testirano točko.
Parnostni vsebnostni test znatno pospeši, če namesto znane zunanje točke uporabimo točko v neskončnosti.
Ko je območ je rasterizirano. je rastrski test zelo hiter, vendar nenatančen, če se rob poligona ne ujema /. rastrom. Najbolj zanesljive rezultate smo dosegli z ovijalnim vsebnostnim testom, kije pri danem poligonu hkrati hitrejši od parnostnega in orientacijskega testa. V več kot dvorazsežnem prostoru uporabljamo druge metode.
4.	Optimizacijske metode za cil jne funkcije
Prostor parametrov, v katerem iščemo minimum ciljne funkcije, je n ■ m razsežen, zato je zelo pomembno uporabiti učinkovito optimizacijsko metodo, ki se bo minimu funkcije pri velikem številu mrežnih točk karseda hitro in poljubno približala. Za vsako izmed funkcij si oglejmo, kako smo iskali njen globalni minimum.
Prva ciljna funkcija doseže svoj minimum pri poljubni (sto-hastični) porazdelitvi mrežnih točk po območju. Ce je bila neka mrežna točka od roba območja oddaljena za manj kot ~im, smo jo prestavili v najbližjo točko na robu, sicer bi pri rešitvenem postopku v fizikalnem modelu primanjkovalo robnih pogojev.
Druga ciljna funkcija je znotraj območja fi gladka funkcija koordinat, zato izkoristimo informacijo o njenih parcialnih
ObjFunc;,( a )
A""-v' p p

(6)
-I
xh-xk

djj
To minimizacijo smo primerjali s posebej konstruiranim, na Boltzmannovcm simuliranem ohlajanju temelječim postopkom, ki statistično zagotavlja konvergenco h globalnemu minimumu funkcije. V prvi fazi postopka po območju stohastično posejemo mrežne točke in jim predpišemo skupno umetno temperaturo kot mero fluktuacij njihovih leg v koordinatnem prostoru. V korakih nižamo temperaturo, pri tem premikamo lege točk in spremenjena stan ja sprejmemo ali zavržemo glede na spremembo vrednosti ciljne funkcije. Če je nova vrednost ciljne funkcije nižja kot stara, potem novo stanje brezpogojno sprejmemo, sicer pa le sorazmerno z Boltzmannovo verjetnostjo, odvisno od vrednosti umetno predpisane temperature. Od zahtevane kvalitete rezultatov je odvisno, kako dolgo postopek simuliranega ohlajanja ponavljamo. Točnemu globalnemu minimumu se sčasoma poljubno približamo.
5. Str jevanje neskončnega vogala. Re/.ultati generatorja mreže
Ocenimo najprej učinkovitost optimizacijske metode. Optimizacija prvih dveh funkcij ni problematična. Prva funkcija doseže svoj minimum pri katerikoli porazdelitvi točk po Cl. Minimumu druge funkcije se pri m točkah približamo na relativno natančnost koordinat 10° v približno 0.1 ■ m zaporednih korakih metode konjugiranega gradienta.
Izkaže sc. da je minimalna vrednost tret je ciljne funkcije na »-razsežni kocki pri pogoju m'" e N enaka
ObjFunc, (a ) =
-I
m -1)
(7)
Z metodo konjugiranega gradienta smo se tej vrednosti pri poljubni začetni porazdelitvi m točk pri m> 10 približali, vendar je nismo dosegli, ker optimizacijski algoritem ni našel globalnega minimuma. Z metodo simuliranega ohlajanja smo se uspešneje približali globalnemu minimumu, kar kaže slika 1.
Generator mreže smo preskusili na modelu strjevanja neskončnega pravokotnega vogala. Model je bolj obširno opisan v ref. 1,2. Neskončno območje smo aproksimirali s končnim kvadratom s stranico, dolgo 1,5 enote. Ob začetnem času je vogal pri konstantni začetni temperaturi, višji od temperature tališča. Po začetnem času temperatura na robu območja skokovito pade na temperaturo, nižjo od evtektične, in vogal začne zmrzo-vati.
-10"
-10-
1.	stohastična A
2.	stohastična x konjugirani gradient 0 simulirano ohlajanje +
enačba (7) o

A

A

A
A
A '
-1
6 7 8 9 10 število mrežnih točk
16
Slika 1: Primerjava dveh optimizacijskih metod za minimizacijo na enotskem kvadratu. Najprej so predstavljene vrednosti dveh stohastičnih porazdelitev točk. kot smo jih uporabili za začetni približek v optimizacijskem postopku. Nato so predstavljene minimalne vrednosti, dobljene z metodo konjugiranega gradienta, z metodo simuliranega ohlajanja in z enačbo (7). Slika kaže. da se pri številu mrežnih točk, enakemu kvadratu celega števila, s konjugiranim gradientom in simuliranim ohlajanjem dobro približamo teoretičnemu
globalnemu minimu ciljne funkcije Figure 1: A comparison of two optimization methods for minimizing on a unit square. First, values were given for two stochasticaIly
distributed points vvhich vvere used as starting points for an optimization algorithm. Then the minimal values given by the conjugated gradient method, the simulated annealing method and by equation (7) vvere added. The figure shovvs that if the number of grid points is equal to the square of an integer. the theoretical minimum value of the object function can be closely approximated vvith the conjugated gradient and simulated annealing methods
V	delu2 smo raziskovali, kako na temperaturno polje vpliva širina talilnega intervala in kako se vede numerična metoda pri različnih vrednostih Štefanovega števila in pri različnih diskretizacijah območja. Numerično izračunano temperaturno polje v mrežnih točkah smo primerjali s semianalitičnimi vrednostmi. dobljenimi s standardno referenčno metodo,:.
V	pričujočem delu smo te raziskave nadaljevali s poudarkom na optimalni porazdelitvi danega števila mrežnih točk. tako daje napaka temperaturnega polja kar najmanjša. Podrobneje smo primerjali povprečno in maksimalno vrednost napake temperaturnega polja, dobljenega z novimi diskretizacijami območja, z napakami, izračunanimi z enakomerno, neenakomerno in premaknjeno diskretizacijo območja1. Diskretizacije so shematično predstavljene na sliki 2 z značilnimi mrežami, ki jih tvorimo s šestintridesetimi mrežnimi točkami.
Ugotovili smo, da je najmanj primerna stohastična mreža. Pri njej so napake temperaturnega polja povprečno za velikostni red večje od napak temperaturnega polja na mrežah, dobljenih z minimizacijo ostalih dveh ciljnih funkcij. Napake temperaturnega polja na mrežah, dobljenih z minimizacijo druge in tretje ciljne funkcije, so si po velikosti podobne. V izračunanih primerih se napake razlikujejo za največ 10 %, pri čemer so pri drugi ciljni funkciji vedno večje.
Če je število mrežnih točk kvadrat celega števila, smo z minimizacijo tretje ciljne funkcije dobili po kvadratu enako porazdeljene točke in dosegli enako napako, kot pri enakomerni mreži, cf.l, s čimer smo preverili uspešnost minimizacijskega postopka.
Porazdelitev točk po območju 1.5
1
0.5 0
0.5
1.5
(a)
Porazdelitev točk po območju 1.5
1
0.5 0
0.5 1 1.5
(c)
Porazdelitev točk po območju 1.5
1
0.5 0
0.5 1 1.5
(d)
Porazdelitev točk po območju 1.5
1
0.5 0
0.5
1.5
(b)
Porazdelitev točk po območju
1.5
1
0.5 0
0.5 1 1.5
(e)
Porazdelitev točk po območju 1.5
1
0.5 0
0.5 1 1.5
(e)
Slika 2: Primerjava porazdelitev mrežnih točk, kot jih da minimizacija prve (a), druge (c) in tretje (d) ciljne funkcije, z enakomerno (b), neenakomerno (č) in premaknjeno (e) mrežo Figure 2: Comparison of grid point distribution given by the minimization of the first (a), second (c) and third (d) object function. vvith equidistant (b). non-equidistant (č) and displaced (e) grids
Slika 3: Območje v obliki fraktalne Kochove snežinke smo z minimizacijo druge ciljne funkcije diskretizirali s 36 točkami. Ni znano, ali obstaja za poljubno število mrežnih točk na tem območju ena sama porazdelitev točk. ki da globalni minimum druge oziroma tretje ciljne funkcije, vsekakor pa je v njegovi okolici veliko lokalnih minimumov, od katerih smo enega dosegli. S slike razberemo, da je
notranjih točk 8 in robnih 28 Figure 3: Domain in the form of a fractal Koch snovvflake vvhich vvas discretized over 36 points vvith the second object function. It is not knovvn vvhether there is a unique distribution for an arbitrary number of grid points vvhich gives a global minimum value for the second or
third object functions, hovvever many local minima exist in their vicinity, vvhich can be used to approximate the global minimum. The figure indicates 8 interior and 28 boundary points
6. Zaključek
Po dosedanjih preskusih ocenjujemo, da je predlagana optimizacijska tehnika s simuliranim ohlajanjem primerna in učinkovita metoda za generacijo optimalne diskretizacije poljubnega povezanega območja, na katerem rešujemo nelinearni problem s taljenjem in strjevanjem z metodo robnih elementov z dvojno recipročnostjo. Kot primer uspešne realizacije mreže naj navedemo diskretizacijo Kochove fraktalne snežinke na sliki 3.
Čeprav smo v dosedanji raziskavi preskusili več mini-mizacijskih tehnik, je pri velikem številu mrežnih točk (/;;>104) minimizacija še vedno dolgotrajna. Nadaljevanje raziskav bo zato kazalo predvsem v smer, kako pospešiti minimizacijo.
Zahvala
Predstavljena raziskava je del temeljnih raziskav pri projektu Dvofluidno modeliranje sistemov s trdno-kapljevinskimi faznimi prehodi. Avtorja se za podporo zahvaljujeta Ministrstvu za znanost in tehnologijo Republike Slovenije in Mednarodnemu biroju Raziskovalnega centra Jiilich, Nemčija.
Literatura
Sarler, B., Košir. A., Kuhn, G.: Dual reciprocitv boundary element method for Štefan problems. (v recenziji) Int../. Numer. Methods En<>. - Sarler, B., Košir, A.: Solution of melting and solidification problems by the dual reciprocitv boundary element method, Lewis. R. W.: Numerical Methods in Thermal Problems, Pineridge Press, Svvansea 1993, 139-150
Košir, A.. Šarier, B.: Influence of mesh on dual reciprocity boundary element method for Štefan problems, GAMM Abstracts 1994 4 Knuth. D. E.: Axioms and hulls. Springer, Berlin 1992 J Arvo, J.: Graphics gems. Academic Press. New York 1992 " Harrington, S.: Computer graphics, McGraw-Hill, New York 19X7 Fletcher, R., Powell, M. J. D.: A rapidly convergent method for mini-mization, Comp. J. 6. 1963, 163-168
11	Fletcher, R.. Reeves, C. M.: Function minimization by conjugate gra-dients, Comp..1. 7. 1964, 149-154
'' Polak, E.: Computational methods in optimization, Academic Press, New York 1971
10 Shewchuk, J. R.: An introduction to the conjugate gradient method vvithout the agonizing pain, Carnegie Mellon University report 1994 1 Szu, H„ Hartley, R.: Fast simulated annealing. Plivs. Lett. A 122, 1987, 157-162
12	Rathjen, K. A., Jiji, L. M.: Heat conduction with melting or freezing in a corner../. HeatTransfer9Z. 1971, 101-109