3. Szimulációs vizsgálatok
3.1. A folyamatos izoterm kristályosítók dinamikus folyamatai
3.1.1. Az MSMPR izoterm kristályosító dinamikai vizsgálata
Az MSMPR izoterm kristályosítóban kialakuló kristályok méreteloszlása a (2.1.1.12) vagy a (2.1.3.1) egyenlettel megadott populációs mérlegegyenlettel írható le a h
szelekciós függvény alkalmazásával. A két populációs mérlegegyenlet a gócképződési forrástag kezelésében különbözik egymástól, amelyet az alkalmazott matematikai megoldási módszerek tesznek szükségessé, de a két egyenletrendszer egyenértékűnek tekinthető. A (2.1.1.12) egyenlettel megadott populációs mérlegegyenletet a momentumok módszerével oldottam meg, amellyel a parciális differenciálegyenletet közönséges differenciálegyenlet-rendszerré transzformáljuk. A momentum transzformáció után kapott dimenziómentesített momentum modellt a (2.1.2.11-2.1.2.17) egyenletek írják le. A differenciálegyenlet-rendszert a MATLAB 5.3 szoftver ode23tb közönséges differenciálegyenlet rendszer megoldójával oldottam meg. A differenciálegyenlet rendszer megoldásával a méreteloszlás momentumainak, valamint az oldott komponens és az oldószer koncentráció időbeli változását kapjuk meg. A 2.1.3 és 2.1.4. fejezetben bemutatott adaptív végeselemű ortogonális kollokáció alkalmazásával a (2.1.3.1) egyenletből származtatott (2.1.3.9) dimenziómentesített populációs mérlegegyenletet oldottam meg. A kollokációs módszer alkalmazásával a méreteloszlás értékeit nyerjük a kollokációs pontokban, amelyek egyben a szakaszonkénti közelítő függvények együtthatói. Az együtthatókra felírt (2.1.3.40-2.1.3.43) differenciál-algebrai egyenletrendszer kiegészül a (2.1.3.17-2.1.3.21) egyenletek által megadott oldott komponens és oldószer tömegmérleggel, valamint a (2.1.3.22-2.1.3.23) kinetikai összefüggésekkel. A differenciál-algebrai egyenletrendszer megoldására a DASSL programkódot (Petzold, 1983) használtam.
( )
L =1A (2.1.2.11-2.1.2.17) momentum modell alkalmazásával a kinetikai paraméterek hatását vizsgáltam az MSMPR izoterm kristályosító dinamikai viselkedésére az , és paramétereken keresztül. A modell további paramétereit, valamint a kezdeti feltételeket a 3.1. táblázat tartalmazza. A szimulációs vizsgálatokban feltételeztem, hogy a kristályosítóba tiszta anyalúg lép be és kezdetben a kristályosító nem tartalmaz kristályokat. A kristályosító működése főbb jellemzőinek az
F Da g
L átlagos kristályméret, az anyalúg koncentráció, valamint az momentum időbeli változása tekinthető. Az
y x3 L
c
átlagos kristályméret számításához a megoldásként kapott , i=0,1,2,3 momentumokat áttranszformáltam a µ , i=0,1,2,3 momentumokra majd az átlagos méretet a µ kifejezéssel adtam meg. Az így kapott átlagos kristályméret az L méretváltozó mentén értelmezett átlagos méretet adja meg. A transzformációt a 2.1.2.
fejezetben megadott s , i=0,1,2,3 tényezők segítségével hajtjuk végre. Az s , i=0,1,2,3 tényezők meghatározásához az s
xi i
0 1/µ
i i
c, kg és kv paraméterek megadása szükséges. Az s tényezőt és egyben a c telítési koncentrációt a ζ dimenziómentes paraméter definíciója alapján határoztam meg a c
s c
in betáplálási koncentráció illetve a ρc
kristálysűrűség specifikációja után:
( )
−1
= −
c c in c s
c c
ζ ρ
ζ és
s in
c c c
s = 1−
(3.1.1) A kg és kv paraméterek értéke szabadon választható, a számításoknál a 3.1.
táblázatban megadott értékeket használtam.
3.1. táblázat. A (2.1.2.11-2.1.2.17) momentum modell paramétereinek megadása.
a 0 st 1/τ x3(0) 0 cin 211.475 kg/m3
ζc 20.0 x0(0) 0 y(0) 0.5 kg 1.0×10-5 (m/h)⋅(m3/kg)-g yin 1.0 x1(0) 0 kv 0.47
τ 2.5 h x2(0) 0 ρc 1760 kg/m3
A 3.1.1.1.-3.1.1.3. ábrák az MSMPR izoterm kristályosító dinamikai viselkedését mutatják különböző értékek mellett. Az F paraméter növekvő értéke mellett a kristályosító fokozatosan instabil állapotba kerül és csillapítatlan oszcilláció kialakulása figyelhető meg mind az átlagméretben, koncentrációban és az x momentumban. A gócképződés sebessége az paraméter értékének csökkenésével növekszik és nagyobb számú, de kisebb átlagméretű kristályok kialakulását vonja maga után, amely megfigyelhető a 3.1.1.1. és 3.1.1.3. ábrán. A kristályosító hozama, amely a harmadik momentummal arányos, az F paraméter csökkenésével növekszik. Ezzel ellentétes tendencia figyelhető meg az átlagos kristályméret értékében, amely egy nagyságrendnyi eltérést mutat a vizsgált F paraméterek esetében. A 3.1.1.4.-3.1.1.6. ábrákon a különböző paraméterek esetén látható az átlagméret, koncentráció és a harmadik momentum időbeli változása. Csökkenő paraméter mellett a kristályosító a stabil működési állapot felé tolódik el és csillapodó oszcilláció figyelhető meg a kristályosító működésében. Az átlagos kristályméret és a kristályosító hozamának nagyságrendje nem változik. A g hatványkitevő értékének növekedésével a kristályosító működése stabilizálódik (3.1.1.7.-3.1.1.9. ábrák). A kristályosító hozama nem változik jelentős mértékben, azonban az átlagos kristályméret a g paraméter növekedésével jelentősen megnövekszik, amely a megnövekedett kristálynövekedési sebességnek köszönhető.
F
3
F
Da
Da
A dinamikai vizsgálatok második csoportjában az ortogonális kollokáció módszerét használtam a (2.1.3.9-2.1.3.23) egyenletekkel megadott modell megoldására h
szelekciós függvény esetében. Az ortogonális kollokáció alkalmazásának eredményeként a momentum modell megoldásával kapott eredmények mellett megkapjuk a képződő kristályok méreteloszlásának dinamikáját is. Abból a célból, hogy a két matematikai módszer által kapott eredmények összehasonlíthatók legyenek a (2.1.3.9-2.1.3.23) egyenletrendszer paramétereit úgy adom meg, hogy azok megfeleljenek a momentum modell paramétereinek. Összehasonlítva a (2.1.2.11-2.1.2.17) és (2.1.3.9-2.1.3.23) modellt látható, hogy az oldott komponens tömegmérlege és az oldószer mérleg megegyezik a két esetben, így a momentum modellben szereplő g, ζ
( )
L =1c, a és yin paraméterek valamint az y(0) kezdeti feltétel megegyeznek a (2.1.3.9-2.1.3.23) modell paramétereivel. Hasonlóan az F paraméter is, amely ebben ez esetben a (2.1.3.22) gócképződési kinetikai egyenlet paramétere. A populációs mérlegegyenlet eltérő kezeléséből adódóan azonban a (2.1.3.22-2.1.3.23) kinetikai egyenletekben szereplő kp és kg sebességi állandókat, valamint a populációs mérlegegyenletben szereplő α paraméter értékét meg kell adni, illetve a momentum modell megfelelő paramétereiből kell meghatározni. A kp gócképződési sebességi állandó, valamint az α
3.1.1.1. ábra. Az átlagos kristályméret időbeli változása (Da=500,g=1,ζc=20,a=0).
3.1.1.2. ábra. Az anyalúg koncentráció időbeli változása (Da=500,g=1,ζc=20,a=0).
3.1.1.3. ábra. Az x3 momentum időbeli változása (Da=500,g=1,ζc=20,a=0).
3.1.1.4. ábra. Az átlagos kristályméret időbeli változása (F=2,g=1,ζc=20,a=0).
3.1.1.5. ábra. Az anyalúg koncentráció időbeli változása (F=2,g=1,ζc=20,a=0).
3.1.1.6. ábra. Az x3 momentum időbeli változása (F=2,g=1,ζc=20,a=0).
3.1.1.7. ábra. Az átlagos kristályméret időbeli változása (F=2,Da=500,ζc=20,a=0).
3.1.1.8. ábra. Az anyalúg koncentráció időbeli változása (F=2,Da=500,ζc=20,a=0).
3.1.1.9. ábra. Az x3 momentum időbeli változása (F=2,Da=500,ζc=20,a=0).
(a) (b) (c)
3.1.1.10. ábra. A méreteloszlás időbeli változása (F=2,Da=500,g=1,ζc=20,a=0,∆x0=2e-5,B0,s=1e12,G0,s=1e-4).
paraméter meghatározható a Da illetve az a paraméter értékéből a 2.1.3. fejezetben megadott definíciók alapján:
3 3 4
6 v g
g c t a
p k k
s s
k = D (3.1.2)
g g c t
k s
= as
α (3.1.3)
Az összefüggésekben szereplő sc tényező értékét a (3.1.1) egyenletnek megfelelően határoztam meg, míg a kv formatényezőt, valamint a kg kristálynövekedési sebességi állandót szabadon választható, amellyel egyben a (2.1.3.23) kristálynövekedési kinetikai egyenletet is specifikáljuk.
A kollokációs módszerrel történő szimulációs vizsgálatokban a momentum modell analízisnél alkalmazott ugyanazon F, g és Da paraméterek mellett vizsgáljuk a kialakuló méreteloszlás dinamikáját. A (2.1.3.9-2.1.3.23) egyenletrendszer további paramétereit, amelyek a szimulációs vizsgálatok mindegyikében azonosak, valamint a kp sebességi állandó és α paraméter meghatározásához szükséges paramétereket a 3.2. táblázat tartalmazza. A kinetikai egyenletekben szereplő B0,s illetve G0,s paraméterek értékét a konkrét eredményeknél ismertetem.
3.2. táblázat. A (2.1.3.9-2.1.3.23) modell paramétereinek megadása.
a 0 yin 1.0 st 1/τ kg 1.0×10-5 (m/h)⋅(m3/kg)-g ρc 1760 kg/m3 ζc 20.0 τ 2.5 h kv 0.47 cin 211.475 kg/m3 y(0) 0.5
A számítási intervallum nagyságának, valamint az osztópontok helyének kezdeti értékeit 50 egyenlő hosszúságú végeselem megadásával definiáltam az alábbi egyenlettel:
, i=0,...,50 (3.1.4) x0
i xi = ∆
ahol ∆x0 a végeselemek kezdeti hosszát jelöli. A végeselemeken minden esetben azonos fokú közelítő polinomokat alkalmaztam, azaz Mi=M, i=0,...,N-1. A közelítő függvény , k=1,...,M, i=0,...,N-1 együtthatóinak kezdeti értékeit, amelyek a kezdeti méreteloszlás értékeit definiálják a kollokációs pontokban, az alábbi egyenlettel adtam meg:
i
fk
, k=1,...,M; i=0,...,N-1 (3.1.5)
( )
0 =0i
fk
azaz a kezdeti számítási intervallumon zérus kezdeti eloszlást definiáltam. Az oldott komponens és oldószer tömegmérlegben szereplő x2 és x3 momentumokat, amelyek megegyeznek a (2.1.2.11-2.1.2.17) momentum modell megfelelő momentumaival, a kollokációs módszer eredményeként kapott közelítő megoldásból határoztam meg két lépésben. Az első lépésben a ς dimenziómentes méretváltozó mentén rendelkezésre álló közelítő megoldásból Gauss-Lobatto kvadratúrával meghatároztam a méreteloszlás momentumait a 2.1.3. fejezetben ismertetett módon. Minden esetben 19 kvadratúra pontot használtam az egyes végeselemeken. A második lépésben a (2.1.3.18) egyenlet felhasználásával a ς méretváltozó mentén megkapott momentumokat az sn, sL és s , i=0,1,2,3 tényezők segítségével az , i=0,1,2,3 momentumokká transzformáltam.
i
xi
A számítási intervallum adaptív végeselem felosztásának paramétereit a 3.3. táblázat tartalmazza.
3.3. táblázat. A számítási intervallum adaptív felosztásának paraméterei.
M 7 r2 0.5 ND 100 I1 1 r1 1.05 NS 20 NI 1 NQ 19
Az ortogonális kollokációval kapott méreteloszlás dinamikai változását a 3.1.10. és M.1.1.-M.1.4. ábrák (1. Melléklet) mutatják különböző paraméter értékek esetében. A 3.1.10. és M.1.1.-M.1.2. ábrákon az F paraméter különböző értékei mellett kapott eredmények láthatók. Az ábrák (a), (b) és (c) esetei különböző időpillanatig mutatják a méreteloszlás kifejlődését, amelyeken megfigyelhető a kristályosítási folyamat korai szakaszában kialakuló méreteloszlás dinamikája is. A méreteloszlásban bekövetkező változások jól értelmezhetők az anyalúg koncentráció változását mutató 3.1.1.2. ábra alapján. F=2 esetben a méreteloszlás kifejlődése egy meredek front kialakulásával kezdődik, amely folyamatosan halad előre a mérettengely mentén. Az anyalúg koncentráció és így a gócképződési sebesség emelkedésével a méreteloszlás peremértéke folyamatosan növekszik és a meredek front fokozatosan ellaposodik (3.1.1.10a. ábra). Az anyalúg koncentrációban ξ=0.6 körül egy maximum figyelhető meg, majd a koncentráció csökkenésével a gócképződés mértéke is csökkenni kezd. A gócképződés csökkenésével csökken a peremérték és egy maximummal rendelkező méreteloszlás alakul ki. ξ=1.0 után az eloszlás peremértéke közelítőleg nullára csökken és az eloszlás maximuma a kristályok növekedésével párhuzamosan halad előre, mialatt a kristályosítóból távozó kristályoknak köszönhetően az eloszlás ellaposodik (3.1.1.10b.
ábra). A kristályok fokozatos eltávozásával az oldott komponens kiválásához rendelkezésre álló kristálytömeg csökkenni kezd, amely megfigyelhető a 3.1.1.3. ábrán is. Ennek következtében az anyalúg koncentráció ξ=2 után ismét emelkedni kezd, és ξ=5 után a peremérték ismételt emelkedése figyelhető meg. A peremérték fokozatosan növekszik az anyalúg koncentrációban megfigyelhető újabb maximum kialakulásáig és a koncentráció ismételt csökkenésével az előzőnél egy kisebb maximummal rendelkező eloszlás alakul ki. A folyamatok ismétlődésével végül egy csillapítatlan oszcilláció kialakulása figyelhető meg a méreteloszlás dinamikájában (3.1.1.10c. ábra).
F=0.1 esetben a méreteloszlás kifejlődése szintén egy meredek front kialakulásával kezdődik, de az eloszlás peremértéke nem változik jelentősen (M.1.1a. ábra). A 3.1.1.2.
ábrán látható, hogy az F=0.1 paraméter esetében az anyalúg koncentráció röviddel a kristályosító indulása után meredeken csökkenni kezd. Ennek következménye jól látható az M.1.1b ábrán, ξ=0.3 után az anyalúg koncentráció csökkenésének hatására a gócképződés sebessége és ennek megfelelően a peremérték is hirtelen lecsökken. Az így kialakuló eloszlásprofilban egy meredeken emelkedő, egy középső egyenletesen csökkenő, valamint egy meredeken csökkenő szakasz figyelhető meg, amely fokozatosan halad előre a mérettengely mentén. A koncentráció és gócképződés további csökkenésével, valamint a kristályok fokozatos kiürülésével a méreteloszlás fokozatosan ellaposodik és a peremérték közelítőleg nulla lesz. ξ=5 után a koncentráció és a gócképződés felfutásával mind az anyalúg koncentráció, mind a méreteloszlás csillapodó oszcillációt mutat, és végül beáll az exponenciális eloszlásnak megfelelő stacionárius állapot (M.1.1c. ábra). A 3.1.1.10c. és M.1.1c. ábrák összehasonlításából megállapítható, hogy a méreteloszlás domináns tartománya jelentősen lecsökken az F=0.1 paraméter mellett.
Az F=0.01 érték esetében a méreteloszlás kezdeti dinamikája hasonló, mint az előző esetben (M.1.2a. ábra), azonban a méreteloszlás peremértéke egy maximum eléréséig folyamatosan növekszik (M.1.2b. ábra). ξ=0.5 után a koncentráció és a gócképződés csökkenésének köszönhetően a peremérték csökkenni kezd, majd egy minimumon
keresztül, a koncentráció növekedésének megindulásával eléri a stacionárius állapotot.
A stacionárius állapotban egy exponenciális eloszlás alakul ki (M.1.2c. ábra).
Az M.1.3a-c. ábrákon a méreteloszlás dinamikai változása látható a Da=500, Da=100 és Da=10 paraméterek mellett. Mint az ábrákból látható, a Da paraméter csökkenésével a kialakuló eloszlásprofil jellege nem változik, mint az a különböző F paraméterek melletti vizsgálatoknál megfigyelhető volt. Az M.1.3a-c. ábrák összehasonlításából látható, hogy a Da paraméter csökkenésével a méreteloszlásban jelentkező csillapítatlan oszcillációt csillapodó oszcilláció váltja fel és a kristályosító működése stabilizálódik.
Da=10 esetben jól látható a méreteloszlásban jelentkező csillapodó oszcilláció, amely egy csillapodóan oszcilláló exponenciális eloszlás. Az ábrákból továbbá megfigyelhető, hogy a Da paraméter csökkenésével a kristályosítóban kialakuló méreteloszlás peremértékének nagyságrendje csökken, amely a csökkenő gócképződés hatását mutatja a Da paraméter csökkenésével.
Az M.1.4a-c. ábrák a g paraméter hatását mutatják a méreteloszlás dinamikájára a g=0.5, g=1.0 és g=1.5 esetekben. Az M.1.4a-c. ábrákból látható, hogy a g paraméter növekedésével a méreteloszlásban megfigyelhető oszcilláció frekvenciája növekszik. A g=1.5 esetben jól látható a méreteloszlásban kialakuló csillapodó oszcilláció, amely egy csillapodóan oszcilláló exponenciális eloszlás. A g=0.5 és g=1.5 esetekben megfigyelhető, hogy a kristályok mérettartománya egy nagyságrenddel kisebb illetve nagyobb, mint a g=1 paraméter esetében, amely a g paraméter által meghatározott kristálynövekedési sebesség hatását mutatja.
A 3.1.1.11.–3.1.1.13. ábrák különböző időpillanatokban mutatják a méreteloszlás profilját, az osztópontok aktuális elhelyezkedését, valamint az eloszlásértékek osztópontbeli értékeit az eloszlásgörbén az F paraméter különböző értékei mellett. Az ábrák jól szemléltetik, hogy az eloszlásgörbe meredekebb szakaszainál az osztópontok sűrűbben helyezkednek el, míg a kisebb meredekséggel rendelkező szakaszoknál az osztópontok távolsága jelentősen nagyobb. A 3.3. táblázat adataiból látható, hogy a számítási intervallum adaptív végeselem felosztásánál az NI=1 paramétert használtam, azaz a méreteloszlás görbéjét az egész hosszában azonos hosszúságú ívszakaszokkal osztottam fel, és mint a 3.1.1.11.–3.1.1.13. ábrákból látható, az eloszlásértékek osztópontbeli értékei az eloszlásgörbén egyenletes távolságokban helyezkednek el.
Összefoglalóan megállapítható, hogy a kidolgozott számítási intervallum felosztási módszer hatékonynak bizonyult az osztópontok megfelelő elhelyezésében a legkülönbözőbb formájú eloszlásprofilok esetében.
A 3.1.1.14.-3.1.1.16. ábrák a részintervallumokat elválasztó osztópontok helyének időbeli változását mutatják az idő függvényében. Mint az előzőkben látható volt, a különböző F paraméterek esetében jelentősen különbözik egymástól a méreteloszlás kialakulása és dinamikája. Az ábrákon jól látható a méreteloszlás meredek szakaszainak mozgása, ahol az osztópontok sűrűbben helyezkednek el. A 3.1.1.14b. ábrán megfigyelhető a méreteloszlásban F=2 esetben kialakuló oszcilláció. Az osztópontok többségének elhelyezkedése a ς=[0,20] intervallumban változik, amely a kialakuló méreteloszlás domináns részét foglalja magában. A legfelső görbe a ςN osztópont mozgását adja meg, amelynek értéke közelítőleg a maximális kristálymérettel egyezik meg. A felsőbb görbék által megadott osztópontoknál az eloszlás értékének nagyságrendje kisebb, mint a méreteloszlás domináns tartományában, de a méreteloszlás harmadik momentumának pontos meghatározásában fontos szerepet töltenek be. A 3.1.1.15a. ábrán jól látható az M.1.1b. ábrán bemutatott méreteloszlás dinamikának megfelelő osztópont mozgás. Az ábrán megfigyelhető az eloszlás meredeken emelkedő és lecsengő szakaszának mozgása, ahol az osztópontok sűrűbben helyezkednek el. A 3.1.1.14b. és 3.1.1.15b. ábrák összehasonlításából látható, hogy az F
(a) (b)
3.1.1.11. ábra. A méreteloszlás profilja és az osztópontok elhelyezkedése (F=2,Da=500,g=1,ζc=20,a=0).
(a) (b)
3.1.1.12. ábra. A méreteloszlás profilja és az osztópontok elhelyezkedése (F=0.1,Da=500,g=1,ζc=20,a=0).
(a) (b)
3.1.1.13. ábra. A méreteloszlás profilja és az osztópontok elhelyezkedése (F=0.01,Da=500,g=1,ζc=20,a=0).
(a) (b) 3.1.1.14. ábra. Az osztópontok helyének időbeli változása
(F=2,Da=500,g=1,ζc=20,a=0).
(a) (b)
3.1.1.15. ábra. Az osztópontok helyének időbeli változása (F=0.1,Da=500,g=1,ζc=20,a=0).
(a) (b)
3.1.1.16. ábra. Az osztópontok helyének időbeli változása (F=0.01,Da=500,g=1,ζc=20,a=0).
paraméter értékének csökkenésével a méreteloszlásban megfigyelhető csillapítatlan oszcilláció megszűnik és csillapodó oszcilláció alakul ki. A 3.1.1.16b. ábrán jól látható az exponenciális eloszlás kialakulása, az osztópontok elhelyezkedése fokozatosasan sűrűsödik a ς=0 méret felé.
A 3.1.1.17. ábrán a momentum modell megoldásával kapott x3 momentum dinamiká-
3.1.1.17. ábra. Az x3 momentum időbeli változása (Da=500,g=1,ζc=20,a=0).
ját hasonlítom össze az adaptív végeselemű ortogonális kollokációval kapott méreteloszlásból meghatározott x3 momentummal. Az összehasonlítást az F paraméter különböző értékei mellett végzett dinamikai vizsgála-toknál végezzük el, mivel ezekben az esetekben kaptuk a legkülönbözőbb lefutású méreteloszlás dinamikát. Az összehasonlítás céljára az x3 mo-mentumokat választottuk, mivel a (2.1.2.11-2.1.2.17) momentum modell az x0-3 momentumokra nézve már zárt, és a (2.1.3.9-2.1.3.23) modellben előforduló legnagyobb rendű momen-tum az x3 momentum. A 3.1.1.17 ábrából látható, hogy az F paraméterek mindegyikénél kitűnően egyezik a két matematikai módszer által kapott x3 momentum dinamikája, amely a kidolgozott adaptív végeselem felosztási módszer hatékonyságát mutatja a kristályosítók populációs mérlegegyenleteinek megoldásában.