In dispersed two-phase flow, the phases are subject not only to interphase drag forces, but also to interphase lift forces. These forces are particularly important for the prediction of phase- separation and phase-distribution phenomenon ( see for example Lahey [1990, 1993] and Prakash [1987] ).
For example, the lateral phase distribution for fully-developed bubbly flow through a circular duct shows that bubbles accumulate near the walls for upward flow, and near the centre of the duct for downward flow ( see Prakash [1987] and Serizawa et al [1992] ).
The modelling of the interphase-lift forces has been considered by a large number of workers, including: Beyerlein [1985], Drew and Lahey [1987], Prakash [1987], Ciccone et al [1989], Huang [1989], Lee and Chang [1991], Petersen [1992] and Lahey [1990, 1993].
All of these researchers introduce the following source terms into the momentum equations:
Sc = Cl * rhoc * rd * (Ud-Uc) x curl(Uc) (1.1)
Sd = -Cl * rhoc * rd * (Ud-Uc) x curl(Uc) (1.2)
wherein:
x denotes the vector cross product;
Sc and Sd are the volumetric force vectors for the continuous
and dispersed phases, respectively;
Ud and Uc are the corresponding velocity vectors;
Cl is the interphase lift coefficient;
rhoc is the density of the continuous phase; and
rd is the volume fraction of the dispersed phase.
These sources, which have the accumulation effects mentioned above, are represented in PHOENICS.
Drew and Lahey [1987], Huang [1989] and Petersen [1992] recommend that the lift coefficient Clshould be given the same value as Cvm, the virtual mass coefficient ( see the Encyclopedia entry VIRTUAL- MASS MOMENTUM-SOURCE TERMS ).
PHOENICS provides not only for a constant value of Cl, but also for values which very with rd, as follows:
Cl = Cla * rc ; and (1.3)
Cl = Cla * [1.-2.78*min(0.2,rd)] (1.4)
where Cla is a constant, usually 0.5, and rc is the volume fraction of the continuous phase.
Equations (1.3) and (1.4) are identical to the relations provided in PHOENICS for Cvm in the virtual-mass momentum sources.
Lahey [1990, 1993] reports that a value of Cl=0.5 is appropriate only for a single bubble in inviscid flow, and that Cl may be as small as 0.01 for highly viscous flows.
It should be mentioned that other workers ( Huang [1989], Petersen [1992] and Svensen et al [1992] ) have employed negative values of Cl for the simulation of bubble-column reactors.
The interfacial-lift source terms for the momentum equations are activated by insertion in the Q1 file of the PIL statement:
INTSOR(LIFT,CLIFT,CLIFTA)
which is equivalent to:
PATCH(LIFT,CELL,1,NX,1,NY,1,NZ,1,LSTEP)
COVAL(LIFT,U1,TINY,GRND4);COVAL(LIFT,U2,TINY,GRND4)
COVAL(LIFT,V1,TINY,GRND4);COVAL(LIFT,V2,TINY,GRND4)
COVAL(LIFT,W1,TINY,GRND4);COVAL(LIFT,W2,TINY,GRND4)
The value ascribed to the PIL variable CLIFT determines how the lift coefficient Cl is to be calculated, as follows:
| CLIFT | = 0.0 cuts out the interphase-lift terms entirely. | |
| CLIFT | = constant, Cl, say, sets the lift coefficient to Cl. | |
| CLIFT | = GRND1 selects: Cl = CLIFTA * rc where CLIFTA is a positive constant and rc is the volume fraction of the continuous phase. | |
| CLIFT | = GRND2 selects: Cl = CLIFTA * [1. - 2.78 * min(0.2,rd)] where CLIFTA is a positive constant and rd is the volume fraction of the dispersed phase. |
If CLIFT and CLIFTA are not entered in the argument list of INTSOR, then the INTSOR command sets CLIFT=0.5 and CLIFTA=0.5.
The interphase-lift coding presumes phase 1 is the continuous phase and phase 2 the dispersed phase. However, if CLIFTB > 0.0, phase 2 is taken as the continuous phase and phase 1 as the dispersed phase.
When STORE(LISU,LISV,LISW) appears in the Q1 file, the interfacial lift forces for each cell, as given by equation (1) and integrated over the control volume, may be printed in the RESULT file or viewed via PHOTON and AUTOPLOT. In addition, storage of LISU, LISV and LISW permits linear under-relaxation of the interfacial-lift sources via the RELAX command. This facility can also be activated by extending the INTSOR argument list to include RELAX and RLXFAC as the last two arguments, as follows:
INTSOR(LIFT,CLIFT,CLIFTA,RELAX,RLXFAC)
The argument RELAX activates the facility with a default linear- relaxation factor of 0.1. The argument RLXFAC is optional, and it serves to change the relaxation factor from its default value.
The FORTRAN coding resides Subroutine GXLIFT, which is called from the following groups of GREX3: Group 1, Section 1: Group 13, Section 16; and Group 19, Section 11. The structure and purpose of the various coding sequences is described in the comment lines given at the head of subroutine GXLIFT. The functionality of the model does not extend to BFC=T, although this limitation will be removed in the near future.
Examples of the use of the interphase-lift feature may be found in the advanced multi-phase section of the Input File libraries.
S.W.Beyerlein, R.K.Crossmann and H.J.Richter, 'Prediction of bubble concentration profiles in vertical turbulent two-phase flow', Int.J.Multiphase Flow, Vol.11, No.5, p629, (1985).
A.D.Ciccone, J.G.Kawall and J.F.Keffer, 'On the determination of particle-erosion rates within a turbulent boundary layer', 7th Symposium on Turbulent Shear Flows, 15.4, Stanford University, USA, August, (1989).
D.A.Drew and T.J.Lahey, 'The virtual mass and lift force on a sphere in rotating and straining inviscid flow',Int.J.Multiphase Flow, Vol.13, No.1, p113, (1987).
R.T.Lahey, 'The analysis of phase separation and phase distribution phenomena using two-fluid models', Nuclear Engng. & Design, Vol.122, p17, (1990).
R.T.Lahey, M.Lopez de Bertodano & O.C.Jones, 'Phase distribution in complex geometry ducts', Nuclear Engng. & Design, Vol.141, p177, (1993).
G.J.Lee and S.H.Chang, 'Physical modelling and finite-element method for the analysis of lateral phase distribution phenomena', Int.Comm.Heat Mass Transfer, Vol.18, p333, (1991).
B.Huang, 'Modelisation numerique d'ecoulements diphasiques a bulles dans des reacteurs chimiques', PhD Thesis, L'Universite Claude Bernard - Lyon, (1989).
K.O.Petersen, 'Etude experimentale et numerique des ecoulements diphasiques dans les reacteurs chimiques', PhD Thesis, L'Universite Claude Bernard - Lyon, (1992).
C.Prakash, 'Prediction of some complex multi-dimensional two- phase flow phenomena using the PHOENICS code', Proc. 2nd Int. PHOENICS Users Conference, Heathrow, London, (1987).
A.Serizawa, I.Kataoka and I.Michiyoshi, 'Phase distribution in bubbly flow: Data Set No.24', in Multiphase Science and Technology, Vol.6, p257, Ed. G.F.Hewitt, J.M.Delhaye and N.Zuber, Hemisphere Publishing Corporation, (1992).
H.F.Svensen, H.A.Jakobsen and R.Torvik, 'Local flow structures in internal loop and bubble-column recators', Chem.Eng.Sci., Vol.47, No.13-14, pp3297, (1992).
1. Introduction 2. Activation 3. Exemplification and References
The basic form of the pressure terms in the momentum equations for two-phase flow is as follows ( see Drew [1983] ):
Sk,i = - rk grad(Pk) + ( Pki - Pk ) grad(rk) (1.1)
where: Sk,i is the volumetric source for phase k in direction i; rk is the volume fraction of phase k; Pk is the pressure of phase k in the bulk; and Pki is the pressure of phase k at the interface.
When ONEPHS=F, the default in PHOENICS is to assume that there are no pressure differences between the phases, i.e.
Pk = Pki = P (1.2)
so that equation (1.1) reduces to:
Sk,i = - rk grad(P) (1.3)
This assumption is adequate in applications which do not involve acoustic effects or bubble expansion or contraction.
The interfacial-pressure source terms allow for the possibility of momentum transfer due to pressure discontinuities between the bulk phases and the interface. The option exists for these effects to be included in the PHOENICS momentum equations for two-phase flow via the source terms (1.1) above.
PHOENICS employs two alternative formulations of the interfacial- pressure terms. The first is the one described by Huang [1989], Lahey et al [1993] and Petersen [1992] for dispersed bubbly two- phase flow. The second is the one described by Stuhmiller [1977] and Prakash [1987] for the same applications.
The formulation of Lahey et al [1993] employs the following relationships for the pressure differences:
Pdi - Pd = 0 (1.4)
Pci - Pc = - Cp * rhoc * |Ur.Ur| (1.5)
Pdi - Pci = 2. * sigma * k (1.6)
where: c denotes the continuous phase; d the dispersed phase; Ur is the relative velocity between the two phases; rho is the density; sigma is the surface tension; k is the mean curvature; and Cp is an empirical coefficient defined below.
Equation (1.5) represents a Bernoulli effect in the continuous- phase flow field whereby the pressure on the liquid side of the interface is less than the bulk liquid pressure due to the velocity increase associated with the deflection of the flow around the dispersed bubbles.
If the surface tension is presumed uniform, then when equations (1.4) to (1.6) are substituted into equation (1.1), there results:
Sc,i = - rc grad(Pc) + Cp*rhoc*|Ur.Ur| * grad(rd) (1.7)
and
Sd,i = - rd grad(Pc) + rd grad(Cp*rhoc*|Ur.Ur|) (1.8)
Consequently, the following additional pressure-source terms appear in the momentum equations to describe interfacial-pressure effects:
Sc,i = + Cp*rhoc*|Ur.Ur| * grad(rd) (1.9)
and
Sd,i = + rd grad(Cp*rhoc*|Ur.Ur|) (1.10)
The formulation of Stuhmiller [1987] employs relations (1.5) and (1.6) above, together with:
Pd - Pc = 2. * sigma * k (1.11)
which implies (Pdi-Pd)=(Pci-Pc) rather than (Pdi-Pd)=0 as used in the Lahey formulation.
If the surface tension is presumed uniform, then, when equations (1.5), (1.6) and (1.11) are substituted into equation (1.1), the following interfacial-pressure source terms appear in the momentum equations:
Sc,i = + Cp*rhoc*|Ur.Ur| * grad(rd) (1.12)
and
Sd,i = - Cp*rhoc*|Ur.Ur| * grad(rd) (1.13)
For Cp, Lahey et al [1990] and Prakash [1987] take Cp = 0.25, whereas Huang [1989] and Petersen [1992] take:
Cp = 0.25 * (1.+rd) * rc**2 (1.14)
and Antal et al [1991] take:
Cp = 0.25 * rc (1.15)
where rd and rc are the volume fractions of the continuous and dispersed phases respectively. Although most workers employ Cp=0.25, Lahey et al [1993] suggest larger values for distorted bubbles, such as for example Cp=1.0.
The interfacial-pressure source terms for the momentum equations are activated by the PIL statement:
INTSOR(INTP*,CPIP,CPIPA)
which is equivalent to:
PATCH(INTP*,CELL,1,NX,1,NY,1,NZ,1,LSTEP)
COVAL(INTP*,U1,TINY,GRND4);COVAL(INTP*,U2,TINY,GRND4)
COVAL(INTP*,V1,TINY,GRND4);COVAL(INTP*,V2,TINY,GRND4)
COVAL(INTP*,W1,TINY,GRND4);COVAL(INTP*,W2,TINY,GRND4)
where *=L for the formulation of Lahey et al [1993], and *=S for the formulation of Stuhmiller [1977].
The value ascribed to the PIL variable CPIP determines how the pressure coefficient Cp is to be calculated, as follows:
| CPIP | = 0.0 cuts out the interfacial-pressure terms entirely. | |
| CPIP | = positive constant, Cp, say, sets the pressure coefficient to Cp. | |
| CPIP | = GRND1 selects: Cp = CPIPA * rc | |
| CPIP | = GRND2 selects: Cp = CPIPA * (1.+rd) * rc**2 |
where CPIPA=0.25 for spherical bubbles.
If CPIP and CPIPA are not entered in the argument list of INTSOR, then the INTSOR command sets CPIP=0.25 and CPIPA=0.25.
The interfacial-pressure coding presumes phase 1 is the continuous phase and phase 2 the dispersed phase. However, if CPIP is set to a negative value, phase 2 is taken as the continuous phase and phase 1 as the dispersed phase.
When STORE(IPSU,IPSV,IPSW) appears in the Q1 file, the Stuhmiller interfacial pressure forces for each cell, as given by equation (1.12) and integrated over the control volume, may be printed in the RESULT file or viewed via PHOTON and AUTOPLOT.
In addition, storage of IPSU, IPSV and IPSW permits linear under- relaxation of the Stuhmiller sources via the RELAX command.
If Lahey's formulation is employed then the foregoing facility is available for the dispersed-phase only, in which case the forces are those corresponding to equation (1.10).
The facility can also be activated by extending the INTSOR argument list to include RELAX and RLXFAC as follows:
INTSOR(INTP*,CPIP,CPIPA,RELAX,RLXFAC)
The argument RELAX activates the facility with a default linear- relaxation factor of 0.1. The argument RLXFAC is optional, and it serves to change the relaxation factor from its default value.
The FORTRAN coding resides in Subroutine GXINTP, which is called from the following groups of GREX3: Group 1, Section 1: Group 13, Section 16; and Group 19, Section 4.
The structure and purpose of the various coding sequences is described in the comment lines given at the head of subroutine GXINTP.
For parabolic flows, the interfacial pressure terms in the W1 & W2 momentum equations are neglected.
Examples of the use of the interfacial-pressure feature may be found in the advanced multi-phase section of the Input File libraries.
S.P.Antal, R.T.Lahey and J.E.Flaherty, 'Analysis of phase distribution in fully-developed laminar bubbly two-phase flow', Int. J. Multiphase Flow, Vol.17, No.5, p635, (1991).
D.A.Drew, 'Mathematical modelling of two-phase flow', Ann. Rev. Fluid Mech., Vol.15, p261, (1983).
R.T.Lahey, M.Lopez de Bertodano & O.C.Jones, 'Phase distribution in complex geometry ducts', Nuclear Engng. & Design, Vol.141, p177, (1993).
B.Huang, 'Modelisation numerique d'ecoulements diphasiques a bulles dans des reacteurs chimiques', PhD Thesis, L'Universite Claude Bernard - Lyon, (1989).
K.O.Petersen, 'Etude experimentale et numerique des ecoulements diphasiques dans les reacteurs chimiques', PhD Thesis, L'Universite Claude Bernard - Lyon, (1992).
C.Prakash, 'Prediction of some complex multi-dimensional two- phase flow phenomena using the PHOENICS code', Proc. 2nd Int. PHOENICS Users Conference, Heathrow, London, (1987).
J.H.Stuhmiller, 'The influence of interfacial pressure forces on the character of two-phase flow model equations', Int. J. Multiphase Flow, Vol.3, p551, (1977).
wbs