VII.3. A TREX-Euler modell

Az Euler-típusú modellek a légkör meghatározott részét rácshálóval fedik le, és ennek pontjaira oldják meg a fizikai folyamatokat leíró egyenletrendszereket úgy, hogy valamilyen állandó, vagy változó időlépésenként kapják meg a megoldást. A TREX-Euler modell (Mészáros et al., 2010) egy Közép-Európát lefedő rácson számolja a különböző szennyezőanyagok diszperzióját (VII.5. ábra).

A modell a terjedés leírásához használt légköri transzportegyenletekben az advekciót, a függőleges és vízszintes diffúziót, az ülepedést, a kémiai reakciókat és az emissziót veszi figyelembe (lásd: I.41. egyenletet is):

,

(VII.1.)

ahol az adott anyagfajta átlagos koncentrációja, az átlagos háromdimenziós szélmező, kd a száraz és kw a nedves ülepedési együttható, a turbulens diffúziós együtthatók vektora, melynek az egyes komponensei a horizontális és a vertikális diffúziós együttható, E az adott anyagfajta emissziós értéke, R a kémiai folyamatok hatására bekövetkező koncentrációváltozás sebessége.

A modell kvázi-háromdimenziós, mint a mai gyakorlatban leginkább használt modellek többsége. A modellben a légkör vizsgált részét függőleges irányban rétegekre bontjuk, a rétegekben a koncentráció-változást külön-külön kétdimenziós modellek írják le, a rétegek közötti anyagszállítást (transzportot) pedig az erre alkalmas fizikai modellekkel szimuláljuk. A vertikális átkeveredés pontos leírása érdekében 32 magassági szintet használunk. A felszíntől 200 méteres magasságig 20, a 200 métertől a 3000 méteres magasságig pedig további 12 szintet definiáltunk, olyan módon, hogy 200 méter alatt és felett is külön-külön egyenlő legyen az egyes szintek közötti nyomáskülönbség. (VII.6. ábra).

A TREX-Euler terjedési modell felépítése

VII.5. ábra. A TREX-Euler terjedési modell felépítése

A TREX-Euler a modell vertikális szintjei

VII.6. ábra. A TREX-Euler a modell vertikális szintjei

Az időlépcső és a rácsfelbontás megválasztása a megoldás pontossága szempontjából döntő fontosságú. Az időlépcső nem megfelelő megadása numerikus hibát, konvergencia és stabilitási problémákat eredményezhet.

A diffúzió számításánál stabil megoldást kapunk, ha a Kx horizontális (x irányú) turbulens diffúziós állandó, a Δt időlépés és a Δx rácsfelbontás között az alábbi összefüggés áll fenn:

.

(VII.2)

Advekció számolásakor stabil megoldást kapunk, ha teljesül az ún. CFL-feltétel (Courant–Friedrichs–Levy-feltétel). Eszerint a parciális differenciálegyenletek (így a transzport egyenlet) megoldásában a Δx alkalmazott rácsfelbontás illetve a modellben használt Δt időlépcső hányadosa nagyobb kell, hogy legyen, mint az adott irányba leggyorsabban terjedő hullám fázissebessége (szemi implicit séma alkalmazásával, lásd a II.3. fejezetet is). Az x irányba terjedő hullám fázissebessége legyen uc. Ekkor:

.

(VII.3)

Látható, hogy adott diffúziós állandó és a fázissebesség mellett a rácsfelbontás növelésével, illetve az időlépés csökkentésével lehet a megoldás stabilitását biztosítani. Azonban ha durva rácsfelosztást használunk, akkor a kibocsátás azonnal nagy területre átlagolódik, ami ’szétkeni’ a meredek gradienst, és nagy numerikus diffúziót okoz. Ennek következtében a csóvában alábecsüljük a maximális koncentrációt és túlbecsüljük a csóva szélességét. Az időlépcső csökkentésével – és kis rácsfelbontás esetén – a számítási idő nő meg jelentősen. Ezek együttes figyelembevételével kell valamilyen kompromisszumos idő- és rácsfelbontást választani. Az általunk kifejlesztett modell egy Közép-Európát lefedő területen, 0,0375 × 0,025 fokos (~2,5 km × ~2,5 km) térbeli felbontással, 10 másodperces időlépéssel számítja a forrásokból kiinduló szennyezőanyag koncentrációját és ülepedését.

VII.3.1. A program felépítése

A programkód több részből tevődik össze. A főprogram az adatok beolvasását, a különböző függvények meghívását és ciklusba szervezését, végül az eredmények kiíratását végzi.

Az első almodul a horizontális és vertikális határfeltételeket adja meg. A tartomány peremén ’no-flux’ határfeltételt használtunk, vagyis azt feltételezzük, hogy a határon nincs anyagáramlás. Külön rutin végzi az advekció, a vertikális és a horizontális diffúzió számítását, illetve a magassági szintek meghatározását. A számításához szükséges Monin–Obukhov-hosszat (L) és a vertikális turbulens diffúziós együtthatót (Kz) is külön függvény számítja. A különböző anyagtranszportok (advekció, diffúzió), illetve a kémiai reakció és ülepedés különálló számítására a későbbiekben leírt operátor-splitting (operátor-szeletelés) módszer ad lehetőséget.

A horizontális diffúziós együtthatót állandónak vettük a modellben. A vertikális turbulens diffúzió a K-elmélet alapján kerül számításra és a magasságfüggő Kz diffúziós együtthatóval vesszük figyelembe. A modell futási idejének lerövidítése érdekében a Kz kiszámítását sztochasztikus, véletlen módszerrel végezzük (Mészáros et al., 2010). Az egyes anyagfajták vertikális eloszlását (profilját) a turbulens diffúziós egyenlettel adjuk meg:

.

(VII.4.)

A vertikális turbulens diffúziós együtthatót a Monin–Obukhov-féle hasonlósági elmélet felhasználásával (lád az I.2.4. fejezetet is) parametrizáltuk a következő módon:

.

(VII.5.)

Eszerint a turbulens diffúziós együttható adott z szinten a keveredési réteg magassága (Hz) a súrlódási sebesség (u*), a stabilitási függvény (Ψ), a Kármán-állandó (κ) és a Monin–Obukhov-féle hossz (L) függvényeként írható fel (Brandt, 1998). A Monin–Obukhov-féle hosszat és a súrlódási sebességet Ács (2004) leírása alapján határoztuk meg (lásd az I.2.4.1. fejezetet is). A száraz ülepedés számításakor állandó ülepedési együtthatót vettünk figyelembe. Nedves ülepedést akkor számítottunk, amikor a relatív nedvesség 80% feletti (Brandt, 1998). Továbbá feltételeztük, hogy ülepedés csak az első, földközeli rétegből történhet.

 A program az adatok beolvasása, a magassági szintek, továbbá a kezdeti- és peremfeltételek megadása után minden időlépcsőben szintenként kiszámítja az advekciót, ezután légoszloponként meghatározza a vertikális keveredést, a turbulens diffúziós együtthatót és a hozzá szükséges Monin–Obukhov-féle hosszat. Ezután következik – radioaktív szennyezők esetén – a radioaktív bomlás számítása. Végül a földközeli (vagy más néven a felszínközeli) rétegben az ülepedés meghatározása történik. A következő időlépésben újból elölről ismétlődik a fent leírt folyamat.

VII.3.2. Numerikus megoldás

Az elfogadható pontosságú 3D modellek hatalmas számítási kapacitást és kifinomult numerikus megoldási technikát igényelnek. A TREX-Euler modellben az egyenletek megoldására az operátor-splitting módszert használtuk, vagyis a parciális differenciálegyenletekben szereplő tagokat külön-külön oldottuk meg. A térbeli transzport tagokat véges differencia sémával diszkretizáltuk. Az első lépésben csak az advekciós tagot (advekció hatását) vettük figyelembe és így határoztuk meg a cadv koncentrációt (az advekció hatására kialakuló új koncentráció eloszlás) a régi cold koncentráció-értékből:

.

(VII.6.)

Ezután az előzőekben kapott cadv koncentráció felhasználásával határoztuk meg a diffúzió hatására kialakuló cdiff koncentrációt (külön a vízszintes és függőleges irányú diffúzió számításával):

,

(VII.7.)

Végül a kémiai reakció és a száraz és a nedves ülepedés hatását az előző lépésben meghatározott koncentrációból számítottuk ki az alábbi egyenlet alapján:

.

(VII.8.)

Így a harmadik lépés után kapott cnew koncentráció az adott időlépés után mindhárom tényező hatását tartalmazza. Az egyenletekben Aadv az advekciós operátor Adiff a diffúziós operátor, míg Achem a kémiai reakciót és az ülepedést leíró operátor. Ezek megoldására különböző módszereket alkalmaztunk.

A parciális differenciálegyenletek egyik hatékony megoldási módszere a „method of lines” technika. A módszer lényege a transzport tagok térbeli diszketizációja után keletkezett közönséges differenciálegyenlet-rendszer időbeli integrálása megfelelő kezdeti- és peremfeltételek alkalmazásával. Az advekció térbeli diszkretizálására ’second upwind’ módszert, a vertikális diffúzió számítására ’first upwind’ módszert használtunk. Az első és másodrendű upwind módszerek az advekció és diffúzió megoldásának stabilitását biztosító sémák. A kémiai reakció, a száraz és a nedves ülepedés esetén nem szerepel térbeli derivált, ott csak az időbeli integrálást kell elvégezni. A diszkretizált tagok időbeli integrálására explicit Euler-sémát használtunk.

VII.3.3. Alkalmazások

A háromdimenziós TREX-Euler terjedési modellel első lépésben pontforrásokból a légkörbe jutó szennyezők terjedését vizsgáltuk. A VII.7. ábra a Paksi Atomerőmű területéről egy feltételezett baleset során a légkörbe kerülő jódizotóp felszínközeli eloszlását mutatja 6 órával egy feltételezett baleset után.

Példa a TREX-Euler terjedésszámító modell alkalmazására

VII.7. ábra.A TREX-Euler terjedésszámító modell alkalmazása: 131I koncentráció eloszlása egy Közép-Európát lefedő rácson, 6 órával egy feltételezett baleset után, folyamatos kibocsátás mellett, az ALADIN modell meteorológiai mezőinek felhasználásával.

A bemenő meteorológiai adatokat az Országos Meteorológiai Szolgálatnál operatívan futtatott ALADIN előrejelzési modell 0–48 órás előrejelzései szolgáltatták. Különböző időjárási helyzetekre végzett esettanulmányok (Mészáros et al., 2010) és részletes érzékenységi vizsgálatok mellett egy hosszabb időszakra (1 év) folyamatosan végeztünk szimulációkat (az év minden egyes napján egy baleseti kibocsátást feltételezve és a diszperziót szimulálva), ami már statisztikai elemzéseket is lehetővé tett (Mészáros et al., 2012a).

Az Euler-modellt alkalmaztuk egyéb gázok terjedésének szimulációjára is, folyamatos kibocsátás mellett (Mészáros et al., 2012b). Az ehhez szükséges emissziós adatokat az EMEP adatbázisból származtattuk. A VII.8. ábrán a szárazon, 48 óra alatt (2006. 03. 24–25.) kiülepedett kén-dioxid mennyiségét mutatjuk be példaként.

TREX-Euler modellel számított kén-dioxid ülepedés

VII.8. ábra. A 48 óra alatt szárazon kiülepedett kén-dioxid. (Az értékek mértékegysége cm2.)