AUTOPLOT USE
  file
  phi 5
 
 
  da 1 u1; da 1 uliq
  col3 1;blb4 2;scale y 0 .015
  redr
  msg 1st-phase velocity
  msg blue PHOENICS   + analytical
  msg press  to continue
  pause
  cl
  screen
  da 1 r1; da 1 rliq
  col3 1;blb4 2;scale y 0 1.
  redr
  msg 1st-phase volume fraction
  msg blue PHOENICS   + analytical
  msg press  to continue
  pause
  cl
  msg press e to END
  enduse
 
TEXT(1D SOLIDS TRANSPORT IN WATER
TITLE
  DISPLAY
  The problem considered is the one-dimensional vertical flow of
  spherical particles by water in a duct. The inlet void fraction
  is 0.4 and the ratio of the solids mass inflow rate to that of
  the gas is 0.5. The duct area is 1 m**2 and the height is 1m.
  The particles are 100 microns in diameter and the density ratio
  is 2.6. The inlet liquid velocity is taken as twice the particle
  terminal velocity. The task is to predict the volume fractions
  and phase velocities throughout the duct. The calculations employ
  the built-in particle-fluidisation drag model (CFIPD=7.), which
  uses the Ergun correlation for void fractions less than 0.8
  (dense fluidisation), and a modified spherical-particle drag
  correlation for void fraction greater than 0.8 (dilute
  fluidisation). The test case provides a test of both dense and
  dilute suspensions according to the setting of DENSE. In both
  cases the PHOENICS predictions show excellent agreement with the
  approximate analytical solutions.
  ENDDIS
 
  * CONPHS=1 selects 1st phase as the continuous ( gas ) phase
             and thus CFIPS=GRND7
  *       =2 selects 2nd phase as the continuous ( gas ) phase
             and thus CFIPS=GRND8
  * VARLAM=T tests coding in GXCFIP when ENUL is variable
  * CFIPD=7. selects the particle-fluidisation drag correlation
  * DENSE=T selects test of Ergun interphase drag model
         =F   ''     ''  of Dilute particle drag model
BOOLEAN(VARLAM,DENSE);INTEGER(CONPHS);VARLAM=T;CONPHS=1
DENSE=T
  * for 1dx calculation set: CH1=X;CH2=U1;CH3=U2
  * for 1dy calculation set: CH1=Y;CH2=V1;CH3=V2
  * for 1dz calculation set: CH1=Z;CH2=W1;CH3=W2
CHAR(CH1,CH2,CH3);CH1=X;CH2=U1;CH3=U2
   CH1=Y;CH2=V1;CH3=V2
   CH1=Z;CH2=W1;CH3=W2
REAL(MRATIO,XLEN,UINL,UINP,UIN1,UIN2,R1IN,R2IN,REYP,DIAMP,DENRAT)
REAL(UTERM,RINL,RINP,FLOW1,FLOW2,dens,RHOL,UMF,DELRHO)
XLEN=1.0;RHOL=1000.0;dens=2600.0;DENRAT=dens/RHOL
IF(DENSE) THEN
+ RINL=0.4;MRATIO=0.5
ELSE
+ RINL=0.6;MRATIO=0.2
ENDIF
DELRHO=dens-RHOL;DIAMP=100.E-6;ENULA=1.0E-3/RHOL;RINP=1.0-RINL
  ** compute minimum fluidisation velocity & Re
UMF=DIAMP**2*DELRHO*9.81*RINL**3/(180.*ENULA*RHOL*RINP)
umf
  ** compute terminal velocity of particles
REYP=UMF*DIAMP/ENULA
UTERM=9.81*DELRHO*DIAMP**2/(18.*RHOL*ENULA)
uterm
REYP=UTERM*DIAMP/ENULA
UINL=2.*UTERM;UINP=MRATIO*(RHOL/dens)*(RINL/RINP)*UINL
uinl
uinp
IF(CONPHS.EQ.1) THEN
+ UIN1=UINL;UIN2=UINP;R2IN=RINP;R1IN=1-R2IN
ELSE
+ UIN1=UINP;UIN2=UINL;R1IN=RINP;R2IN=1-R1IN
ENDIF
  ** compute analytical liquid volume fraction & velocity
REAL(AA,BB,CC,QL,QP,ROUTL,DELR,UOUTL,UOUTP)
QL=UINL*RINL;QP=UINP*RINP
IF(DENSE) THEN
  ** assume viscous term of Ergun correlation dominates
AA=150.*RHOL*ENULA/(DELRHO*9.81*DIAMP*DIAMP)
BB=AA*QL;AA=AA*(QL+QP);ROUTL=RINL
DO JJ=1,6
+ DELR=(ROUTL**3+AA*ROUTL-BB)/(3.*ROUTL**2+BB)
+ ROUTL = ROUTL - DELR
ENDDO
ELSE
  ** assume Stokes flow
ROUTL=0.9;AA=DELRHO*9.81*DIAMP*DIAMP/(18.*RHOL*ENULA)
DO JJ=1,6
+ BB = (1.-ROUTL)*(AA*ROUTL**4.65-QL)+QP*ROUTL
+ CC=4.65*AA*ROUTL**3.65-5.65*AA*ROUTL**4.65+(QL+QP)
+ DELR=BB/CC
+ ROUTL = ROUTL - DELR
ENDDO
ENDIF
routl
UOUTL=QL/ROUTL;UOUTP=QP/(1.-ROUTL)
uoutl
uoutp
    GROUP 1. Run title and other preliminaries
    GROUP 2. Transience; time-step specification
    GROUP 3. X-direction grid specification
GRDPWR(:CH1:,20,XLEN,1.0)
    GROUP 4. Y-direction grid specification
    GROUP 5. Z-direction grid specification
    GROUP 6. Body-fitted coordinates or grid distortion
    GROUP 7. Variables stored, solved & named
ONEPHS=F;SOLVE(P1,:CH2:,:CH3:,R1,R2)
  Activate storage for printout of interphase drag properties
STORE(REYN,VREL,CFIP,SIZE,ULIQ,RLIQ)
    GROUP 8. Terms (in differential equations) & devices
    GROUP 9. Properties of the medium (or media)
IF(CONPHS.EQ.1) THEN
+ RHO1=RHOL;RHO2=dens
ELSE
+ RHO2=dens;RHO1=RHOL
ENDIF
IF(VARLAM) THEN
    Use ENUL=ENULA+ENULB*TMP1
+ STORE(TMP1);ENUL=LINTEM;ENULA=1.E-6;ENULB=0.0;TMP1=CONST
ELSE
+ ENUL=ENULA
ENDIF
    GROUP 10. Inter-phase-transfer processes and properties
IF(CONPHS.EQ.1) THEN
+ CFIPS=GRND7;VARMIN(R2)=1.E-10
ELSE
+ CFIPS=GRND8;VARMIN(R1)=1.E-10
ENDIF
  ** CFIPA = minimum slip velocity  CFIPB = bubble size
CFIPA=1.E-4;CFIPB=DIAMP;CFIPD=7.
    GROUP 11. Initialization of variable or porosity fields
FIINIT(U1)=UIN1;FIINIT(U2)=UIN2;FIINIT(R2)=R1IN;FIINIT(R1)=R2IN
FIINIT(RLIQ)=ROUTL;FIINIT(ULIQ)=UOUTL
    GROUP 12. Unused
    GROUP 13. Boundary conditions and special sources
FLOW1=RHO1*UIN1*R1IN;FLOW2=RHO2*UIN2*R2IN
INLET(IN,CELL,$1,$1,$1,$1,$1,$1,1,1)
VALUE(IN,P1,FLOW1);VALUE(IN,:CH2:,UIN1)
VALUE(IN,P2,FLOW2);VALUE(IN,:CH3:,UIN2)
OUTLET(OUT,CELL,%1,%1,%1,%1,%1,%1,1,1)
  ** use expected outflow velocities to allow for slip
     in the exit-cell outflow
IF(CONPHS.EQ.1) THEN
+ COVAL(OUT,P1,RHO1*UOUTL*1.E2,0)
+ COVAL(OUT,P2,RHO2*UOUTP*1.E2,0)
ELSE
+ COVAL(OUT,P1,RHO1*UOUTP*1.E2,0)
+ COVAL(OUT,P2,RHO2*UOUTL*1.E2,0)
ENDIF
PATCH(GRAVITY,PHASEM,1,NX,1,NY,1,NZ,1,1)
IF(CONPHS.EQ.1) THEN
+ COVAL(GRAVITY,:CH3:,FIXFLU,-9.81*(1.-RHOL/dens))
ELSE
+ COVAL(GRAVITY,:CH2:,FIXFLU,-9.81*(1.-RHOL/dens))
ENDIF
    GROUP 15. Termination of sweeps
LSWEEP=250
RESREF(P1)=1.E-12*(FLOW1+FLOW2)
RESREF(R1)=1.E-12*FLOW1;RESREF(R2)=1.E-12*FLOW2
RESREF(:CH2:)=RESREF(R1)*UIN1;RESREF(:CH3:)=RESREF(R2)*UIN2
    GROUP 16. Termination of iterations
    GROUP 17. Under-relaxation devices
REAL(DTF);DTF=0.1*XLEN/UINL/N:CH1:
RELAX(R1,LINRLX,0.3);RELAX(R2,LINRLX,0.3)
RELAX(:CH2:,FALSDT,DTF);RELAX(:CH3:,FALSDT,DTF)
RELAX(CFIP,LINRLX,0.3)
    GROUP 18. Limits on variables or increments to them
    GROUP 19. Data communicated by satellite to GROUND
    GROUP 20. Preliminary print-out
    GROUP 21. Print-out of variables
OUTPUT(CFIP,P,P,P,P,Y,P);NXPRIN=1;NYPRIN=1;NZPRIN=1
    GROUP 22. Spot-value print-out
IF(:CH1:.EQ.X) THEN
+ IXMON=NX/2;IZMON=1;IYMON=1
ELSE
+ IZMON=NZ/2;IXMON=1;IYMON=1
ENDIF
IF(:CH1:.EQ.Y) THEN
+ IYMON=NY/2;IXMON=1;IZMON=1
ENDIF
TSTSWP=-1
    GROUP 23. Field print-out and plot control
    GROUP 24. Dumps for restarts