TALK=T;RUN( 1, 1)
** LOAD(x110) from the x Input Library
GROUP 1. Run title and other preliminaries
TEXT(K-W_BKWRD FACING STEP Y-Z :T110
TITLE
DISPLAY
This simulation concerns 2d incompressible, turbulent flow
over a backward-facing step in a closed channel. The edge of the
step provides a fixed point of flow separation.The case is
similar to library case 290, except the calculation is performed
in the Y-Z plane.
The geometry considered has a expansion ratio of 1.5 and an
exit width of 3h, where h is the step height. The calculations
are started at x= -4h and terminated at x=16h, where a fixed-
pressure boundary condition is applied. The Reynolds number based
on the step height is 45,000 and uniform profiles of u, k and e
are specified at the inlet to the calculation domain. Scalable
equilibrium wall functions are used in the simulations.
For this case, the standard form of the k-e model is known to
underpredict the reattachment length XR by about 14%. Calculations
may be made with the high-Re forms of the standard, Chen-Kim, RNG
and realisable k-e models; the Wilcox 1988 & 2008 k-w models, and
the Menter k-w and k-w SST models.
ENDDIS
The calculation employs a non-uniform mesh of NY=60 by NZ=90 for
which the standard k-e model predicts XR/H=5.9 and the measurements
indicate XR/H=7.1. The mesh sensitivity of the computed results
shown below has not been investigated, but certainly it isn't fine
enough around the vicinity of the step to resolve the experimentally-
observed secondary-separation region in the corner just downstream
of the step.
K-E RKE CK RNG KW KWR KWM KW-SST EXPT
XR/H = 5.9 7.0 7.6 6.7 5.3 6.3 5.2 6.2 7.1
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
phida 3
da 1 W1 y 1
divide x .0381 1
shift x -4 1
col9 1
level y 0;level x 0
scale x 0 15
msg Velocity (W1) profile
msg Press e to END
ENDUSE
CHAR(CTURB,TLSC)
REAL(HEIGHT,WIDTH,CLEN,SLEN,REYNO,UIN,TKEIN,EPSIN,MIXL,FRIC,OMEGIN)
INTEGER(NYS,NZS,NZDS,NYCH)
BOOLEAN(KWMOD);KWMOD=F
** Calculation of domain specifications
HEIGHT=0.0381;WIDTH=3.*HEIGHT;SLEN=4.*HEIGHT;CLEN=20.*HEIGHT
REYNO=4.5E4;UIN=13.
FRIC=0.018;TKEIN=0.25*UIN*UIN*FRIC;MIXL=0.09*HEIGHT
EPSIN=0.1643*TKEIN**1.5/MIXL
OMEGIN=EPSIN/(0.09*TKEIN)
GROUP 3. X-direction grid specification
GROUP 4. Y-direction grid specification
** channel length = 0.762 & channel width = 0.1143
NZS=15;NZDS=75
NREGZ=2
IREGZ=1;GRDPWR(Z,NZS,SLEN,1.0)
IREGZ=2;GRDPWR(Z,NZDS,-(CLEN-SLEN),1.02)
NYS=25;NYCH=35
NREGY=2
IREGY=1;GRDPWR(Y,-NYS,-HEIGHT,1.05)
IREGY=2;GRDPWR(Y,-NYCH,-(WIDTH-HEIGHT),1.05)
GROUP 7. Variables stored, solved & named
SOLVE(P1,W1,V1);SOLUTN(P1,Y,Y,Y,N,N,N);STORE(ENUT)
SOLUTN(W1,P,P,P,P,P,N);SOLUTN(V1,P,P,P,P,P,N)
MESG( Enter the required turbulence model:
MESG( CK - Chen-Kim k-e model
MESG( KE - Standard k-e model
MESG( RNG - RNG k-e model
MESG( RKE - Realisable k-e model
MESG( KW - Wilcox 1988 k-w model (default)
MESG( KWR - Wilcox 2008 k-w model
MESG( KWM - Menter 1992 k-w model
MESG( KWS - k-w SST model
MESG(
READVDU(CTURB,CHAR,KOW)
CASE :CTURB: OF
WHEN CK,2
+ TEXT(CHEN KIM K-E_BKWRD FACING STEP Y-Z :T110
+ MESG(Chen-Kim k-e model
+ TURMOD(KECHEN);TLSC=EP
WHEN KE,2
+ TEXT(K-E_BKWRD FACING STEP Y-Z :T110
+ MESG(Standard k-e model
+ TURMOD(KEMODL);TLSC=EP
WHEN RNG,3
+ TEXT(RNG K-E_BKWRD FACING STEP Y-Z :T110
+ MESG(RNG k-e model
+ TURMOD(KERNG);TLSC=EP
+ STORE(ETA,ALF,GEN1)
+ OUTPUT(ALF,Y,N,P,Y,Y,Y);OUTPUT(ETA,Y,N,P,Y,Y,Y)
WHEN RKE,3
+ TEXT(Realisable K-E_BKWRD FACING STEP Y-Z :T110
+ MESG(Realisable k-e model
+ TURMOD(KEREAL);TLSC=EP;STORE(C1E)
+ OUTPUT(CMU,P,P,P,P,Y,Y);OUTPUT(C1E,P,P,P,P,Y,Y)
WHEN KW,2
+ TEXT(K-W_1988_BKWRD FACING STEP Y-Z :T110
+ MESG(Wilcox 1988 k-w model (default)
+ TURMOD(KWMODL);TLSC=OMEG
+ OMEGIN=EPSIN/(0.09*TKEIN);KWMOD=T
WHEN KWR,3
+ TEXT(K-W_2008_BKWRD FACING STEP Y-Z :T110
+ MESG(Wilcox 2008 k-w model
+ TURMOD(KWMODLR);TLSC=OMEG;STORE(FBP);FIINIT(FBP)=1.0
+ OMEGIN=EPSIN/(0.09*TKEIN);KWMOD=T
** Solve whole-field to accelerate convergence
+ SOLUTN(V1,Y,Y,Y,P,P,P);SOLUTN(W1,Y,Y,Y,P,P,P)
+ SOLUTN(KE,Y,Y,Y,P,P,P);SOLUTN(OMEG,Y,Y,Y,P,P,P)
WHEN KWM,3
TEXT(Menter k-w BKWRD FACING STEP Y-Z :T110
+ MESG(Menter 1992 k-w model
+ TURMOD(KWMENTER);TLSC=OMEG
+ FIINIT(BF1)=1.0
+ OMEGIN=EPSIN/(0.09*TKEIN);KWMOD=T
+ STORE(BF1,SIGK,SIGW,DKDY,DKDZ,DFDY,DFDZ,CDWS)
+ STORE(CWAL,CWBE)
WHEN KWS,3
TEXT(k-w SST BKWRD FACING STEP Y-Z :T110
+ MESG(Menter k-w SST model
+ TURMOD(KWSST);TLSC=OMEG
+ OMEGIN=EPSIN/(0.09*TKEIN);KWMOD=T
+ STORE(BF1,BF2,GEN1,SIGK,SIGW,CDWS)
+ STORE(CWAL,CWBE)
+ FIINIT(BF1)=1.0;FIINIT(BF2)=1.0
ENDCASE
** for output purposes only
STORE(SKIN,YPLS,VABS,STRS)
GROUP 8. Terms (in differential equations) & devices
GROUP 9. Properties of the medium (or media)
RHO1=1.0;ENUL=UIN*HEIGHT/REYNO
GROUP 11. Initialization of variable or porosity fields
FRIC=0.018;TKEIN=0.25*UIN*UIN*FRIC;MIXL=0.09*HEIGHT
EPSIN=0.1643*TKEIN**1.5/MIXL;FIINIT(W1)=UIN;FIINIT(P1)=1.3E-4
FIINIT(KE)=TKEIN;FIINIT(EP)=EPSIN;FIINIT(V1)=0.001*UIN
IF(KWMOD) THEN
+ FIINIT(OMEG)=OMEGIN
ENDIF
** Initialization of variables in blocked region
STORE(PRPS)
PATCH(STEP,INIVAL,1,NX,1,NYS,1,NZS,1,1)
INIT(STEP,PRPS,0.,198)
EGWF=T
GROUP 13. Boundary conditions and special sources
INLET(INLET,LOW,#1,#1,#2,#NREGY,#1,#1,1,1)
VALUE(INLET,P1,UIN);VALUE(INLET,W1,UIN)
VALUE(INLET,KE,TKEIN);VALUE(INLET,EP,EPSIN)
IF(KWMOD) THEN
+VALUE(INLET,OMEG,OMEGIN)
ENDIF
PATCH(OUTLET,HIGH,#1,#1,#1,#NREGY,#NREGZ,#NREGZ,1,1)
COVAL(OUTLET,P1,1.0E5,0.0)
COVAL(OUTLET,W1,ONLYMS,0.0);COVAL(OUTLET,V1,ONLYMS,0.0)
COVAL(OUTLET,KE,ONLYMS,0.0);COVAL(OUTLET,:TLSC:,ONLYMS,0.0)
** scalable wall functions
SCALWF=T
** Use default equilibrium wall functions, although ideally
non-equilibrium wall functions should be used, which
are available only for k-e models.
WALL (WFNN,NORTH,#1,#1,#NREGY,#NREGY,#1,#NREGZ,1,1)
WALL (WFNS,SOUTH,#1,#1,#1,#1,#2,#NREGZ,1,1)
IF(IENUTA.EQ.10.OR.IENUTA.EQ.17.OR.IENUTA.EQ.19) THEN
+ COVAL(WFNN,OMEG,GRND2,GRND2);COVAL(WFNN,EP,0.0,0.0)
+ COVAL(WFNS,OMEG,GRND2,GRND2);COVAL(WFNS,EP,0.0,0.0)
ENDIF
GROUP 15. Termination of sweeps
LSWEEP=1000
GROUP 16. Termination of iterations
REAL(MASIN,DTF);MASIN=2.*HEIGHT*UIN*RHO1; RESREF(P1)=1.E-12*MASIN
RESREF(W1)=RESREF(P1)*UIN; RESREF(V1)=RESREF(W1)
RESREF(KE)=RESREF(P1)*TKEIN; RESREF(EP)=RESREF(P1)*EPSIN
IF(KWMOD) THEN
+ RESREF(OMEG)=RESREF(P1)*OMEGIN
ENDIF
GROUP 17. Under-relaxation devices
DTF=ZWLAST/UIN; RELAX(W1,FALSDT,DTF); RELAX(V1,FALSDT,DTF)
IF(:TLSC:.EQ.EP) THEN
+ KELIN=3
+ RELAX(KE,LINRLX,0.3);RELAX(EP,LINRLX,0.3)
+ DTF=DTF/NZ
+ RELAX(W1,FALSDT,DTF); RELAX(V1,FALSDT,DTF)
ELSE
+ RELAX(KE,FALSDT,DTF/2.);RELAX(:TLSC:,FALSDT,DTF/2.)
+ DTF=5.*DTF/NZ
+ RELAX(V1,FALSDT,DTF);RELAX(W1,FALSDT,DTF)
ENDIF
CASE :CTURB: OF
WHEN KWR,3
** deactivate scalable wall functions & relax the
vortex-stretching factor
+ SCALWF=F;RELAX(FBP,LINRLX,0.01)
+ RELAX(W1,FALSDT,DTF);RELAX(V1,FALSDT,DTF)
+ RELAX(KE,FALSDT,2.*DTF);RELAX(OMEG,FALSDT,2.*DTF)
+ LSWEEP=1000
WHEN KWM,3
+ OUTPUT(BF1,Y,N,Y,N,Y,Y)
+ RELAX(BF1,LINRLX,0.01)
WHEN KWS,3
+ RELAX(BF1,LINRLX,0.01)
+ OUTPUT(BF1,Y,N,Y,N,Y,Y);OUTPUT(BF2,Y,N,Y,N,Y,Y)
WHEN RKE,3
+ RELAX(ENUT,LINRLX,0.3)
ENDCASE
LITER(P1)=50
GROUP 23. Field print-out and plot control
ITABL=3;NPLT=10;IPLTL=LSWEEP;NZPRIN=2;NYPRIN=2
IYMON=NYS-2;IZMON=NZS+2;NPRMON=100
TSTSWP=-1;LITER(KE)=5;LITER(:TLSC:)=5
SPEDAT(SET,GXMONI,PLOTALL,L,T)
SPEDAT(SET,OUTPUT,NOFIELD,L,T)
DISTIL=T
CASE :CTURB: OF
WHEN CK,2
+EX(P1 )=1.951E+01;EX(V1 )=2.056E-01
+EX(W1 )=8.090E+00;EX(KE )=1.351E+00
+EX(EP )=1.400E+02;EX(EPKE)=7.087E+01
+EX(PRPS)=9.306E-01;EX(STRS)=8.007E-03
+EX(VABS)=8.103E+00;EX(YPLS)=9.739E-01
+EX(SKIN)=2.415E-04;EX(ENUT)=2.056E-03
WHEN KE,2
+EX(P1 )=1.742E+01;EX(V1 )=2.317E-01
+EX(W1 )=8.101E+00;EX(KE )=1.625E+00
+EX(EP )=1.525E+02;EX(EPKE)=6.773E+01
+EX(PRPS)=9.306E-01;EX(STRS)=8.001E-03
+EX(VABS)=8.118E+00;EX(YPLS)=9.765E-01
+EX(SKIN)=2.394E-04;EX(ENUT)=2.553E-03
WHEN RNG,3
+EX(P1 )=1.881E+01;EX(V1 )=2.198E-01
+EX(W1 )=8.088E+00;EX(KE )=1.486E+00
+EX(EP )=1.406E+02;EX(EPKE)=7.062E+01
+EX(PRPS)=9.306E-01;EX(STRS)=7.877E-03
+EX(VABS)=8.103E+00;EX(YPLS)=9.695E-01
+EX(SKIN)=2.404E-04;EX(GEN1)=3.809E+05
+EX(ALF )=2.865E-01;EX(ETA )=1.990E+00
+EX(ENUT)=2.177E-03
WHEN RKE,3
+EX(P1 )=1.914E+01;EX(V1 )=2.093E-01
+EX(W1 )=8.073E+00;EX(KE )=1.413E+00
+EX(EP )=1.381E+02;EX(PRPS)=9.306E-01
+EX(STRS)=7.789E-03;EX(VABS)=8.087E+00
+EX(YPLS)=9.630E-01;EX(SKIN)=2.418E-04
+EX(C1E )=4.105E-01;EX(DWDZ)=9.498E+00
+EX(DWDY)=1.747E+02;EX(DVDZ)=1.911E+00
+EX(DVDY)=9.507E+00;EX(EPKE)=7.196E+01
+EX(CMU )=1.387E-01;EX(ENUT)=2.810E-03
WHEN KW,2
+EX(P1 )=1.620E+01;EX(V1 )=2.393E-01
+EX(W1 )=8.111E+00;EX(KE )=1.849E+00
+EX(EP )=1.731E+02;EX(PRPS)=9.306E-01
+EX(STRS)=8.366E-03;EX(VABS)=8.130E+00
+EX(YPLS)=1.023E+00;EX(SKIN)=2.331E-04
+EX(OMEG)=7.219E+02;EX(ENUT)=2.870E-03
WHEN KWR,3
+EX(P1 )=1.748E+01;EX(V1 )=2.219E-01
+EX(W1 )=8.105E+00;EX(KE )=1.624E+00
+EX(EP )=1.623E+02;EX(PRPS)=9.306E-01
+EX(STRS)=8.453E-03;EX(VABS)=8.122E+00
+EX(YPLS)=9.802E-01;EX(SKIN)=4.398E-04
+EX(DWDZ)=1.002E+01;EX(DWDY)=1.635E+02
+EX(DVDZ)=2.605E+00;EX(DVDY)=9.993E+00
+EX(GEN1)=6.680E+05;EX(FBP )=9.286E-01
+EX(XWP )=6.408E+01;EX(OMEG)=7.446E+02
+EX(ENUT)=2.362E-03
WHEN KWM,3
+EX(P1 )=1.622E+01;EX(V1 )=2.365E-01
+EX(W1 )=8.104E+00;EX(KE )=1.729E+00
+EX(EP )=1.649E+02;EX(PRPS)=9.306E-01
+EX(STRS)=8.089E-03;EX(VABS)=8.124E+00
+EX(YPLS)=9.934E-01;EX(SKIN)=2.364E-04
+EX(CWBE)=7.245E-02;EX(CWAL)=4.763E-01
+EX(CDWS)=1.540E+05;EX(DFDZ)=4.050E+03
+EX(DFDY)=1.781E+05;EX(DKDZ)=1.210E+01
+EX(DKDY)=1.377E+02;EX(SIGW)=1.578E+00
+EX(SIGK)=1.520E+00;EX(LTLS)=8.345E-04
+EX(WDIS)=2.257E-02;EX(BF1 )=5.898E-01
+EX(OMEG)=7.505E+02;EX(ENUT)=2.558E-03
WHEN KWS,3
+EX(P1 )=1.743E+01;EX(V1 )=2.194E-01
+EX(W1 )=8.096E+00;EX(KE )=1.574E+00
+EX(EP )=1.544E+02;EX(PRPS)=9.306E-01
+EX(STRS)=7.760E-03;EX(VABS)=8.113E+00
+EX(YPLS)=9.628E-01;EX(SKIN)=2.412E-04
+EX(CWBE)=7.263E-02;EX(CWAL)=4.737E-01
+EX(CDWS)=1.435E+05;EX(SIGW)=1.559E+00
+EX(SIGK)=1.498E+00;EX(LTLS)=8.345E-04
+EX(WDIS)=2.257E-02;EX(GEN1)=8.740E+05
+EX(BF2 )=8.577E-01;EX(BF1 )=5.672E-01
+EX(OMEG)=7.492E+02;EX(ENUT)=2.264E-03
ENDCASE