GROUP 1. Run title and other preliminaries
TEXT(REALISABLE KE_2D FLOW PAST A SQUARE RIB: T307
TITLE
  DISPLAY
  The case considered is 2D, steady, incompressible, turbulent flow
  past a surface-mounted square rib in a channel. This case has
  been studied experimentally by D.Crabb et al, Proc. 4th Brazilian
  Congress on Mech. Engng., Florianopolis, Brazil, p415,(1997).
 
  The height H of the rib is 8.5% of that of the channel. The flow
  Reynolds number based on channel bulk velocity and rib height
  H is 300,000. The inlet plane is located 6H upstream of the rib,
  and the outlet plane 20H downstream of the rib. A fixed-pressure
  boundary condition is applied at the outlet, and uniform flow
  profiles are specified at the inlet.
 
  The case is set up to run one of six versions of the k-e model
  with scalable wall functions, namely the standard model, the MMK
  model, the Kato-Launder model, the Chen-Kim & RNG models, or
  the realisable k-e model. An option is provided to also run the 
  standard k-w model.
 
  For this case, the main parameter characterising separation is
  the length of the separation zone behind the rib. The experimental
  and computed results are:
 
             K-E  K-O  KL    MMK   RKE   Chen  RNG   EXPT
 
     Lr/H =  5.2  6.3 11.4   11.7  12.8  14.1  13.0  12.3

  where the separation length Lr is measured from the front of the
  rib. These results are not grid independent, and the mesh is not 
  fine enough to resolve the separation regions around the rib. 
  The standard k-e & k-w models seriously underpredict the length
  of the separation behind the rib, and the realisable k-e model
  gives closest agreement with the data. 
  ENDDIS

    The AUTOPLOT sequence below provides a plot of the axial
    velocity W1 along the bottom surface of the solution domain 
	versus normalised axial distance. The axial coordinate 0.0 
	corresponds to the rear surface of the rib. The reattachment 
	point behind the cube corresponds the axial location where W1 
	changes from negative to positive.
   AUTOPLOT USE
    AUTOPLOT
    FILE                                                                            
    PHIDA 3                                                                         
                                                                                                                                                                                                                                    
    D 1 W1 Y 1                                                                      
    PLOT                                                                            
    LEVEL Y 0                                                                       
    SHIFT X -7 1                                                                    
    REDR                                                                            
    SCALE X 0 15    
    msg Press e to END
  ENDUSE
CHAR(CTURB)
REAL(HCUBE,CLUP,CLDOWN,CHIGHT,CWIDTH)
REAL(REYNO,UIN,TKEIN,EPSIN,MIXL,FRIC,OMIN)
INTEGER(NYC,NZC,NZUP,NZDOWN,NYUP,JKO)
     ** Calculation of domain specifications
HCUBE=1.0;UIN=1.0
CHIGHT=11.75*HCUBE
CLUP=6.*HCUBE;CLDOWN=20.*HCUBE;REYNO=3.E5
 
FRIC=0.014;TKEIN=0.25*UIN*UIN*FRIC
MIXL=0.09*CHIGHT;EPSIN=0.1643*TKEIN**1.5/MIXL
 
NYC=20;NYUP=44
NZUP=34;NZC=12;NZDOWN=64
    GROUP 3. X-direction grid specification
    GROUP 4. Y-direction grid specification
NREGY=2
IREGY=1;GRDPWR(Y,NYC,HCUBE,1.0)
IREGY=2;GRDPWR(Y,NYUP,-(CHIGHT-HCUBE),1.07)
    GROUP 5. Z-direction grid specification
NREGZ=3
IREGZ=1;GRDPWR(Z,NZUP,-CLUP,-1.05)
IREGZ=2;GRDPWR(Z,NZC,HCUBE,1.0)
IREGZ=3;GRDPWR(Z,NZDOWN,-CLDOWN,1.04)
    GROUP 7. Variables stored, solved & named
SOLVE(P1,V1,W1);SOLUTN(P1,Y,Y,Y,N,N,N);STORE(ENUT)
SOLUTN(V1,P,P,P,P,P,N);SOLUTN(W1,P,P,P,P,P,N)
JKO=0
MESG( Enter the required turbulence model:
MESG(  KEM  -  Standard k-e model  
MESG(  MMK  -  MMK  k-e model      
MESG(  KLM  -  KL   k-e model      
MESG(  CKM  -  Chen-Kim  k-e model
MESG(  RNG  -  RNG k-e model
MESG(  KOM  -  k-w model
MESG(  RKE  -  realisable k-e model (default)
MESG(
READVDU(CTURB,CHAR,RKE)
CASE :CTURB: OF
WHEN KEM,3
+ TEXT(K-E 2D SQUARE RIB FLOW  :T307
+ MESG(Standard k-e model
+ TURMOD(KEMODL)
WHEN MMK,3
+ TEXT(MMK K-E SQUARE RIB FLOW :T307
+ MESG(MMK k-e model
+ TURMOD(KEMMK)
WHEN KLM,3
+ TEXT(KLM K-E SQUARE RIB FLOW :T307
+ MESG(KLM k-e model
+ TURMOD(KEKL)
WHEN CKM,3
+ TEXT(Chen-Kim K-E SQUARE RIB FLOW :T307
+ MESG(Chen Kim k-e model
+ TURMOD(KECHEN)
WHEN RNG,3
+ TEXT(RNG K-E SQUARE RIB FLOW :T307
+ MESG(RNG k-e model
+ TURMOD(KERNG)
WHEN KOM,3
+ TEXT(K-W SQUARE RIB FLOW :T307
+ MESG(k-w model
+ TURMOD(KOMODL);STORE(EP)
+ JKO=1;OMIN=EPSIN/(0.09*TKEIN)
WHEN RKE,3
+ TEXT(RK K-E SQUARE RIB FLOW :T307
+ MESG(RK k-e model
+ TURMOD(KEREAL);STORE(C1E)
ENDCASE
STORE(YPLS)
SCALWF=T  ! Scalable wall functions 
    GROUP 8. Terms (in differential equations) & devices
    GROUP 9. Properties of the medium (or media)
RHO1=1.0;ENUL=UIN*HCUBE/REYNO
    GROUP 11. Initialization of variable or porosity fields
FIINIT(W1)=UIN;FIINIT(P1)=1.3E-4
FIINIT(KE)=TKEIN;FIINIT(EP)=EPSIN;FIINIT(V1)=0.001*UIN
IF(JKO.EQ.1) THEN
+ FIINIT(OMEG)=OMIN
ENDIF
     ** Initialization of variables in blocked region
CONPOR(RIB,0.0,CELL,-#1,-#1,-#1,-#1,-#2,-#2)
    GROUP 13. Boundary conditions and special sources
INLET(INLET,LOW,#1,#NREGX,#1,#NREGY,#1,#1,1,1)
VALUE(INLET,P1,UIN);VALUE(INLET,W1,UIN)
VALUE(INLET,KE,TKEIN);VALUE(INLET,EP,EPSIN)
IF(JKO.GT.0) THEN
+VALUE(INLET,OMEG,OMIN)
ENDIF
 
PATCH(OUTL,HIGH,#1,#NREGX,#1,#NREGY,#NREGZ,#NREGZ,1,1)
COVAL(OUTL,P1,1.0E3,0.0)
COVAL(OUTL,W1,ONLYMS,0.0);COVAL(OUTL,V1,ONLYMS,0.0)
COVAL(OUTL,KE,ONLYMS,0.0);COVAL(OUTL,EP,ONLYMS,0.0)
 
WALL(WALLN,NORTH,#1,#NREGX,#NREGY,#NREGY,#1,#NREGZ,1,1)
WALL(WALLS,SOUTH,#1,#NREGX,#1,#1,#1,#NREGZ,1,1)
 
    GROUP 15. Termination of sweeps
LSWEEP=1500
    GROUP 16. Termination of iterations
SELREF=T
LITER(P1)=50;LITER(KE)=5;LITER(EP)=5
    GROUP 17. Under-relaxation devices
KELIN=3
REAL(DTF);DTF=ZWLAST/UIN/NZ/2
RELAX(W1,FALSDT,DTF);RELAX(V1,FALSDT,DTF)
RELAX(KE,FALSDT,DTF); RELAX(EP,FALSDT,DTF)
IF(JKO.EQ.1) THEN
+ RELAX(OMEG,FALSDT,DTF)
ENDIF  
IYMON=NYC-4;IXMON=1;IZMON=NZUP+NZC+10;NPRMON=100
    GROUP 23. Field print-out and plot control
ITABL=3;NPLT=10;IPLTL=LSWEEP;NZPRIN=2;NYPRIN=2
TSTSWP=-1
SPEDAT(SET,GXMONI,PLOTALL,L,T)
 
DISTIL=T
STORE(PRPS); EX(PRPS)= 9.659E-01;EX(VPOR)=9.659E-01
CASE :CTURB: OF
WHEN KLM,3
+EX(P1  )=1.776E-01;EX(V1  )=9.811E-02
+EX(W1  )=7.631E-01;EX(KE  )=1.623E-02
+EX(EP  )=2.893E-03;EX(EPKE)=9.602E-02
+EX(YPLS)=7.306E+01
+EX(DWDY)=2.863E-01;EX(DVDZ)=6.170E-02
+EX(FOMG)=6.110E-01;EX(ENUT)=2.569E-02
WHEN KEM,3
+EX(P1  )=1.205E-01;EX(V1  )=7.313E-02 
+EX(W1  )=8.293E-01;EX(KE  )=2.615E-02 
+EX(EP  )=6.198E-03;EX(EPKE)=1.013E-01 
+EX(YPLS)=7.179E+01;EX(ENUT)=3.082E-02 
WHEN MMK,3
+EX(P1  )=1.787E-01;EX(V1  )=9.850E-02
+EX(W1  )=7.610E-01;EX(KE  )=1.618E-02
+EX(EP  )=2.788E-03;EX(EPKE)=9.374E-02
+EX(YPLS)=7.315E+01
+EX(DWDY)=2.871E-01;EX(DVDZ)=6.213E-02
+EX(FOMG)=5.503E-01;EX(ENUT)=1.460E-02
WHEN RNG,3
+EX(P1  )=1.919E-01;EX(V1  )=1.039E-01
+EX(W1  )=7.595E-01;EX(KE  )=1.480E-02
+EX(EP  )=2.418E-03;EX(EPKE)=1.066E-01
+EX(YPLS)=7.333E+01
+EX(ENUT)=2.054E-02
WHEN CKM,3
+EX(P1  )=1.960E-01;EX(V1  )=1.047E-01
+EX(W1  )=7.581E-01;EX(KE  )=1.333E-02
+EX(EP  )=2.434E-03;EX(EPKE)=1.114E-01
+EX(YPLS)=7.349E+01
 EX(ENUT)=1.981E-02;
WHEN KOM,3
+EX(P1  )=1.393E-01;EX(V1  )=8.191E-02 
+EX(W1  )=8.115E-01;EX(KE  )=2.906E-02 
+EX(EP  )=6.036E-03;EX(EPKE)=9.659E-11 
+EX(YPLS)=7.200E+01 
+EX(OMEG)=1.065E+00;EX(ENUT)=3.745E-02 
WHEN RKE,3
+EX(P1  )=1.866E-01;EX(V1  )=1.027E-01
+EX(W1  )=7.593E-01;EX(KE  )=1.526E-02
+EX(EP  )=2.468E-03;EX(VPOR)=9.659E-01
+EX(YPLS)=7.337E+01;EX(C1E )=5.111E-01
+EX(DWDZ)=8.975E-02;EX(DWDY)=2.800E-01
+EX(DVDZ)=6.251E-02;EX(DVDY)=9.097E-02
+EX(EPKE)=9.428E-02;EX(CMU )=8.361E-02
+EX(ENUT)=2.273E-02
ENDCASE