TEXT(2-PHASE TURBULENT AIR-SOLIDS PIPE FLOW
TITLE
  DISPLAY
    The case considered is 2-phase solids-laden turbulent vertical
    flow of air in a pipe, as studied experimentally by Tsuji et
    al [J.Fluid Mech, Vol.139, p417, 1984]. 
    The calculation is performed with the parabolic option.
    For comparison with data the calculation should be carried out
    until the flow is fully developed, i.e. 70 diameters downstream.
    The calculation may be made with the standard k-e model, or
    alternatively with the Chen-Wood or Mostafa-Mongia k-e variants
    which allow for gas turbulence modulation due to the presence 
    of particles. 
    The pipe Reynolds number is 3E4, the density ratio
    is 866, the particle mass-flow loading may be 1.0 or 2.1, and
    the particle diameter is 200 microns. The task is to calculate
    the fully-developed vertical velocity profiles for comparison
    with the measured profiles.
    For the results of a study made with a 1995 version of 
    PHOENICS, click
  here.
  ENDDIS
  PHOTON USE
  p                                                                               
  parphi  
  20 1                                                                        
                                                                                
  msg flow is from right to left; wall is at top, axis at bottom
  msg radial dimension enlarged 20-fold                                                                              
  msg axial velocity contours for phase 1                                                                         
  con w1 x 1 fi ;.001
  pause;con off; red                                                             
  msg axial velocity contours for phase 2                                                                         
  con w2 x 1 fi;0.001                                                             
  pause;con off;red                                                             
  msg volume-fraction contours for phase 2                                                                         
  con r2 x 1 fi;0.001                                                             
  pause;con off;red                                                             
  msg interphase-friction factors                                                                         
  con cfip x 1 fi;0.001                                                             
  pause; con off; red                                                             
  msg effective (i.e. turbulent) viscosity                                                                       
  con enut x 1 fi;0.001                                                             
  pause; con off; red                                                             
  msg axial velocity contours for phase 1                                                                       
  vec x 1 sh                                                                        
  enduse

CHAR(CTURB);INTEGER(JRUN,NSTEP)
MESG( Enter the required turbulence-modulation option;
MESG( default NONE
MESG( The options are:
MESG(  MOST  - Mostafa-Mongia k-e modulation sources
MESG(  CHEN  - Chen-Wood      k-e modulation sources
MESG(  NONE  - Standard k-e model
MESG(
READVDU(CTURB,CHAR,NONE)
REAL(RHOS,RHOG,EMUGAS,RGAS,RSOL,TKEIN,EPSIN,FLOWG,FLOWS)
REAL(VGAS,VSOL,DIAM,GRAD,PLEN,DIAMS,REY,GRAVAC,MLOAD)
RHOS=1020.;GRAVAC=-9.81;DIAM=0.0305;RHOG=1.178;GRAD=0.5*DIAM
MESG( Enter JRUN: default =1
MESG(             =1 solids-to-air mass-flow ratio=1.0
MESG(             =2 solids-to-air mass-flow ratio=2.1
READVDU(JRUN,INT,1)
IF(JRUN.EQ.1) THEN
  ** forward step-size is 0.175*diameter
+ MLOAD=1.0;REY=3.E4;PLEN=70.*DIAM;NSTEP=400
ELSE
  ** forward step-size is 0.2*diameter
+ MLOAD=2.1;REY=3.3E4;PLEN=100.*DIAM;NSTEP=500
ENDIF
plen
plen=plen*4
RSOL=MLOAD*RHOG/RHOS;RGAS=1.0-RSOL;EMUGAS=1.868E-5
VGAS=REY*EMUGAS/(RHOG*DIAM);VSOL=VGAS
   *** specify mean particle diameter
DIAMS=200.E-6
    GROUP 1. Run title and other preliminaries
    GROUP 2. Transience; time-step specification
PARAB=T
    GROUP 3. X-direction grid specification
CARTES=F;XULAST=0.1;NX=1
    GROUP 4. Y-direction grid specification
GRDPWR(Y,25,GRAD,1.0)
    GROUP 5. Z-direction grid specification
NZ=10
nz=200
   NZ=NSTEP
GRDPWR(Z,NZ,NZ*PLEN/NSTEP,1.0)
GRDPWR(Z,NZ,PLEN,1.0)
    GROUP 6. Body-fitted coordinates or grid distortion
    GROUP 7. Variables stored, solved & named
ONEPHS=F;SOLVE(P1,V1,V2,W1,W2,R1,R2);TURMOD(KEMODL)
STORE(CFIP,ENUT,CD,REYN)
    GROUP 8. Terms (in differential equations) & devices
    GROUP 9. Properties of the medium (or media)
RHO1=RHOG;RHO2=RHOS
FLOWG=RHOG*RGAS*VGAS;FLOWS=RHOS*RSOL*VSOL;ENUL=EMUGAS/RHOG
TKEIN=(0.05*VGAS)**2;EPSIN=0.1643*TKEIN**1.5/(0.1*GRAD)
PRNDTL(V2)=1.E10;PRNDTL(W2)=1.E10;PRT(V2)=1.E10;PRT(W2)=1.E10
    GROUP 10. Inter-phase-transfer processes and properties
CFIPS=GRND7;CFIPA=1.E-3;RLOLIM=1.E-6;CFIPB=DIAMS
    GROUP 11. Initialization of variable or porosity fields
FIINIT(W1)=VGAS;FIINIT(W2)=VSOL;FIINIT(R1)=RGAS;FIINIT(R2)=RSOL
    GROUP 12. Unused
    GROUP 13. Boundary conditions and special sources
  ** inlet boundary
INLET(IN,LOW,1,NX,1,NY,1,1,1,1)
VALUE(IN,P1,FLOWG);VALUE(IN,W1,VGAS)
VALUE(IN,P2,FLOWS);VALUE(IN,W2,VSOL)
VALUE(IN,KE,TKEIN);VALUE(IN,EP,EPSIN)
  ** gravity
PATCH(GRAVITY,PHASEM,1,NX,1,NY,1,NZ,1,1)
COVAL(GRAVITY,W1,FIXFLU,GRAVAC);COVAL(GRAVITY,W2,FIXFLU,GRAVAC)
  ** wall boundary
WALL(NWALL,NORTH,1,NX,NY,NY,1,NZ,1,1)
  ** turbulence-modulation sources
CASE :CTURB: OF
WHEN NONE,4
+ MESG(No turbulence modulation sources
+ TEXT(AIR-SOLIDS FLOW- No turb. Mods
WHEN MOST,4
+ MESG(Mostafa-Mongia turbulence-modulation sources
+ PATCH(KEDISP,CELL,1,NX,1,NY,1,NZ,1,LSTEP);STORE(TE,TP)
+ COVAL(KEDISP,KE,GRND1,ZERO);COVAL(KEDISP,EP,GRND1,ZERO)
+ TEXT(AIR-SOLIDS FLOW- Mostafa-Mongia turb. mods
WHEN CHEN,4
+ TEXT(AIR-SOLIDS FLOW- Chen-Wood turb. mods
+ MESG(Chen-Wood turbulence-modulation sources
+ PATCH(KEDISP,CELL,1,NX,1,NY,1,NZ,1,LSTEP);STORE(TE,TP)
+ COVAL(KEDISP,KE,GRND2,ZERO);COVAL(KEDISP,EP,GRND2,ZERO)
ENDCASE
    GROUP 15. Termination of sweeps
LITHYD=10;SELREF=T;RESFAC=0.01
    GROUP 16. Termination of iterations
    GROUP 17. Under-relaxation devices
RELAX(V1,LINRLX,0.4);RELAX(V2,LINRLX,0.4)
RELAX(W1,LINRLX,0.5);RELAX(W2,LINRLX,0.5)
RELAX(R1,LINRLX,0.3);RELAX(R2,LINRLX,0.3)
RELAX(KE,LINRLX,0.3);RELAX(EP,LINRLX,0.3)
    GROUP 18. Limits on variables or increments to them
VARMIN(W1)=1.E-10;VARMIN(W2)=1.E-10
VARMIN(R1)=1.E-10;VARMIN(R2)=1.E-10
    GROUP 19. Data communicated by satellite to GROUND
    GROUP 20. Preliminary print-out
    GROUP 21. Print-out of variables
  NYPRIN=1
NZPRIN=NZ/5
    GROUP 22. Spot-value print-out
IYMON=NY-2;TSTSWP=LITHYD;NPLT=4;IDISPA=NSTEP/10
    GROUP 23. Field print-out and plot control
patch(endprof,profil,1,1,1,ny,1,nz,1,1)    
coval(endprof,w1,0,0)    
coval(endprof,w2,0,0) 
idispa=10
tstswp=-10
#maxabs  
    GROUP 24. Dumps for restarts