Group Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Title: 14.4.5 - Two Phase Cross Jet in a Fuel Rod Assembly Using DNS/Level-Set Method
ALL VOLUMES CITATION THUMBNAILS PAGE IMAGE ZOOMABLE
Full Citation
STANDARD VIEW MARC VIEW
Permanent Link: http://ufdc.ufl.edu/UF00102023/00354
 Material Information
Title: 14.4.5 - Two Phase Cross Jet in a Fuel Rod Assembly Using DNS/Level-Set Method Computational Techniques for Multiphase Flows
Series Title: 7th International Conference on Multiphase Flow - ICMF 2010 Proceedings
Physical Description: Conference Papers
Creator: Behafarid, F.
Bolotnov, I.A.
Podowski, M.Z.
Jansen, K.E.
Publisher: International Conference on Multiphase Flow (ICMF)
Publication Date: June 4, 2010
 Subjects
Subject: two-phase flow
stabilized finite element
level set method
DNS
 Notes
Abstract: One possible accident scenario in a hexagonal fuel assembly of a Gen. IV sodium-cooled fast reactor (SFR) is associated with local heat degradation around the hottest fuel element due to the coolant flow blockage by foreign objects and/or mechanical failure of the hottest fuel rods. As a result, high pressure fission gases would be injected into the liquid sodium coolant channels as a high speed two-phase cross-flow jet. A model of gas-jet/liquid-sodium flow evolution has been modeled using the Parallel Hierarchic Adaptive Stabilized Transient Analysis (PHASTA) software. This new model of turbulent, unsteady, three-dimensional two-phase flow uses unstructured grids and multiple (more than 8000) processors. The PHASTA code uses a direct numerical simulation (DNS) method and is capable of capturing the shape of gas/liquid interfaces using the level set method. An intrinsic flaw was found in the standard approach to the modeling of gas/liquid/solid interfaces using the level set method, and appropriate steps to correct it are discussed in the paper. It has also been demonstrated that the pressure outflow boundary conditions (BCs) commonly used in single-phase flow simulations are not applicable to two-phase flows with resolved gas/liquid interfaces, since they cannot handle the surface-tension-induced sharp pressure gradients across the interface. A new approach to treating the outflow BC in two-phase flow problems has been proposed.
General Note: The International Conference on Multiphase Flow (ICMF) first was held in Tsukuba, Japan in 1991 and the second ICMF took place in Kyoto, Japan in 1995. During this conference, it was decided to establish an International Governing Board which oversees the major aspects of the conference and makes decisions about future conference locations. Due to the great importance of the field, it was furthermore decided to hold the conference every three years successively in Asia including Australia, Europe including Africa, Russia and the Near East and America. Hence, ICMF 1998 was held in Lyon, France, ICMF 2001 in New Orleans, USA, ICMF 2004 in Yokohama, Japan, and ICMF 2007 in Leipzig, Germany. ICMF-2010 is devoted to all aspects of Multiphase Flow. Researchers from all over the world gathered in order to introduce their recent advances in the field and thereby promote the exchange of new ideas, results and techniques. The conference is a key event in Multiphase Flow and supports the advancement of science in this very important field. The major research topics relevant for the conference are as follows: Bio-Fluid Dynamics; Boiling; Bubbly Flows; Cavitation; Colloidal and Suspension Dynamics; Collision, Agglomeration and Breakup; Computational Techniques for Multiphase Flows; Droplet Flows; Environmental and Geophysical Flows; Experimental Methods for Multiphase Flows; Fluidized and Circulating Fluidized Beds; Fluid Structure Interactions; Granular Media; Industrial Applications; Instabilities; Interfacial Flows; Micro and Nano-Scale Multiphase Flows; Microgravity in Two-Phase Flow; Multiphase Flows with Heat and Mass Transfer; Non-Newtonian Multiphase Flows; Particle-Laden Flows; Particle, Bubble and Drop Dynamics; Reactive Multiphase Flows
 Record Information
Bibliographic ID: UF00102023
Volume ID: VID00354
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: 1445-Behafarid-ICMF2010.pdf

Full Text



7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010

Two Phase Cross Jet in a Fuel Rod Assembly Using DNS/Level-Set Method

F. Behafarid, I.A. Bolotnov, M.Z. Podowski and K.E. Jansen

Center for Multiphase Research,
Rensselaer Polytechnic Institute,
110 8th St., Troy, NY 12180, USA
behaff@rpi.edu


Keywords: Two phase flow, stabilized finite element, level set, DNS



Abstract

One possible accident scenario in a hexagonal fuel assembly of a Gen. IV sodium-cooled fast reactor (SFR) is associated with
local heat degradation around the hottest fuel element due to the coolant flow blockage by foreign objects and/or mechanical
failure of the hottest fuel rods. As a result, high pressure fission gases would be injected into the liquid sodium coolant
channels as a high speed two-phase cross-flow jet.
A model of gas-jet/1iquid-sodium flow evolution has been modeled using the Parallel Hierarchic Adaptive Stabilized
Transient Analysis (PHASTA) software. This new model of turbulent, unsteady, three-dimensional two-phase flow uses
unstructured grids and multiple (more than 8000) processors. The PHASTA code uses a direct numerical simulation (DNS)
method and is capable of capturing the shape of gas/1iquid interfaces using the level set method.
An intrinsic flaw was found in the standard approach to the modeling of gas/1iquid/solid interfaces using the level set method,
and appropriate steps to correct it are discussed in the paper. It has also been demonstrated that the pressure outflow boundary
conditions (BCs) commonly used in single-phase flow simulations are not applicable to two-phase flows with resolved
gas/1iquid interfaces, since they cannot handle the surface-tension-induced sharp pressure gradients across the interface. A
new approach to treating the outflow BC in two-phase flow problems has been proposed.


V, Volume of each phase within a cell (m )
w Pseudo velocity in re-distancing

Greek letters
e Half of the thickness of the interface (m)
pu Dynamic viscosity (kg/m s)
p Density (kg/m )
r Pseudo time in re-distancing

E.Stress tensor (N/m )

Subscripts
1 first phase: Liquid
2 Second phase: Gas
i or j Coordinates (Defines spatial components)
,i or j Spatial derivative in that direction
,t Time derivative


PHASTA Description
PHASTA is a parallel, hierarchic (higher order accurate,
from the 2nd to 5th order accuracy, depending on function
choice), unstructured/adaptive, stabilized finite element
transient analysis flow solver for both incompressible and


Introduction

Development of next generation nuclear reactors, in
particular Gen. IV reactors, will require a new generation
of advanced modeling tools. Thus, it is important that
advanced computational capabilities be developed to
capture the multi-scale physics phenomena governing
reactor transients and accidents.
The purpose of this paper is to discuss a novel
mechanistic approach to 3-D simulations of the injection
of gaseous fission products into a partially blocked SFR
coolant channel following locahized cladding overheat and
failure. The proposed methodology has been developed
for, and implemented in, the PHASTA DNS code,
coupled with the Level Set Method (LSM) solver.
The computational mesh used in the simulation is highly
refined around the injection zone, to let PHASTA capture
the effect of turbulence and accurately predict the
evolution of interface shape.

Nomenclature
d Corrected distance field (m)
SBody forces (N)
H, Heaviside kernel function
p Pressure (N/m )
S, Strain rate tensor (s')
ui Velocity (m/s)






7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010

incompressible flow of a Newtonian fluid, the stress
tensor is related to the strain rate tensor, SU, as:
.i. = 2pSij.= p(ut.j+ ujp)




PHASTA
domain












Figure 2. Top view of the PHASTA computational
domain.

Using the Continuum Surface Tension (CST) model of
Brackbill et al. (1992), the surface tension force is
computed as a local interfacial force density.
The finite element formulation is based on the so called
weak torm of equations which are generated by dotting
continuity and momentum equations on the left by the
weight functions and integrating over the domain.
Following the stabilized formulation of Taylor et al
(1998) the residuals of the continuity and momentum
equations are added together.
The stabilized finite element formulation is described in
detail in Rodriguez (2009) and Jansen et al (1999). The
stabilization and its parameters were described by Franca
and Frey (1992) and by Whiting and Jansen (2001).

Level Set Method

The level set method of Sussman (Sussman et al, 1998,
Sussman and Fatemi, 1999 and Sussman et al., 1999) and
Sethian (1999) involves modeling the interface as the zero
level set of a smooth function, 9, which represents the
signed distance from the interface. Hence, the interface is
defined by p = 0. The scalar 9, the so-called first scalar, is
convected within a moving fluid according to,

D =fv +.V= 0

Phase-1, typically the liquid phase, is indicated by a
positive level set, p > 0, and phase-2 by a negative level
set, p < 0. Evaluating the jump in physical properties
using a step change across the interface leads to poor
computational results. Therefore, the properties near an
interface are defined using a smoothed Heaviside kernel
function, H,, given by:




HE ( ) = 1+ + i



The fluid properties are thus defined as:


compressible flows (Jansen et al., 2000, Whiting and
Jansen, 2001)
PHASTA was the first unstructured grid LES code
(Jansen, 1993) and has been applied to turbulent flows
ranging from validation benchmarks (channel flow, decay
of isotropic turbulence) to complex flows (airfoils at
maximum lift, flow over a cavity, near lip jet engine flows
and fin-tube heat exchangers).
PHASTA uses advanced anisotropic adaptive algorithms
(Sahni et al., 2006) and the most advanced LES/DES
models (Tejada-Martinez and Jansen, 2005). Note that
DES, LES, and DNS are computationally intensive even
for single phase flows. The two-phase version of
PHASTA utilizes the Level Set method to define the
interface between the gas and liquid phases (Sethian,
1999). Here, it has been used to simulate gas jet injection
into the coolant region.

Simulation Domain

PHASTA simulations were performed for a small region
in front of the crack. Figures 1 and 2 show the
computational domain sharing in the multi-scale approach
used in the simulations.
PHASTA outlet two-phase
mixture outflow


Inlet
velocity
profile


Figure 1. Computational domain (side view).


PHASTA Formulation
Governing Equations
PHASTA is a parallel, hierarchic, unstructured adaptive,
stabilized (finite element) transient analysis flow solver
(both incompressible and compressible). The spatial and
temporal discretization of the Incompressible NS
equations within PHASTA has been described in Whiting
(1999) and Nagrath !1lI -H~ This section follows their
discussion.
The strong form of the INS equations is given by:
Continuity: ui j 0

momentum: pui~t +pu ui, j ,i +9r, j i

where ui is ith component of velocity and fi represents
body forces along the ith coordinate. For the






7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010

where H' is aHE (d) Bd and is given by:



HE()=11 1 a


Hence, to minimize the deviation in volume each new di~
is constrained to satisfy the following:
J] H'/(do)(dk Odo x= 0

This is accomplished by projecting the value of the
current level set, dy onto the new level set field, dio ,

where di is computed from:

di = di + f fk () )HE (d )
The scaling function /1U is assumed constant within each
cell and is computed from:

HE (d ) ~k ()d
k O


p(#) = PlHE ~ P2(1- HE ())

pU( ) = pl1HE(~ 2 (1- HE ())

Although the solution may be reasonably good in the
immediate vicinity of the interface, the distance field may
not be correct throughout the domain since the varying
fluid velocities throughout the flow field distorts the level
set contours. Thus the level set is corrected with a re-
distancing operation by solving the following PDE:
ad
= S( ) 1- Vd

where d is a scalar that represents the corrected distance
field and r is the pseudo time over which the PDE is
solved to steady-state. This may be more clearly
expressed as the following transport equation:


+ E:~Vd = S(4)


The second scalar, d, is originally assigned the level set
field, p, and is convected with a pseudo velocity Ei,
where
Vd
G = S( ) Vd

and S(0) is defined as:



S() )= + -]1si ,<



Note that the zeroth level set or interface, p = 0, does not
move since its convecting velocity Ei is zero. Solving the
second scalar to steady-state restores the distance field
to Vd = f-1but does not alter the location of the interface.
First scalar, 9, is updated using the steady solution of d-

Volume Constraint on Level Set Re-distance Step
Sussman proposed an additional constraint to be applied
during the re-distance step to help ensure the interface
(p=0) does not move (Sussman and Fatemi, 1999 and
Susman et al., 1999). They found imposing this constraint
improves the convergence of the re-distance step while
maintaining the original zero level set. The essence of the
constraint is to preserve the original volume of each phase
during the re-distance step. The volume of each phase
within a cell, V -, is defined as:



where dk is the value of the distance field, d, at iteration k.
Small changes in the zero level set from pseudo time step
O to Zk may be approximated by:
k () k () He(d )

k a
= J~ H:(d ) d -d )dv
e e


"j HE(d ~) dx

Note when the distance field has converged,


dk = ,) and


/1U= 0. In incompressible PHASTA, the volume constraint
correction occurs at the global level after each re-
distancing iteration. The calculation of /1U occurs at the
element level and is projected back to the global level.
The Level set method was implemented into PHASTA by
Nagrath !II II -1, and Rodriguez !II II **,s For the simulations
presented herein the authors used the version of the code
provided by Dr. Rodriguez with several modifications
described below.


LSM Shadowing Problem
Level set method models the exact shape of the interface
as the zero level set of a smooth function which represents
the signed distance from the interface. The convection of
the distance field by the liquid flow moves the interface
properly but distorts the distance field. To recover it, a
few re-distancing (or re-initiating) iterations should be
applied after each time step. Re-distancing has an intrinsic
flaw which can cause errors and introduce a false phase in
some parts of the domain.
This problem is related to the lack of boundary conditions
for distance field and has not been discussed in the
literature before. In re-distancing, the information flows in
the direction normal and away from both sides of the
interface. If the interface is not attached to the boundary
or is attached in the normal direction to it (see Figure 3),
the level set method does not cause errors.
However, when the interface is attached to the wall
(creating a triple point) at an arbitrary angle, there are
some regions that cannot get the correct information about
the distance field due to the location of those regions in
relation to the interface and to the wall. We call these
regions "shadowed". One example of a shadowed region
is shown in Figure 4.









If angle is 900 or less, there is

O ~no shadowing problem






Figure 3. Schematic of conditions when LSM has no
problem even if the interface touches the boundary
In the shadowed regions, the information should come
from the boundary which has no boundary condition for
the distance field; therefore, meaningless information may
propagate into, and pollute, the domain. In this case, more
re-distancing iterations make the distance field worse-
















Figure 4. Schematic of two typical conditions when
interface: touches the boundary (Solid) and LSM
encounters the shadowing problem.


7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010


\h,,,1. -as. J region


Figure 6. Distance filed after re-distancing. The distance
field in shadowed region is getting worse.

In most cases it is difficult to set a boundary condition for
the distance field in order to prevent the shadowing
problem. The current approach is to minimize the impact
of this error by using the limited number of re-distancing
operations.
Future work will include a level set method enhancement
which will allow the enforcement of a physically based
value of a contact angle at the triple point.

New BC for Two Phase Flow Outlet

Three types of BCs commonly used for outflow for
incompressible fluids are "Natural Pressure BC",
"Convective Pressure BC" and "No Pressure BC". None
of the three can properly handle the lateral pressure
gradient due to surface tension force in two phase flow
simulations.
The common way to avoid problems in outflow is to use
periodic boundary conditions. The issue can also be
avoided in cases where the pressure jump through
interface is low; either due to a large interface curvature
or a small value of the surface tension. In such cases, the
"Natural Pressure Boundary Condition" can handle this
small pressure difference with little or no effect on the
interface far from the outlet. In other similar cases, the
main stream velocity is high enough to convect the
unphysical data out of the domain. In the case of interest
to this work, none of the above mentioned conditions are
satisfied. Since no solutions are available in the open
literature a new method has been developed to solve the
problem.
A simple explanation of the problem is as follows. The
surface tension force acts on, and is balanced by the
pressure gradient across, the interface. At the outlet, the
"Natural Pressure BC" tries to make the pressure uniform
throughout the outlet and eliminates the pressure
difference. Therefore, the surface tension cannot be
balanced by pressure difference and moves the interface
in a non-physical manner. Figures 7 9 show a co-flow
case with the "natural pressure outlet BC In order to
solve this problem, we added an artificial body force at
each nodal point in a very thin layer adjacent to the outlet,
to maintain the surface tension/pressure gradient balance.
The simulation results mn this region can be considered as
a numerical extension to the real domain and have no
physical meaning. Smnce this region is very thmn (of the
order of 0.1 mm), the additional computational expenses
are negligible. To properly apply the body forces, this thin
layer should have at least three to five elements across.


Figure 5. Distance field after convection by flow and
before re-distancing. The arrows show the direction of
information, normal to the interface used to re-build the
distance field.

Figures 5 and 6 show a shadowing problem observed in
the current PHASTA simulations. In Figure 5, the
distance field is shown after it has been convected by the
fluid flow solution. At this stage the interface has moved
and the distance field is distorted. A series of re-
distancing iterations is applied to reconstruct the distance
field.
Figure 6 shows the distance field after re-distancing
iterations have been performed. As discussed above, the
distance field degrades after re-distancing in the
shadowed regions only, where the distance field
information flows from the boundary to the interior of the
domain. In this region, no matter how many re-distancing
iterations are used, the distance field cannot be recovered.































I tI I I I


Figure 10. Schematic of the artificial body forces applied
to modify the outlet BC when theof plume of gas exits the
computational domain.

Outlet: Natural Pressure
Inlet: Velocity poi
0s~i 0200ir

00150"










Figure 11. Test case geometry and velocity profile.

velocity Magnitude
000 0.0500 0.100 0.150 0.200


7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010


SPressume gradient
Surface tension
SInterface
SBody forces inside of the interface
SBody forces in the gas region


Artificial body
forces are added
after this line


SOutlet


Figure 7. Geometry colored by velocity magnitude. The
small opening in the middle reflects gas injection region.


Outlet
Natural Pressure


Liquid bIlet 7


Gas bIlet


Figure 8. Pressure field for the co-flow case before
gas/1iquid interface reaches domain outlet.


Figure 9. Effect of surface tension on pressure imbalance
when the interface reaches domain outlet.

Figure 10 shows a schematic of the plume of gas reaching
the outlet and shows the size and direction of the artificial
body forces. The violet arrows show the body forces
applied inside the interface to balance the difference
between surface tension and pressure gradient. On the
other hand, in the gas region the pressure goes down in
the direction normal to the outlet. To avoid gas
acceleration, another body force has applied which is
shown as black arrows.
Figure 11 shows the geometry, velocity profile and the
boundary conditions of a simple case used to test the new
outlet BC.
Figures 12 and 13 present the results of simulations using
the new modified boundary condition for the case of co-
flow with a cylindrical interface.
Since the initial interface is defined as a perfect steady
cylinder, even very small changes in the interface shape at
the outlet can be readily noticed. As shown in Figure 13,
pressure inside the gas plume goes down to zero but the
artificial body forces balance the surface tension.


Figure 12. Cross section of the velocity field in a co-flow
test case. The initial interface is a cylindrical plume of
gas (white lines) all the way to the outlet.

In test calculations, the interface maintained its original
shape after 2000 time steps. Similar tests were performed
for other interfaces including bubbles and more complex
shapes, and all of them shown desired results.


t t i t 1 i t

c
c
c
c
c
c

1 lit f f 1






7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010

Different gas injection velocities have been used to test
their effect on the interface shape. A parabolic velocity
profile was used at the crack inlet, corresponding to three
different values of maximum velocity: 5 m/s, 15 m/s and
120 m/s.


070s 10


Figure 13. Pressure field after an artificial body force has
been applied to the right of the white line. A uniform
pressure is observed at the outlet, whereas the interface
shape stays unchanged.

Results

Gas injection

The first geometry under consideration, shown in Figure
14, included two fuel rods. Figure 15 shows the 3-D
geometry colored by distance field. The crack has a
rectangular 0.1 mm by 0.1 mm shape, and it is initially
filled with gas (as shown by a negative sign of the
distance field variable).

/ conywasutroal







r


/ *


Figure 14. Computational domain (top view).
The following fluid properties have been used in the
simulations presented below:
kg ; 1431nkg
p, = 754.9 ,:p .8x0 -


p, = 4.481 k ,=4051 g
n2 n2 s
where the liquid sodium properties are taken at 800oC and
the fission gas is assumed to have the properties of
nitrogen. The surface tension coefficient used in the
simulations was, a = 0.19 N/2.


sea cano e
33 01 06


scrl

MO.OO89
0.00430

0.00270

0.00110
-0 000501




Figure 15. Computational domain including a small
crack colored by distance field (having negative values in
the gas region. Initially, the interface was at the
intersection of the crack and the main domain.


Figures 16 18 show the velocity and pressure fields, as
well as the 3-D interface shape for the maximum gas
injection velocity of 5 m/s after three different numbers of
time steps. The left part of each figure shows the velocity
magnitude distribution in a 2-D cut of the computational
domain. The cut plane was chosen to include the
directions of mean jet velocity and mean liquid velocity.
The solid white line shows the gas liquid interface. The
center part of each figure presents the pressure
distribution in the same plane cut. The interface is shown
with a black solid line. Higher pressure values can be
observed in the gas region due to the surface tension
effect. Finally, the right part of each figure presents a 3-D
view of the gas/1iquid interface. The color of the interface
represents the local velocity magnitude.
As discussed before, in these simulations a parabolic
velocity profile has been used at the gas injection
boundary condition. It should be noticed that the
maximum liquid velocity before the gas injection is
almost 10 times lower than the maximum gas injection
velocity, but the effect of gravity is strong enough to shed
bubbles from the inlet gas stream. In this simulation and
the following ones, the gas density was almost 168 times
less than the liquid sodium.
In figures 19 21, the maximum gas injection velocity
increased from 5 m/s to 15 m/s. In this case enough gas is
provided to avoid bubble shedding. Also, the gas flow is
less stable and exhibits somewhat turbulent behavior
compared to the previous case.
Figures 22 and 23 show the results for a maximum gas
injection velocity of 120 n2/s. The gas flow is very
turbulent and starts to expand the interface close to the
wall. For this special geometry and gas injection velocity,
the interface is expected to create a lot of small bubbles.






7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010


essure
1 500+03
938.

375-
-188

750


11.2

7.50


3 7

2 50


11.2


Figure 16. Simulation of two-phase cross jet injection.
The gas and liquid velocity fields are colored by pressure,
via; = 5 m/s, of time step number: 2,000.


Figure 19. Simulation of two-phase cross jet injection.
The gas and liquid velocity fields are colored by pressure,
vi4j= 15 m/s, of time step number: 5,000.









11.2 938. 11.2

7.50 1 = 376./ 7.50











Figure 20. Simulation of two-phase cross jet injection.
The gas and liquid velocity fields are colored by pressure,
vi4j= 15 m/s, time step number: 8,500.









11.2 938. 11.2

7.50 376. 7.50

H .75 -188 37









Figure 21. Simulation of two-phase cross jet injection.
The gas and liquid velocity fields are colored by pressure,
vi,j = 15 m/s, time step number: 11,500.


Figure 17. Simulation of two-phase cross jet injection.
The gas and liquid velocity fields are colored by pressure,
v~i = 5 m/s, time step number: 4,000.


3mi

2 50

25:


Figure 18. Simulation of two-phase cross jet injection.
The gas and liquid velocity fields are colored by pressure,
vin; = 5 m/s, time step number: 5,500.






7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010

coarse region is eight times bigger than in the refined
region, the computational cost of the newly added domain
is about 500 times (8 ) smaller per unit volume.


Figure 24. Top view of the velocity field for the
computational domain including 2 fuel rods, vmax=0.462


Tegsadlqiveoiyfedarcooebypressure,
v;,,y= 120m/s, time step onume:100

a/ I









Figure 23. Simulation of two-phase cross jet injection.
The gas and liquid velocity fields are colored by pressure,
vi,j = 120 m/s, time step number: 2,000.









.


mesh resoltio which cain be tpaffred inos jea pareticula


Inthe prevn ious aeapaaoi a velocity fed r ord proie was ,

s et y asa oundar condition capt the ga tint To obtain a
motrephyical vhnoeloit profle winside the crack, the
comutatisonald doain wtas late extmenddto inclusdiae gas


(oresult 00) Tu, the rc oehsbeent etene areosse the

Iths allso bubeenosre in the simulations thateen uo the
syme etry faces which are cloe tfo the interace patcanaffct


dmai ie h odified geometry ishonnFgue

Iet shoul be menti oned here tat the compultationa ostai of
thre adde reion i veriy smafll since the mehcrank be

loars in the regionnls farfrom the interfacdte whchpertin .


Ito hsingle-hs ben fowFr exape if the smeltos h szei the


Figure 25. Top view of the velocity field for the
computational domain including several fuel rods
surrounding one intact fuel element, vmax= 0.008 m/s.

Azimuthal flow instability
The inclusion of a larger domain also yielded an
interesting finding during the stage of simulation where
the fully developed liquid flow was occurred prior to the
start of gas injection. Since the geometry is quite
complex, no analytical methods could be used to
formulate the initial conditions for the liquid coolant
inflow velocity profile. Instead, a stagnant flow was used
at the starting point, combined with periodic boundary
conditions in the stream-wise coolant direction and a
proper body force to drive the flow and to create the fully
developed velocity profile. The resultant velocity profile
was then used as the initial and coolant inflow conditions.
One of the major differences between the geometries
shown in Figure 24 and Figure 25 is that the latter
geometry does not impose the symmetry boundary
condition which was an inherent part of the former,
As shown in Figures 25 and 26, after the symmetry lines
have been removed and a complete fuel rod was includes
inside the domain, an azimuthal instability have occurred,


0.00592~






7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010

quite a stable behavior. Due to the relatively slow liquid
flow, no bubble shedding is expected for this case. It can
be noticed from the pressure plots in Figures 28 and 29,
that a higher gas pressure pushes the jet out of the crack
(the left part of the central pressure contour plot). Also,
as the size of the gas region in the channel grows, the
incoming gas jet starts to show a more turbulent behavior,
as it is no more constrained by the relatively heavy fluid.
Figure 30 shows the velocity distribution and the interface
evolution at various downstream distances from the
injection zone. It is interesting to notice two nearly
symmetrical recirculation regions in the middle of the gas
bubble in Figure 30 (b, c).


enhancing turbulence development at a much lower mass
flow rate than before. This instability first starts rotating
the fluid in one direction, and then in the opposite
direction, so that the flow eventually becomes turbulent,


10 .0 2


Figure 26. Unstable velocity field around the failed rod,
vmax = 0.07 m/s. A counter-clockwise rotation can be
noticed.

The maximum velocity to create a laminar velocity profile
in Figure 25 is almost 57 times smaller than that of the
Figure 24. This can be explained by the fact that the
effective characteristic diameter of the flow shown in
Figure 25 is much larger due to the enclosed liquid loop
around the fuel rod. This fact results in reaching the
critical Reynolds number at much lower liquid velocities
compared to the case shown in Figure 24. Figure 27
shows the mesh corresponding to the case of a larger
computational domain with the ring. This mesh has
around 15 million elements and more than 8000
processors were used in the simulations.


Figure 28. Side-view cuts of the pressure and velocity
fields and a top view of the velocity field: vinj, max= 8.4
m/s: crack size: 1 by 2mm, time step number: 12,000:
transient duration (physical time): 0.021 s.


Figure 29: Side-view cuts of the pressure and velocity
fields and a top view of the velocity field: vinj, max = 8.4
m/s: crack size: 1 by 2 mm, time step number: 24,000:
transient duration (physical time): 0.043 s.

Summary and Conclusions
A DNS approach with Level Set Method for interface
tracking has been used to numerically simulate a
hypothetical nuclear reactor accident. It was shown that
the standard Level Set Method has an intrinsic problem on
the boundaries near the gas/1iquid/solid interfaces during
the re-distancing stage, a ways to overcome this problem
were discussed. It was also was shown that no common
outlet boundary condition can handle the sharp pressure
jump across the interface due to surface tension force. A


Figure 27. Cross-section of the mesh. The mesh is refined
inside and around the crack, and is coarse far from it.

Figures 28 and 29 show the evolution of gas flow for the
new case, based on the modified geometry and the new
velocity profile at the crack. The gas/1iquid interface is
shown as a white solid line in the left and center part of
the figures. The right parts of Figures 28 and 29 show the
domain slices across the coolant flow direction, with the
interface shown again as a white line.
It has been observed that with a larger size crack and a
relatively small jet velocity, the gas/1iquid interface is
well controlled by the surface tension force and it exhibits


velocity Magnitude


velocity Magn tude


4.2~113


velocity Magiitude
.169


.056


U






7th International Conference on Multiphase Flow,
ICMF 2010, Tampa, FL, May 30 June 4, 2010

Jansen, K. E., Collis, S. S., Whiting, C., Shakib, F., "A
better consistency for low-order stabilized finite element
methods" Computer methods in applied mechanics and
engineering, Vol. 174, 153-170 (1999).
Jansen, K. E., Whiting, C.H. and Hulbert, G.M., "A
generalized-alpha method for integrating the filtered
Navier-Stokes equations with a stabilized finite element
method", Computer Medrods in Applied Mechcaries card
Engineering, Vol. 190, 3-4, 305-319 (2000).
Nagrath, S. "Adaptive stabilized finite element analysis of
multi-phase flows using level set approach", Ph.D. Thesis,
MANE, Rensselaer Polytechnic Institute (2004).
Nagrath, S., Jansen, K.E. and Lahey, Jr., R.T.,
"Computation of incompressible bubble dynamics with a
stabilized finite element level set method", Computer
Medrods in Applied Mechcaries and Engineering, Vol.
194(42-44) 4565-4587, (2005).
Rodriguez J. M., "Numerical simulation of two-phase
annular flow", Ph.D. Thesis, MANE, Troy, NY,
Rensselaer Polytechnic Institute (2009).
Sahni, O., Mueller, J., Jansen, K.E., Shephard, M.S. and
Taylor, C.A., "Efficient Anisotropic Adaptive
Discretization of the Cardiovascular System", Computer
Medrods in Applied Mechcaries and Engineering, Vol.
195, 41-43, 5634-5655 (2006).
Sethian, J. A., "Level Set Methods and Fast Marching
Methods", Cambridge University Press. (1999)
Sussman, M., E. Fatemi, P. Smereka and S. Osher., "An
improved level set method for incompressible two-phase
flows", Journal of Computers & Fluids, Vol. 27, Issue 5-
6, 663-680 (1998).
Sussman, M. and E. Fatemi., "An efficient, interface-
preserving level set re-distancing algorithm and its
application to interfacial incompressible fluid flow."
Siam Journal on Scientific Computing, Vol. 20, Issue 4,
Pages 1165-1191 (1999).
Sussman, M., A. S. Almgren, J. B. Bell, P. Colella, L. H.
Howell and M. L. Welcome., "An adaptive level set
approach for incompressible two-phase flows.", Journal of
Computational Physics, Vol. 148, Issue 1, 81-124 (1999).
Taylor, C. A., T. J. R. Hughes and C. K. Zarins, "Finite
element modeling of blood flow in arteries", Computer
Methods in Applied Mechanics and Engineering, Vol.
158, Issues 1-2, 155-196 (1998).
Tejada-Martinez, A.E. and Jansen, K.E., "A parameter-
free dynamic subgrid-scale model for large-eddy
simulation", Computer Medrods in Applied Mechanics
and Engineering, Vol. 194, No. 9, 1225-1248 (2005).
Whiting, C. H., "Stabilized finite element methods for
fluid dynamics using a hierarchical basis", Ph.D. Thesis,
MEAE, Troy, NY, Rensselaer Polytechnic Institute
(1999).
Whiting, C.H. and Jansen, K.E., "A Stabilized Finite
Element Formulation For The Incompressible Navier-
Stokes Equations Using A Hierarchical Basis,"
International Journal of Numerical Medrods in Fluids,
Vol. 35, 93-116 (2001).


new method was presented to provide a proper outlet BC
for two-phase flow simulations.
The results of the gas injection into the liquid sodium for
different crack sizes and gas velocities were shown and
discussed. The selection of the main domain geometry
and its effect on the simulation were investigated.
Future work will meclude the simulations for extended
domains with time-dependent gas mnflow rates and a time-
varying crack size.


;eocty Magntuc~
PO 22 43 65 8.6


Figure 30: Top views of velocity field: (a) at the center of
the crack, (b) 2 mm, (c) 4 mm, and (d) 6 mm downstream
from the center of the crack. Time step number: 13,000.

Acknowledgments
The authors would like to acknowledge the financial
support provided to this study by the U.S. Department of
Energy and the computation resources provided by the
Computational Center for Nanotechnology Innovations
(CCNI) at Rensselaer Polytechnic Institute (RPI).

References
Brackbill, J. U., D. B. Kothe and C. Zemach, "A
continuum method for modeling surface tension", Journal
of Computational Physics, Vol. 100, Issue 2, Pages 335-
354 (1992).
Franca, L. P. and S. L. Frey. "Stabilized Finite-Element
Methods. II. The Incompressible Navier-Stokes
Equations", Computer Methods in Applied Mechanics
and Engineering, Vol. 99, Issues 2-3, Pages 209-233
(1992).
Jansen, K.E., "Unstructured Grid Large Eddy Simulations
of Wall Bounded Flows", Annual Research Briefs, Center
for Turbulence Research, NASA Anzes/Stcozford
University, 151 (1993).


U U




University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - Version 2.9.7 - mvs