Encyclopaedia IndexBack to start of article

7.5.Modifications for Bubble-Induced Turbulence

There have been a number of proposals for extending turbulence model to account for the additional production of turbulence due to the presence of bubbles.

PHOENICS provides three representative options, which are described below. All of these options require specification of the bubble diameter and/or calculation of the local slip-velocity vector, and are thus designed to be used in conjunction with one of the interphase drag laws selected by CFIPS=GRND7 (see the Encyclopaedia entry 'INTERPHASE DRAG MODELS').

Modifications of the KE-EP model

Following Lopez de Bertodano et al [1990] and Svendsen et al [1992], the first option introduces additional volumetric source terms into the KE-EP model:

SKE,b = R1*Pb (7.5.1)
SEP,b = R1*C1*Pb*EP/KE (7.5.2)

where Pb is the production rate of KE due to the drag work that the bubbles do as they move through the continuous phase:

Pb = Cb*(Fi . Vsi) (7.5.3)

where: Cb (=EL1A) is an empirical constant; Vsi is the slip velocity vector; and Fi is the drag-force vector per unit volume, i.e. Fi = {FIP/Vol}*Vslip, where Vslip is the absolute magnitude of the slip velocity vector, Vol is the cell volume, and FIP in the interphase drag coefficient defined by equation (2.9) under the Encyclopaedia entry 'INTERPHASE DRAG MODELS'.

Thus, equation (5.3) reduces to:

Pb = Cb*0.75*Cd*RHO1*R2*R1*Vslip**3/Db (7.5.4)

where Cd is a dimensionless drag coefficient and Db is the bubble diameter ( defined by CFIPB ). There is much uncertainty concerning value of Cb: for example, Lopez de Bertodano et al [1993] use Cb=0.0 whereas Johansen et al [1988] and Svendsen et al [1992] suggest values ranging between Cb=0.1 and 0.75. PHOENICS sets a representative default value of Cb=0.05, which is suitable for low gas fractions, say << 5%.

The FORTRAN coding for the foregoing sources is contained in subroutine GXDISP of the file GX2PHS.FOR. They can be activated by inserting the following PATCH in the Q1 file:


The PIL variable EL1A permits the user to modify Cb; and when STORE(PRKB) appears in the Q1 file, the product (Pb*Vol) may be printed in the RESULT file or viewed via PHOTON and AUTOPLOT.

Augmentation of the turbulent viscosity

Following Lopez de Bertodano et al [1994], an alternative modification is to calculate ENUT by linear superposition of the shear-induced and bubble-induced viscosities, which for the KE-EP model, leads to:

ENUT = CMUCD*KE**2/EP + CMUB*Db*R2*Vslip .... (7.5.5)

where: CMUB(=ENUTC) is an empirical constant; Db(=CFIPB) is the bubble diameter; and R2 is the volume fraction of phase 2 (the dispersed phase). This augmentation of ENUT is activated by setting ENUTB=1.0 in the Q1 file, and ENUTC permits the user to modify CMUB from the default value of 0.6 recommended by Lopez de Bertodano et al [1994].

The modification also allows implicitly for the appearance of bubble-induced turbulent stresses in the phasic momentum equations ( see the discussion given by Lopez de Bertodano et al [1994] ).

The FORTRAN coding can be found in subroutine GXENUT of the file GXPROP.FOR. The modification can be used with the high-Re KE-LM and KE-EP models, and also with any of the high-Re zero-equation ENUT options in subroutine GXENUT, excepting the option ENUT=GRND1 for which ENUTB is already in use. For the KE-LM and KE-EP models, the user must select KELIN=3 for linearisation of the KE-EP source terms.

Modifications of the mixing-length formula

For bubbly two-phase flows, Rice and Geary [1990] proposed the use of the following mixing-length formula in Prandtl's mixing-length model (see the Encyclopaedia entry 'ENUT'):

LM = (Db/h) * Rd/{Rd,av} (7.5.6)

    where: LM is the mixing length;
           Db (=CFIPB) is the bubble diameter;
           h (=EL1A) is a correction factor for bubble deformation,
              which is taken as unity by default;
           Rd is the void fraction; and
           {Rd,av} is the area-averaged void fraction at the current
           flow section, given by:

{Rd,av} = sum(Rd*dA)/sum(dA) (7.5.7)

where sum denotes a radial summation over all cells at the current cross-section, and dA denotes the cell flow area normal to the radial direction.

It can be seen that the model scales the mixing length with the bubble diameter only, and therefore presumes that bubble-induced turbulence dominates shear-induced turbulence. The model appears to be designed specifically for bubble-column reactors, in which liquid flows up the centre of the column and down near the walls.

Therefore, the PHOENICS implementation presumes that the flow geometry has no draught tube to separate the two liquid streams. The implementation is also restricted to bubble-column geometries in which z is the main flow direction.

The model is activated by TURMOD(MIXLEN-RICE), which is equivalent to
ENUT=GRND2; EL1=GRND10; EL1A=1.0; EL1E=1.0; GENK=T

The bubble diameter is set by the PIL variable CFIPB (see the PHENC entry CFIPS), and the bubble-deformation correction factor h is set by the PIL variable EL1A.