GROUP 1. Run title and other preliminaries
TEXT(LAM-BRE K-E BKWRD-FACING STEP Y-X  :T208
TITLE
  DISPLAY
   The problem concerns 2d incompressible, turbulent flow and heat
   transfer over a backward-facing step. The case is similar to that
   described for library case T103, excepting that heat transfer is
   considered and the calculation is performed with any one of the
   following low-Re models: the Lam-Bremhorst k-e model, the Chen-
   Kim k-e model, the 2-layer k-e model, and the Wilcox k-omega
   model. The thermal field is calculated for a constant heat flux
   through the south wall downstream of the step and an adiabatic
   condition on all other walls. The laminar Prandtl number is taken
   as 1.0 and the turbulent Prandtl as 0.86.
 
  ENDDIS
 
  For this case, both the low- and high-Re forms of the k-e model
  are known to underpredict the reattachment length XR by about 14%,
  yielding a value of XR/H=6.1 rather XR/H=7.1. The fine-mesh
  calculations employ a relatively coarse non-uniform mesh of NX=80
  by NY=80 which concentrate grid cells close to the wall. About
  1000 sweeps are required for a converged solution with this mesh.
  However, for speed of computation the calculation is performed on
  a coarse mesh of NX=25 by NY=20. Consequently, the solution is
  unlikely to be of sufficient quality because most of the near-wall
  cells are located too far from the wall for laminar conditions to
  apply on KE and EP. PHOENICS will actually activate log-law wall
  functions for velocities if the near-wall point lies outside the
  viscous sublayer. For demonstration purposes the calculation is
  set up for 30 sweeps, but 200 sweeps are required for complete
  convergence.
 
    This AUTOPLOT sequence provides a plot of the axial
    velocity U1 along the bottom surface of the solution
    domain versus normalised axial distance X. The axial
    coordinate 0.0 corresponds to the step location. The
    reattachment point corresponds the axial location where
    U1 changes from negative to positive.
 
   AUTOPLOT USE
   file
   phi 5
 
   da 1 u1 y 1
   divide x .0381 1
   shift x -4 1
   col9 1
   level y 0;level x 0
   msg Velocity (U1) profile
   msg Press e to END
   ENDUSE
 
CHAR(CTURB,TLSC)
REAL(FRIC);INTEGER(TMODEL,NYS,NXS,NXS2);BOOLEAN(FINEG);FINEG=T
REAL(HGHT,WIDTH,CLEN,SLEN,SLEN2,REYNO,UIN,TKEIN,EPSIN,QIN,AREAO)
     ** Calculation of domain specifications
REYNO=4.5E4;UIN=13.;HGHT=0.0381;WIDTH=3.*HGHT
SLEN=4.*HGHT;SLEN2=5.5*HGHT;CLEN=20.*HGHT;FRIC=0.018
TKEIN=0.25*UIN*UIN*FRIC;EPSIN=0.1643*TKEIN**1.5/(0.09*HGHT)
    GROUP 2. Transience; time-step specification
    GROUP 3. X-direction grid specification
    GROUP 4. Y-direction grid specification
  ** channel length=0.762m & channel width=0.1143m
IF(FINEG) THEN
+ NXS=14;NXS2=61;NYS=40
+ SUBGRD(X,1,NXS,SLEN,-2.0);SUBGRD(X,NX+1,NXS2,SLEN2,2.0)
+ SUBGRD(X,NX+1,80,CLEN-SLEN-SLEN2,1.2)
+ SUBGRD(Y,1,-NYS,HGHT,1.7);SUBGRD(Y,NY+1,-80,2.0*HGHT,2.0)
ELSE
+ NXS=7;NXS2=17;NYS=10
+ SUBGRD(X,1,NXS,SLEN,-2.0);SUBGRD(X,NX+1,NXS2,SLEN2,2.0)
+ SUBGRD(X,NX+1,25,CLEN-SLEN-SLEN2,1.2)
+ SUBGRD(Y,1,-NYS,HGHT,1.7);SUBGRD(Y,NY+1,-20,2.0*HGHT,2.0)
ENDIF
    GROUP 5. Z-direction grid specification
    GROUP 6. Body-fitted coordinates or grid distortion
    GROUP 7. Variables stored, solved & named
SOLVE(P1,U1,V1,H1);SOLUTN(P1,Y,Y,Y,N,N,N);STORE(ENUT,LEN1)
SOLUTN(U1,P,P,P,P,P,N);SOLUTN(V1,P,P,P,P,P,N);STORE(FMU,FTWO)
STORE(PRPS)
MESG( Enter the required turbulence model:
MESG(  LAMB - Low-Re Lam-Brem. k-e model (default)
MESG(  CHEN - Low-Re Lam-Brem. k-e model
MESG(  TWOL - Low-Re 2-layer   k-e model
MESG(  KO   - Low-Re k-omega model
MESG(
READVDU(CTURB,CHAR,LAMB)
CASE :CTURB: OF
WHEN LAMB,4
+ MESG(Low-Re Lam-Bem. k-e turbulence model
+ TMODEL=1;KELIN=1;TLSC=EP
+ TURMOD(KEMODL-LOWRE);STORE(REYN)
WHEN CHEN,4
+ TEXT(CHEN-KIM K-E BKWRD-FACING STEP Y-X  :T208
+ MESG(Low-Re Chen-Kim k-e turbulence model
+ TMODEL=2;KELIN=1;TLSC=EP
+ TURMOD(KECHEN-LOWRE);STORE(REYN)
WHEN TWOL,4
+ TEXT(2-LAYER K-E BKWRD-FACING STEP Y-X  :T208
+ MESG(Low-Re 2-layer k-e turbulence model
+ TMODEL=3;KELIN=1;TLSC=EP
+ TURMOD(KEMODL-2L);STORE(REYN)
WHEN KO,2
+ TEXT(K-OMEGA BKWRD-FACING STEP Y-X  :T208
+ MESG(Low-Re k-omega turbulence model
+ TMODEL=4;TLSC=OMEG
+ TURMOD(KOMODL-LOWRE);STORE(REYT)
+ STORE(EP);EPSIN=EPSIN/(0.09*TKEIN)
ENDCASE
    GROUP 8. Terms (in differential equations) & devices
TERMS(H1,N,P,P,P,P,P)
    GROUP 9. Properties of the medium (or media)
RHO1=1.0;ENUL=UIN*HGHT/REYNO;PRT(H1)=0.86
    GROUP 10. Inter-phase-transfer processes and properties
    GROUP 11. Initialization of variable or porosity fields
CONPOR(0.0,VOLUME,-1,-NXS,-1,-NYS,1,1)
FIINIT(U1)=0.1*UIN;FIINIT(V1)=0.0;FIINIT(P1)=1.3E-4
FIINIT(:TLSC:)=EPSIN;FIINIT(KE)=0.04*UIN*UIN;FIINIT(H1)=0.5
    GROUP 13. Boundary conditions and special sources
INLET(INLET,WEST,1,1,NYS+1,NY,1,1,1,1);VALUE(INLET,P1,RHO1*UIN)
VALUE(INLET,U1,UIN);VALUE(INLET,KE,TKEIN);VALUE(INLET,:TLSC:,EPSIN)
PATCH(OUTLET,EAST,NX,NX,1,NY,1,1,1,1);COVAL(OUTLET,P1,1.E5,0.0)
COVAL(OUTLET,H1,ONLYMS,SAME)
WALL(WALLN,NORTH,1,NX,NY,NY,1,1,1,1)
WALL(WALLS,SOUTH,NXS+1,NX,1,1,1,1,1,1)
   ** Heat input; take cp=1 & set Qin for unit temperature rise
QIN=RHO1*2.*HGHT*UIN;AREAO=ZWLAST*(CLEN-SLEN)
PATCH(HEATIN,SOUTH,NXS+1,NX,1,1,1,NZ,1,1)
COVAL(HEATIN,H1,FIXFLU,QIN/AREAO)
EGWF=F
    GROUP 14. Downstream pressure for PARAB=.TRUE.
    GROUP 15. Termination of sweeps
LSWEEP=1000;TSTSWP=-1
    GROUP 16. Termination of iterations
LITER(U1)=2;LITER(V1)=2;LITER(KE)=5;LITER(:TLSC:)=5;LITER(H1)=10
    GROUP 17. Under-relaxation devices
REAL(DTF);DTF=10.*CLEN/UIN/NX; RELAX(P1,LINRLX,0.5)
CASE (TMODEL) OF
WHEN 1
  ** lam-bremhorst k-e model
+ KELIN=3; RELAX(KE,LINRLX,0.3); RELAX(EP,LINRLX,0.5)
 EX(P1  )=1.976E+01; EX(U1  )=6.848E+00; EX(V1  )=2.304E-01
 EX(KE  )=1.992E+00; EX(EP  )=4.018E+02; EX(H1  )=4.132E+00
 EX(PRPS)=9.125E-01; EX(EPKE)=2.510E+02; EX(VPOR)=9.125E-01
 EX(REYN)=2.129E+03; EX(LTLS)=6.673E-04; EX(WDIS)=1.796E-02
 EX(FTWO)=9.083E-01; EX(FMU )=8.250E-01; EX(LEN1)=3.382E-03
 EX(ENUT)=2.477E-03
WHEN 2
  ** chen-kim k-e model
+ KELIN=3; RELAX(KE,LINRLX,0.3); RELAX(EP,LINRLX,0.5)
 EX(P1  )=2.161E+01; EX(U1  )=6.836E+00; EX(V1  )=1.881E-01
 EX(KE  )=1.543E+00; EX(EP  )=3.550E+02; EX(H1  )=5.842E+00
 EX(PRPS)=9.125E-01; EX(EPKE)=1.349E+03; EX(VPOR)=9.125E-01
 EX(REYN)=1.926E+03; EX(LTLS)=6.673E-04; EX(WDIS)=1.796E-02
 EX(FTWO)=9.079E-01; EX(FMU )=8.176E-01; EX(LEN1)=3.057E-03
 EX(ENUT)=1.869E-03
WHEN 3
  ** two-layer k-e model
  ** this model does not work on V2.2 when SELREF=T
+ SELREF=F;KELIN=3; RELAX(KE,LINRLX,0.3); RELAX(EP,LINRLX,0.5)
 EX(P1  )=2.133E+01; EX(U1  )=6.940E+00; EX(V1  )=2.420E-01
 EX(KE  )=1.780E+00; EX(EP  )=3.938E+02; EX(H1  )=5.698E+00
 EX(PRPS)=9.125E-01; EX(EPKE)=5.739E+02; EX(VPOR)=9.125E-01
 EX(REYN)=2.106E+03; EX(LTLS)=6.673E-04; EX(WDIS)=1.796E-02
 EX(FTWO)=4.452E+00; EX(FMU )=8.314E-01; EX(LEN1)=2.855E-03
 EX(ENUT)=2.038E-03
WHEN 4
  ** k-omega model
+ DTF=DTF/10
+ RELAX(KE,FALSDT,DTF); RELAX(:TLSC:,FALSDT,DTF)
 EX(P1  )=2.246E+01; EX(U1  )=6.814E+00; EX(V1  )=1.981E-01
 EX(KE  )=1.554E+00; EX(EP  )=6.231E+02; EX(H1  )=6.998E+00
 EX(VPOR)=9.125E-01; EX(REYT)=1.699E+02; EX(OMEG)=5.202E+03
 EX(PRPS)=9.125E-01; EX(FTWO)=8.520E-01; EX(FMU )=7.946E-01
 EX(LEN1)=2.840E-03; EX(ENUT)=1.819E-03
ENDCASE
RELAX(U1,FALSDT,DTF); RELAX(V1,FALSDT,DTF); RELAX(H1,FALSDT,DTF)
    GROUP 18. Limits on variables or increments to them
    GROUP 21. Print-out of variables
    GROUP 22. Monitor print-out
ITABL=3;IXMON=NXS+2;IPLTL=2000;IYMON=NYS-2;NPRMON=10000
WALPRN=T
 LIBREF  =       208
DISTIL=T
STOP