© Strojni{ki vestnik 48(2002)8,449-458 © Journal of Mechanical Engineering 48(2002)8,449-458 ISSN 0039-2480 ISSN 0039-2480 UDK 004.94:532.5:536.7 UDC 004.94:532.5:536.7 Izvirni znanstveni ~lanek (1.01) Original scientific paper (1.01) Simuliranje izotermnega QUEOS preskusa me{alne faze eksplozije pare Q08 Simulation of the Isothermal QUEOS Steam-Explosion Premixing Experiment Q08 Matja` Leskovar - Borut Mavko Mešalna faza eksplozije pare obsega interakcijo curka taline z vodo pred izbruhom eksplozije pare. Da bi bolje spoznali hidrodinamične procese med mešalno fazo eksplozije pare, izvajajo poleg “vročih” preskusov, pri katerih je pomembno uparjanje vode, tudi “hladne” izotermne preskuse. Značilnost izotermnih preskusov mešalne faze eksplozije pare je, da se od vseh treh faz, to je vode, zraka in kroglic, z drugima dvema fazama pomešajo le kroglice, medtem ko ostaneta voda in zrak ločena z gladko medfazno ploskvijo. To dejstvo smo upoštevali pri razvoju izvirnega kombiniranega večfaznega modela mešalne faze eksplozije pare. Pri razvitem kombiniranem večfaznem modelu obravnavamo vodo in zrak z modelom proste površine kot eno, združeno fazo z nezveznimi faznimi lastnostmi na medfazni ploskvi voda -zrak, kroglice pa obravnavamo kot običajno z modelom večfaznega toka, pri katerem pomenijo kroglice razpršeno fazo, združena faza voda - zrak pa kontinuum, ki obdaja kroglice. Razviti kombinirani večfazni model smo preverili na izotermnem preskusu Q08, ki so ga izvedli na napravi QUEOS. Simuliranja smo opravili z različnimi numeričnimi obravnavami medfazne ploskve voda -zrak (metoda z nivojsko funkcijo, metoda visoke ločljivosti, metoda “nazaj”) za nestisljiv in stisljiv primer. Rezultate simuliranj smo primerjali med seboj in s podatki preskusa. © 2002 Strojniški vestnik. Vse pravice pridržane. (Ključne besede: tok večfazni, modeli večfazni, eksplozija parne, metode z nivojsko funkcijo, simuliranje) The premixing phase of a steam explosion covers the interaction of the melt jet with the water prior to any steam explosion occurring. To get a better insight into the hydrodynamic processes during the premixing phase, in addition to “hot” premixing experiments, where the water evaporation is significant, “cold” isothermal premixing experiments were also performed. The special feature of isothermal premixing experiments is that three phases are involved - the water, the air and the spheres ’ phase - but only the spheres’ phase mixes with the other two phases, whereas the water and air phases do not mix and remain separated by a free surface. Our idea was to treat the isothermal premixing process with an original, combined multiphase model. In the developed combined multiphase model the water and air phases are treated with a free-surface model as a single, joint phase with discontinuous phase properties at the water-air interface. The spheres are treated, as is usual, with a multiphase flow model, where the spheres represent the dispersed phase and the joint water-air phase represents the continuous phase. The developed combined multiphase model was validated against the QUEOS isothermal premixing experiment Q08. A numerical analysis using different treatments of the water-air interface (level set, high-resolution, upwind) was performed for the incompressible and compressible cases and the results were compared with experimental measurements. © 2002 Journal of Mechanical Engineering. All rights reserved. (Keywords: multiphase flow, multiphase models, steam explosion, level set methods, simulations) 0 UVOD Med resno nezgodo v jedrski elektrarni, ko pride staljena sredica po različnih poteh v stik s hladilno vodo, lahko pride do eksplozije pare. Eksplozijo pare razdelimo na več faz. Prva, mešalna faza eksplozije pare je najbolj pomembna, ker določa začetne pogoje morebitne eksplozije pare in podaja 0 INTRODUCTION During a severe reactor accident following core meltdown when the molten fuel comes into contact with the coolant water a steam explosion may occur. A steam explosion can be divided into more stages. The first stage, called premixing, is the most important since it provides the initial conditions for a possible steam explosion and gfin^OtJJIMISCSD 02-8 stran 449 |^BSSITIMIGC Leskovar M. - Mavko B.: Simuliranje izotermi~nega - Simulation of an Isothermal največjo količino staljenega jedrskega goriva, ki bi lahko sodelovalo pri eksploziji. Da bi bolje spoznali hidrodinamične procese med prodiranjem taline v vodo, izvajajo poleg “vročih” preskusov, pri katerih je pomembno predvsem uparjanje vode, tudi “hladne” izotermne preskuse. Pri teh izotermnih preskusih mešalne faze eksplozije pare spuščajo različne curke hladnih kroglic v posodo, napolnjeno z vodo in opazujejo dogajanja [1]. Značilnost izotermnih preskusov mešalne faze eksplozije pare je, da se od prisotnih treh faz, to je vode, zraka in kroglic, s preostalima dvema fazama pomešajo le kroglice, medtem ko ostaneta voda in zrak ločena z gladko medfazno ploskvijo. To dejstvo smo upoštevali pri razvoju izvirnega kombiniranega večfaznega modela mešalne faze eksplozije pare. Pri razvitem kombiniranem večfaznem modelu obravnavamo vodo in zrak z modelom proste površine kot eno, združeno fazo z nezveznimi faznimi lastnostmi na medfazni ploskvi voda - zrak, kroglice pa obravnavamo kot običajno z modelom večfaznega toka, pri katerem pomenijo kroglice razpršeno fazo, združena faza voda - zrak pa kontinuum, ki obdaja kroglice [2]. 1 GLAVNE ENAČBE determines the maximum quantity of melt, which might then be involved in the explosion. To investigate the mixing process associated with the melt penetration, in addition to “hot” premixing experiments, where the water evaporation is significant, “cold” isothermal premixing experiments were also performed to get a better insight into the hydrodynamic processes during the premixing phase. In these isothermal premixing experiments, different jets of spheres are injected in a water pool [1]. The special feature of isothermal premixing experiments is that three phases are involved – the water, the air and the spheres’ phase – but only the spheres’ phase mixes with the other two phases, whereas the water and air phases do not mix and remain separated by a free surface. Our idea was to treat the isothermal premixing process with an original combined multiphase model. In the developed combined model the water and air phases are treated with a free-surface model as a single, joint phase with discontinuous phase properties at the water-air interface. The spheres are treated, as is usual, with a multiphase flow model, where the spheres represent the dispersed phase and the joint water-air phase represents the continuous phase [2]. 1 GOVERNING EQUATIONS Fazna verjetnost a in gostota rc združene faze voda - zrak c sta definirani kot: The phase presence probability a and the density rc of the joint water-air phase c are defined as: ac=aw+aa, rc= aw rw+aa ra (1), kjer indeksa w in a označujeta fazi vode in zraka. S to definicijo lahko večfazen tok kroglic, vode in zraka obravnavamo kot navidez dvofazen tok razpršenih kroglic d v kontinuumu združene faze voda - zrak c. Vsako fazo p v tako definiranem navidez dvofaznem toku smo opisali s kontinuitetno: where the indices w and a correspond to the water and air phases. With this definition the multiphase flow of spheres, water and air can be described as a quasi two-phase flow of dispersed spheres d in the continuous joint water-air phase c. Each phase p in the so-defined quasi two-phase flow is described using the continuity: Lt (aprp)+V-(aprpvp)=0 in gibalno enačbo: and the momentum equation: dt (2) (3), ki smo ju izpeljali s statističnim povprečenjem enofaznih enačb [3]. Pri tem smo medfazno sklopitveno silo M definirali kot: obtained by ensemble averaging the single-phase r equations [3]. The interfacial friction force M is defined as: Md=M drag + Mvd m+Mldift, acMc=-ad M kjer silo medfaznega trenja računamo z [4]: Mdd r ag = 0, 44 3 r 4 d d silo navidezne mase z [4]: where the drag force is calculated from [4]: v c - v d l(v c - v d ) the virtual mass force from [4]: (4), (5), VBgfFMK stran 450 Leskovar M. - Mavko B.: Simuliranje izotermi~nega - Simulation of an Isothermal M vm =0.5rc|—- — in dvižno silo z [4]: c| Dt Dt and the lift force from [4]: Mldift =0.5rc ( vc-v )x(Vxvc (6) (7). Medfazno trenje kroglic na medfazni ploskvi The spheres’ drag at the water-air surface was voda - zrak smo določevali na poseben način, ki upošteva treated in a special way by considering the nezveznost prehoda zrak - voda, kakor je opisano v [5]. discontinuous air-water transition as described in [5]. Tlačno enačbo za predpostavljeno skupno The pressure equation for the assumed tlačno polje p smo izpeljali iz vsote kontinuitetne common pressure field p was derived from the sum of enačbe za stisljivo tekočino (2) po vseh fazah s the compressible flow continuity equation (2) over projekcijsko metodo [6]: all phases using the projection method [6]: kjer je: I p rp( pn+1-pn ) +v-( apn+1 r np+1/ 2 ) Dt-v- rp dpn an+1 p Vpn+1 vr n+1/ 2 = vr n + Dt where: ( vnp-v) vnp+g+ rp Mnp rp n+1/ 2 Dt J J r np+1/ 2 =r(apn+1 ,p n ) (8), (9). 2 SIMULIRANJA PRESKUSA MEŠALNE FAZE Za analizo izvirnega kombiniranega večfaznega modela smo izbrali preskus mešalne faze Q08, ki so ga izvedli na napravi QUEOS [7]. Pri tem preskusu so spuščali molibdenove kroglice premera 4,2 mm in skupne mase 10 kg z višine 1,3 m skozi cev v posodo v obliki kvadra z osnovno ploskvijo 70 cm x 70 cm in višino 138 cm, napolnjeno z vodo. Premer curka kroglic je bil 18 cm, hitrost kroglic na vodni gladini 5,12 m/s in fazna verjetnost kroglic 0,17. Višina vodne gladine v posodi je bila 100 cm, začetni tlak plina nad vodno gladino pa 0,975 bar. Pri simuliranjih preskusa smo obravnavali le z vodo napolnjeni spodnji del posode in 20 cm debelo plast plina nad vodno gladino (sl. 1). Vse druge značilnosti preskusa smo upoštevali z ustreznimi robnimi pogoji. Predpostavili smo, da je plin nad vodno gladino zrak. V resnici je plin v posodi neznana mešanica vodne pare in neukapljivega argona. Simuliranja smo opravili v valjnem koordinatnem sistemu ob predpostavki osne simetričnosti problema. Ta predpostavka je le delno upravičljiva, saj porazdelitev kroglic predvsem na čelu curka ni osno simetrična zaradi sprostilnega mehanizma kroglic, pri katerem se simetrično odpreta dve ravni premični plošči [8]. Dejanski polmer valjaste posode modela (41 cm) smo določili tako, da sta se ujemali površini osnovnih ploskev modela in preskusa. Padajoči curek kroglic smo v model vključili z robnimi pogoji (RP). Robne vrednosti hitrostnega polja kroglic v in polja fazne verjetnost kroglic a na zgornjem robu simulirnega območja 20 cm nad vodno gladino smo določili posredno iz znanih 2 SIMULATIONS OF THE PREMIXING EXPERIMENT The QUEOS premixing experiment Q08 was chosen for the analysis of the original combined multiphase model [7]. In this experiment molybdenum spheres with a diameter of 4.2 mm and a total mass of 10 kg were discharged from a height of 1.3 m through a tube into a water-filled vessel with an inner cross-section of 70 cm x 70 cm and a height of 138 cm. The spheres’ jet diameter was 18 cm and the spheres entered the water at a velocity of 5.12 m/s and a phase-presence probability of 0.17. The water level in the vessel was 100 cm, and the initial gas pressure over the water surface was 0.975 bar. In the experimental simulations only the part of the vessel filled with water together with a 20-cm-wide gas zone over the water surface was modeled (Figure 1). All other experimental features were taken into account with appropriate boundary conditions. It was assumed that the gas over the water was air. In fact it was an unknown mixture of steam and non-condensable argon. The simulations were performed in the cylindrical coordinate system and assumed an axial symmetry for the problem. This was not really the case in the experiment since the release mechanism with two doors opening in opposite directions caused, especially at the front of the spheres’ cloud, an initially non- axially symmetric distribution of the falling spheres [8]. To match the cross-sectional area of the vessel, an equivalent radius of 41 cm was used in the axially symmetric model. The falling spheres were introduced in the model with boundary conditions (BCs). The appropriate boundary conditions for the spheres’ velocity vz and the spheres’ presence probability as at the upper edge of the simulation region 20 cm above the water surface were | IgfinHŽslbJlIMlIgiCšD I stran 451 glTMDDC Leskovar M. - Mavko B.: Simuliranje izotermi~nega - Simulation of an Isothermal 9 cm curek kroglic spheres jet zrak air vodna gladina water surface voda water simetrijska os symmetry axis robni pogoji / boundary conditions: vz = 4,72 m/s as = 0,184 D t = 0,0443 s p = 0,975 bar 20 cm 100 cm 41 cm Sl. 1. Shema modela QUEOS preskusa Q08 Fig. 1. Schematic model of QUEOS experiment Q08 povprečnih vrednosti na vodni gladini. Čas preleta curka kroglic Dt mimo zgornjega roba simulirnega območja smo določili iz razmerja znane celotne mase kroglic in izračunanega masnega toka kroglic [2]. Polja faznih verjetnosti in hitrostna polja faz smo računali z metodo velike ločljivosti [9], ki je drugega reda natančnosti. Metoda velike ločljivosti je temeljila na privetrni metodi prvega reda natančnosti in metodi Lax-Wendroff, ki je drugega reda natančnosti [10]. Tlačno enačbo smo reševali z metodo CGSTAB [6], ker hitro konvergira tudi za velika razmerja gostot posameznih faz [8]. Nezvezno gostoto združene faze voda - zrak smo določevali z metodo z nivojsko funkcijo [11], katero so razvili za reševanje problemov s prosto površino in se v zadnjih letih veliko uporablja. Simuliranja preskusov smo opravili na mreži velikosti 41 x 120 točk. Tako je bila mrežna razdalja točno 1 cm, kar je, kakor je pokazala konvergenčna analiza, primerno za takšne vrste simuliranj [5]. 3 REZULTATI SIMULIRANJ determined from the known average values at the water’s surface. The time interval Dt during which the spheres were falling through the upper boundary of the simulated region was established from the total spheres’ mass and the calculated spheres’ mass-flow rate [2]. The phases’ presence probabilities and the phases’ velocities were calculated using the second-order-accurate high-resolution method [9]. The highresolution method was based on the first-order-accurate upwind method and the second-order-accurate Lax-Wendroff method [10]. The pressure equation was solved using the CGSTAB method [6], since it also converges quickly for high phase-density ratia [8]. The discontinuous density of the joint water-air phase was determined with the front-capturing level-set method [11], which was developed for free-surface problems and has been widely used in recent years. The experimental simulations were performed on a mesh with 41 x 120 grid points. The grid spacing was exactly 1 cm, which is adequate for this type of simulation, as the convergence analysis showed [5]. 3 RESULTS OF THE SIMULATIONS Da bi bolje spoznali prednosti razvitega To get a better insight into the advantages of kombiniranega večfaznega modela, pri katerem the developed combined multiphase model with the VH^tTPsDDIK stran 452 Leskovar M. - Mavko B.: Simuliranje izotermi~nega - Simulation of an Isothermal obravnavamo medfazno ploskev voda - zrak z metodo z nivojsko funkcijo (NF - LS), ki so jo razvili za reševanje problemov s stično površino, smo medfazno ploskev voda - zrak primerjalno določevali še z metodo velike ločljivosti (VL - HR), ki je drugega reda natančnosti, in privetrno metodo (PM - UW) prvega reda natančnosti. Na sliki 2 so prikazani posnetki preskusa Q08, ki so jih napravili s hitro slikovno kamero ob osvetlitvi z zadnje strani. Čas posnetkov je merjen od sprostitve kroglic. Ker točni podatki niso bili na voljo, smo lahko čas posnetkov določili le približno. Na slikah je vidna referenčna mreža s celicami velikosti 10 cm x 10 cm, ki je pritrjena na notranjo stran okenskih šip posode. Ker so okna široka le 50 cm, je na vsaki strani in tudi na dnu posode zakrit 10 cm širok pas. Poudariti je treba, da pomeni večji del črno senčenega območja na posnetkih plinski steber in ne kroglic. incorporated level-set (LS) method, which was developed for free-surface problems, the water-air surface was, for a comparison, also determined with the second-order-accurate high-resolution (HR) method and the first-order-accurate upwind (UW) method. The images of experiment Q08 taken with a high-speed film camera with backlighting are presented in Figure 2. The time, which is measured after the release of the spheres, could only be determined approximately since no exact data was available. On the images the 10 cm x 10 cm reference grid mounted close to the inside of each of the vessel windows can be seen. Since the windows are only 50-cm wide, there is a 10-cm-wide zone on both sides and also on the bottom of the vessel, which is hidden. It should be stressed that most of the black shaded region on the images represents the gas chimney and not the spheres. -TTT T H----- J t =~520 ms t =~590 ms t =~650 ms t =~725 ms t =~775 ms Sl. 2. Posnetki preskusa Q08, napravljeni ob različnih časih (posnetki vzeti iz [7]) Fig. 2. Images of experiment Q08 taken at different times (images taken from [7]) Na slikah 3 do 5 so predstavljeni rezultati simuliranj, opravljenih z metodami NF, VL in PM. Slike prikazujejo gostoto združene faze voda - zrak (črni obrisi), fazno verjetnost kroglic (sivi obrisi) in hitrostno polje združene faze voda - zrak v valjnem koordinatnem sistemu (r, z) v ekvidistančnih 50 ms dolgih časovnih razmikih Čas je zaradi lažje primerjave merjen tako kakor pri preskusu po sprostitvi kroglic, to je 0,46 s pred začetkom simuliranja. Obrisi gostote združene faze voda - zrak (p = 997 kg/m3, p = 1,2 kg/m3) imajo vrednosti 2, 100 w 200, ..., 800, 900, 996 kg/m3, obrisi fazne verjetnosti kroglic pa imajo vrednosti 0,02; 0,05; 0,1; 0,15; itn. Hitrostno polje je predstavljeno z enotskimi vektorji, ki prikazujejo le smer vektorja hitrosti in ne tudi velikosti hitrosti, saj bi bile sicer slike nepregledne. Kakor smo pričakovali, ostane medfazna ploskev voda - zrak ostra med celotnim simuliranjem le pri izračunu z metodo NF. Pri metodi VL, ki je drugega reda natančnosti, predvsem pa pri metodi PM prvega The results of the simulations performed with the LS, HR and UW methods are presented in Figs. 3 to 5. The figures show the joint water-air phase density (black contours), the spheres’ presence probability (gray contours) and the joint water-air phase velocity field in the cylindrical coordinate system (r, z) in equi-distant 50-ms-long time intervals. The time is measured, as in the experiment, after the release of the spheres, i.e. 0.46 s before the start of the simulation. The joint water-air phase density (rw = 997 kg/m3, ra = 1.2 kg/m3) contours correspond to values 2, 100, 200, ... , 800, 900, 996 kg/m3, and the spheres’ presence probability contours correspond to the values 0.02, 0.05, 0.1, 0.15, etc. The velocity field is represented by unit vectors showing only the direction of the velocity vector but not the speed because the figures would be unclear. As expected, the water-air surface remains sharp during the whole simulation only during the LS-method calculation. For the second-order-accurate HR method, and in particular for the first-order- gfin^OtJJlMlSCSD 02-8 stran 453 |^BSSITIMIGC Leskovar M. - Mavko B.: Simuliranje izotermi~nega - Simulation of an Isothermal reda natančnosti, se združena faza voda - zrak nefizikalno razprši zaradi numerične difuzije. Numerična difuzija je splošna skupna pomanjkljivost metod končnih razlik. Rezultati simuliranj, opravljenih z metodo NF (slika 3), se kakovostno razmeroma dobro ujemajo s posnetki preskusa (sl. 2). t = 560 ms t = 610 ms accurate UW method, the joint water-air phase unphysically mixes due to numerical diffusion, which is a common weakness of finite-differences methods. The simulation results performed with the LS method (Figure 3) agree qualitatively quite well with the images of the experiment (Figure 2). \\\\\wm\\\ 0 10 20 30 40 0 10 20 30 40 0 10 20 30 40 0 10 20 30 40 t = 660 ms t = 710 ms 0 10 20 30 40 t = 760 ms Sl. 3. Gostota združene faze voda - zrak, fazna verjetnost kroglic in hitrostno polje združene faze voda - zrak ob različnih časih pri izračunih z metodo NF Fig. 3. The joint water-air phase density, the spheres’ phase-presence probability and the joint water-air phase velocity field at different times for the LS-method calculation {\\\\\w,m\\ t = 560 ms 0 10 20 30 40 0 10 20 30 40 t = 610 ms 0 10 20 30 40 t = 660 ms t = 710 ms 0 10 20 30 40 0 10 20 30 40 t = 760 ms Sl. 4. Gostota združene faze voda - zrak, fazna verjetnost kroglic in hitrostno polje združene faze voda - zrak ob različnih časih pri izračunih z metodo VL Fig. 4. The joint water-air phase density, the spheres’ phase-presence probability and the joint water-air phase velocity field at different times for the HR-method calculation VH^tTPsDDIK stran 454 120 120 120 120 120 100 100 100 100 100 80 80 80 80 80 60 60 60 60 60 40 40 40 40 40 20 20 20 20 20 0 0 0 0 120 100 80 60 40 20 0 Leskovar M. - Mavko B.: Simuliranje izotermi~nega - Simulation of an Isothermal 0 10 20 30 40 0 10 20 30 40 0 10 20 30 40 0 10 20 30 40 0 10 20 30 40 t = 560 ms t = 610 ms t = 660 ms t = 710 ms t = 760 ms Sl. 5. Gostota združene faze voda - zrak, fazna verjetnost kroglic in hitrostno polje združene faze voda - zrak ob različnih časih pri izračunih z metodo PM Fig. 5. The joint water-air phase density, the spheres’ phase-presence probability and the joint water-air phase velocity field at different times for the UW-method calculation Pri vseh simuliranjih smo spremljali tudi potek tlaka na različnih globinah. Na sliki 6 so prikazane izračunane tlačne krivulje skupaj s preskusnimi podatki. Tlaka p3 in p6 so merili v vodi ob steni posode v točkah 838 mm in 250 mm nad dnom posode. Porast tlaka 0,5 s po izpustu kroglic označuje začetek prodiranja čelnih kroglic v vodo, tlačni vrh po 0,76 s pa zapiranje plinskega stebra. Pri prvem nizu simuliranj smo vse tri faze, to je vodo, zrak in kroglice, obravnavali kot nestisljive. Kakor je razvidno s slike 6, se ob zapiranju plinskega stebra pojavi tlačna konica, ki je najbolj izrazita pri izračunu z metodo NF, pri kateri ostane medfazna ploskev voda - zrak ostra, medtem ko se tlačna nihanja, ki se pojavijo pri preskusu, ne razvijejo. To kakovostno drugačno obnašanje tlaka kakor pri preskusu je posledica nestisljive obravnave zraka, saj je razlog za tlačna nihanja pri preskusu širjenje in krčenje velikega plinskega mehurja, ki se oblikuje po zaprtju plinskega stebra. V drugem nizu simuliranj smo zato fazo zraka obravnavali kot stisljivo. Ker energijske enačbe nismo vključili v model, je bilo treba za stisljivost zraka določiti ustrezno vrednost. Stisljivost zraka je omejena z največjo, izotermno stisljivostjo in najmanjšo, adiabatno stisljivostjo. Preprosta ocena pokaže, da je prevod toplote v plinski mehur v enem nihajnem času zanemarljiv in da je zato adiabatna stisljivost zraka zelo dober približek: During all the simulations the pressure at different levels was also traced. The calculated pressure curves are presented together with the experimentally measured values in Fig. 6. The pressures p3 and p6 were measured in the water at the vessel’s wall, 838 mm and 250 mm from the bottom of the vessel. The pressure rise 0.5 s after the spheres’ release indicates the entry of the first spheres into the water and the pressure peak after 0.76 s is the collapse of the gas chimney. In the first set of simulations, all three phases – the water, the air and the spheres’ phase – were treated as incompressible. As seen in Fig. 6, the pressure produces a spike when the gas chimney collapses, most pronounced for the LS method, where the water-air surface remains sharp, whereas the pressure oscillations, which were observed in the experiment, do not occur. This qualitatively different pressure behavior is due to the incompressible air treatment, since the reason for the pressure oscillations in the experiment is the expansion and compression of the whole, big, gas bubble, which forms after the gas-chimney collapse. In the second series of simulations the air phase was therefore treated as compressible. Since the energy equation was not taken into account the correct air compressibility, which lies somewhere in the range bounded with the maximum isothermal and minimum adiabatic value, had to be determined. A simple estimation shows that the heat conduction in the gas bubble during one oscillation period is negligible and that, consequently, the adiabatic air compressibility is a very good approach: DT la Dto DT a = a Dt osc 10-2W/Km-10-1s 1kg/m3-103J/kgK.10-2m =10-4 <<1 (10), | IgfinHŽslbJlIMlIgiCšD I stran 455 glTMDDC 120 120 120 120 120 100 100 100 100 100 80 80 80 80 80 60 60 60 60 60 40 40 40 40 40 20 20 20 20 20 0 0 0 0 2 Leskovar M. - Mavko B.: Simuliranje izotermi~nega - Simulation of an Isothermal kjer so: AT razlika temperature zraka v mehurju in temperature okolišne vode, AT sprememba temperature zraka v mehurju v enem nihajnem času At zaradi prevoda toplote in Ala linearna izmera mehurja. Kljub temu smo simuliranja opravili tudi z izotermno stisljivostjo zraka, predvsem zato, da bi spoznali vpliv stisljivosti zraka na rezultate simuliranj. Na sliki 6 je prikazan potek tlaka za oba stisljiva primera, pri katerih upoštevamo izotermno oz. adiabatno stisljivost zraka. Kakor smo pričakovali, se v stisljivih primerih pojavijo tlačna nihanja. Tlačna nihanja so najbolj izrazita pri metodi NF, pri kateri ostane medfazna ploskev voda - zrak ostra. nestisljivo / incompressible 1,35 1,3 1,25 1,2 1,15 1,1 1,05 1 Q08-p6 Q08-p3 LS-p6 LS-p3 HR-p6 HR-p3 UW-p6 UW-p3 where AT is the temperature difference between the air in the bubble and the surrounding water, AT is the temperature change of the air in the bubble during one oscillation period At due to heat conduction and Al is the linear dimension of the bubble. Despite that the simulations were also performed with the isothermal air compressibility, primarilly to establish the air compressibility’s influence on the simulation results. The pressure behavior in the compressible cases taking into account the isothermal and adiabatic air compressibility are presented in Fig. 6. As expected, in the compressible cases the pressure oscillations occur and are most pronounced for the LS method, where the water-air surface remains sharp. stisljivo - izotermno / compressible - isothermal --------------------------------- - - - - Q08-p6 Q08-p3 LS-p6 LS-p3 HR-p6 HR-p3 UW-p6 UW-p3 0,5 0,55 0,8 0.5 0,6 0,65 0,7 0, čas / time (s) stisljivo - adiabatno / compressible - adiabatic 1,25 --------------------------------------------------------------------------------------- 0.55 0.6 0.65 čas / time (s) 0.7 0.75 0.8 1,2 1,15 - 1,1 1,05 -¦ 1 -- 0,95------ 0,45 Q08-p6 Q08-p3 LS-p6 LS-p3 HR-p6 HR-p3 UW-p6 UW-p3 A 0,5 0,55 0,7 0,75 0,8 0,6 0,65 čas / time (s) Sl. 6. Primerjava tlakov na različnih globinah (p3 in p6) pri izračunih z metodami NF, VL in PM za nestisljivo, izotermno in adiabatno stisljivost zraka s podatki preskusa Q08 Fig. 6. Comparison of pressures in the vessel at different levels (p3 and p6) between LS-, HR- and UW- methods calculations for the incompressible, isothermal and adiabatic air compressibility and the experimental measurements Q08 Da bi simuliranja lažje primerjali med seboj, smo glavne rezultate simuliranj zbrali v preglednici 1. Za vsa predstavljena simuliranja so navedeni naslednji podatki: največji tlak p6 (p ), čas po katerem tlak doseže vrh (t ) in časovni razmik med prvim tlačnim vrhom in dolom pri tlačnem nihanju (At c/2). Zaradi popolnosti in kot referenca so dodane tudi preskusne vrednosti, ki pa niso namenjene za resno primerjavo, saj so bile preskusne razmere, kar je opisano v 2. poglavju, nekoliko drugačne. Potek tlaka se pri izračunih z metodami PM, VL in NF razlikuje zaradi različno izrazite razpršitve medfazne ploskve voda - zrak, ki jo povzročajo te numerične metode. Metoda PM prvega reda natančnosti povzroči največjo razpršitev (sl. 5) in zato To make the comparison of the simulations easier the main results are assembled in table 1. The maximum pressure p6 (pmax), the time when the pressure peak is reached (tpmax) and the time between the first pressure maximum and minimum during the pressure oscillations (Dtosc/2) are listed for all the presented simulations. The experimental values are also added for completeness and as a reference, but not for a rigorous comparison since the experimental conditions were somewhat different, as explained in section 2. The reason for the various pressure behaviors for the UW-, HR- and LS-method calculations is the different extent of water-air surface-spreading these numerical methods produce. The first-order-accurate UW method produces the highest spreading (Fig. 5), so the VH^tTPsDDIK stran 456 Leskovar M. - Mavko B.: Simuliranje izotermi~nega - Simulation of an Isothermal Preglednica 1. Glavni rezultati nestisljivih in stisljivih simuliranj, opravljenih z metodami NF, VL in PM, skupaj s podatki preskusa Q08 Table 1. Main results of the incompressible and compressible simulations performed with the LS-, HR- and UW-methods together with the experimental measurements Q08 pmax (ba t p max (s) ) — &tosc/2(s) nestisljivo incompressible NF/LS VL/HR PM/UW NF/LS VL/HR PM/UW 1,328 1,098 1,060 1,143 1,108 1,067 J----- Simuliranje Simulation stisljivo - izotermno compressible - isothermal 0,740 0,730 0,716 0,756 0,747 0,725 0,031 0,032 0,036 stisljivo - adiabatno compressible - adiabatic NF/LS 1,155 0,753 0,027 VL/HR 1,111 1,066 1,219 0,028 PM/UW 0,743 0,720 0,760 0,028 preskus experiment Q08 0,030 se začne zračni steber najhitreje zapirati. Ker se gostota medfazne ploskve voda - zrak spreminja najpočasneje, je povišanje tlaka najmanjše. Pri metodi VL je medfazna ploskev voda - zrak ostrejša (sl. 4) in zato je amplituda tlačnih nihanj večja. Tlačne spremembe so največje pri metodi NF, saj ostane medfazna ploskev voda - zrak ostra med celotnim simuliranjem (sl. 3), tako kakor je vidno tudi na posnetkih preskusa (sl. 2). Kakor smo pričakovali, vpliva stisljivost zraka neposredno na frekvenco tlačnih nihanj. Pri manjši, adiabatni stisljivosti zraka je frekvenca tlačnih nihanj višja (preglednica 1). 4 SKLEP Predstavili smo izvirno obravnavo izotermnih preskusov mešalne faze eksplozije pare. Značilnost izotermnih mešalnih preskusov je, da jih ni mogoče ustrezno obravnavati niti izključno z modeli proste površine niti izključno z modeli večfaznega toka, saj se kroglice razpršijo, medtem ko ostane medfazna ploskev voda - zrak ostra. Zato smo se odločili razviti kombinirani večfazni model, pri katerem kroglice kot običajno obravnavamo z modelom večfaznega toka, medfazno ploskev voda - zrak pa z modelom proste površine. Izvirni kombinirani večfazni model smo preverili na izotermnem preskusu mešalne faze eksplozije pare Q08, ki so ga izvedli na napravi QUEOS. Opravili smo številna primerjalna simuliranja, pri katerih smo medfazno ploskev voda - zrak določevali z različnimi numeričnimi postopki (NF, VL in PM), pri tem pa zrak obravnavali kot stisljiv oz. nestisljiv Po pričakovanju so rezultati simuliranj pokazali, da ostane medfazna ploskev voda - zrak ostra le pri izračunih z metodo NF in da je razpršitev medfazne ploskve voda - zrak največja pri izračunih z metodo PM prvega reda natančnosti. Posledično so tlačne spremembe največje pri metodi NF, kjer se gostota na medfazni ploskvi voda - zrak najhitreje spremeni. Tlačna nihanja, ki so se pojavila pri preskusu, smo dobili le v stisljivem primeru, saj so posledica nihanja ujetega plinskega mehurja. Ob upoštevanju predstavljenih pomanjkljivosti modela air chimney starts to close the fastest in this calculation. Since the density of the spread water-air surface changes slowest the pressure rise is the lowest. The HR method produces a sharper water-air interface (Fig. 4); therefore, the amplitude of the pressure oscillations is greater. The pressure changes are the largest for the LS method since it keeps the water-air interface sharp during the whole simulation (Fig. 3), as was also observed in the experiment (Fig. 2). As expected, the air compressibility has a direct influence on the pressure oscillation frequency. At the lower adiabatic air compressibility the pressure oscillation frequency is higher (Table 1). 4 CONCLUSION An original treatment of isothermal steam-explosion premixing experiments was presented. The special feature of isothermal premixing experiments is that they cannot be adequately modeled just with free-surface models or just with multiphase flow models, since the spheres disperse, whereas the water-air surface remains sharp. So we decided to develop a combined multiphase model, where the spheres are treated, as is usual, with a multiphase flow model, whereas the water-air surface is treated with a free-surface model. The original combined multiphase model was validated using the QUEOS isothermal premixing experiment Q08. A number of simulations were performed in a comparison using different numerical methods for the water-air surface determination (LS, HR and UW) treating the air as compressible or incompressible. As expected, the results of the simulations showed that the water-air surface remains sharp only for the LS-method calculations and that the first-order-accurate UW method produces the highest water-air surface spreading. Consequently, the pressure variations are the largest for the LS method, where the density at the water-air interface changes most rapidly. The pressure oscillations that were observed in the experiment could only be reproduced in the compressible case since they are the consequence of the entrapped gas bubble oscillations. The simulation results obtained with the LS gfin^OtJJlMlSCSD 02-8 stran 457 |^BSSITIMIGC Leskovar M. - Mavko B.: Simuliranje izotermi~nega - Simulation of an Isothermal preskusa, se rezultati simuliranj z metodo NF tako kakovostno kakor količinsko dobro ujemajo s podatki preskusa. Zahvala Avtorji se zahvaljujejo Ministrstvu za šolstvo, znanost in šport, ki je finančno podprlo raziskavo v okviru raziskovalnega projekta Z2-3514. method agree, when the deficiencies of the presented experiment model are considered, qualitatively and quantitatively well with the experimental measurements. Acknowledgment The authors gratefully acknowledge the support of the Ministry of Education, Science and Sport of the Republic of Slovenia in the frame of the Research Project Z2-3514. 5 LITERATURA 5 REFERENCES [I] Turland, B.D., G.P. Dobson (1996) Molten fuel coolant interactions: A state of the art report. European Commission, Nuclear Science and Technology, Luxembourg, ISSN 1018-5593. [2] Leskovar, M., J. Marn (1999) A combined single-multiphase flow formulation of the premixing phase using the level set method. Nuclear Energy in Central Europe ‘99, Portorož, Slovenija, Proceedings, 233-240. [3] Marn, J., M. Leskovar (1996) Simulation of steam explosion premixing phase using probabilistic multiphase flow equations. Fluid Mechanics Research, Vol. 22 (1), 44-55. [4] Fletcher, D.F, P.J. Witt (1996) Numerical studies of multiphase mixing with application to some small-scale experiments. Nuclear Engineering and Design, 135-145. [5] Leskovar, M., B. Mavko (2002) Izviren kombiniran večfazni model mešalne faze eksplozije pare. Strojniški vestnik, Vol. 48(2002)8,438-448. [6] Ferziger, J.H., M. Peric (1996) Computational Methods for Fluid Dynamics. Springer-Verlag, Berlin, Heidelberg. [7] Meyer, L., G. Schumacher (1996) QUEOS, a simulation-experiment of the premixing phase of a steam explosion with hot spheres in water, base case experiments. FZKA Report 5612, Forschungszentrum Karlsruhe. [8] Leskovar, M., J. Marn, B. Mavko (2000) The influence of the accuracy of the numerical methods on steam-explosion premixing-phase simulation results. Journal of Mechanical Engineering, Vol. 46, 607-621. [9] Leveque, R.J. (1992) Numerical methods for conservation laws. Lectures in Mathematics, ETH Zuerich, Birkhauser Verlag, Basel. [10] Leskovar, M., J. Marn, B. Mavko (2000) Numerical analysis of multiphase mixing - comparison of first and second order accurate schemes. Fluid Mechanics Research, Vol. 27, 1-32. [II] Sethian, J.A. (1998) Level set methods. Cambridge University Press, Cambridge. Naslov avtorjev: dr. Matjaž Leskovar profdr. Borut Mavko Institut “Jožef Stefan” Jamova 39 1000 Ljubljana matjaz.leskovar@ijs.si borut.mavko@ijs.si Authors’ Address: Dr. Matjaž Leskovar Prof.Dr. Borut Mavko “Jožef Stefan” Institute Jamova 39 1000 Ljubljana, Slovenia matjaz.leskovar@ijs.si borut.mavko@ijs.si Prejeto: Received: 6.8.2001 Sprejeto: Accepted: 20.9.2002 VH^tTPsDDIK stran 458