TALK=T;RUN( 1, 1)
** LOAD(x208) from the x Input Library
GROUP 1. Run title and other preliminaries
TEXT(LRN MODEL BKWRD-FACING STEP Y-X :T208
TITLE
DISPLAY
The problem concerns 2d incompressible, turbulent flow and heat
transfer over a backward-facing step at a Reynolds number of
45,000 based on step height Hs and inlet velocity. The expansion
ratio Hc/Hs = 3.0, where Hc is the downstream channel height 3.0.
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 or Menter k-w models, or the k-w-SST 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 measured reattachment length of
Xr/H=7.1 by about 15%. The reattachments lengths predicted by the
various turbulence models are compared with the measured value below:
LB CK 2-L KW KWM KWR SST EXPT
Xr/H = 5.3 6.8 6.4 7.3 6.6 7.8 7.2 7.1
The calculations employ a relatively coarse non-uniform mesh of NX=80
by NY=80 which concentrates grid cells close to the wall. The solution
is unlikely to produce grid-independent results because most of the
near-wall y+ values cells are too large, typically ranging from 0.05
to 4. The experimental values are taken from:
Kim,J., Klein,S.J. & Johnston,J.P., "Investigation of separation &
reattachment of a turbulent shear layer : Flow over a backward-facing
step" Report No. MD-37. Thermosciences Division, Dept. of Mech Eng.,
Stanford University (1978).
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 u1 y 1
divide x .0381 1
shift x -4 1
col9 1
level y 0;level x 0
scale x -4 15
msg Velocity (U1) profile
msg Press e to END
ENDUSE
CHAR(CTURB,TLSC)
INTEGER(NYS,NXS,NXS2)
REAL(HEIGHT,WIDTH,CLEN,SLEN,SLEN2,REYNO,UIN,TKEIN,EPSIN,FRIC,MIXL)
REAL(QIN,AREAO,OMEGIN)
BOOLEAN(KWMOD);KWMOD=F
** Calculation of domain specifications
** channel length=0.762m & channel width=0.1143m
HEIGHT=0.0381;WIDTH=3.*HEIGHT
SLEN=4.*HEIGHT;SLEN2=5.5*HEIGHT;CLEN=20.*HEIGHT
REYNO=4.5E4;UIN=13.;FRIC=0.018;MIXL=0.09*HEIGHT
TKEIN=0.25*UIN*UIN*FRIC;EPSIN=0.1643*TKEIN**1.5/MIXL
GROUP 2. Transience; time-step specification
GROUP 3. X-direction grid specification
NXS=14;NXS2=61
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)
GROUP 4. Y-direction grid specification
NYS=40
SUBGRD(Y,1,-NYS,HEIGHT,1.7);SUBGRD(Y,NY+1,-80,2.0*HEIGHT,2.0)
GROUP 5. Z-direction grid specification
GROUP 6. Body-fitted coordinates or grid distortion
GROUP 7. Variables stored, solved & named
SOLVE(P1,U1,V1,TEM1);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)
MESG( Enter the required turbulence model:
MESG( LB - Low-Re Lam-Brem. k-e model (default)
MESG( CK - Low-Re Lam-Brem. k-e model
MESG( 2L - Low-Re 2-layer k-e model
MESG( KW - Low-Re Wilcox 1988 k-w model
MESG( KWR - Low-Re Wilcox 2008 k-w model
MESG( KWM - Low-Re Menter 1992 k-w model
MESG( KWS - Low-Re k-w SST model
MESG(
READVDU(CTURB,CHAR,LB)
CASE :CTURB: OF
WHEN LB,2
+ TEXT(Low-Re LB k-e BF Step Y-X :T208
+ MESG(Low-Re LB k-e turbulence model
+ KELIN=1;TLSC=EP
+ TURMOD(KEMODL-LOWRE);STORE(REYN,FMU,FTWO)
WHEN CK,2
+ TEXT(Low-Re CK k-e BF Step Y-X :T208
+ MESG(Low-Re Chen-Kim k-e turbulence model
+ KELIN=1;TLSC=EP
+ TURMOD(KECHEN-LOWRE);STORE(REYN,FMU,FTWO)
WHEN 2L,2
+ TEXT(2-Layer k-e BF Step Y-X :T208
+ MESG(Low-Re 2-layer k-e turbulence model
+ KELIN=1;TLSC=EP
+ TURMOD(KEMODL-2L);STORE(REYN,FMU,FTWO)
WHEN KW,2
+ TEXT(Low-Re Wilcox 1988 k-w BF Step Y-X :T208
+ MESG(Low-Re Wilcox 1988 k-w model
+ TLSC=OMEG
+ TURMOD(KOMODL-LOWRE);STORE(REYT,FMU,FTWO)
+ STORE(EP);EPSIN=EPSIN/(0.09*TKEIN)
WHEN KWR,3
+ TEXT(Low-Re Wilcox 2008 k-w BF Step Y-X :T208
+ MESG(Low-Re Wilcox 2008 k-w model
+ TLSC=OMEG;KWMOD=T
+ TURMOD(KWMODLR-LOWRE);STORE(REYT)
+ EPSIN=EPSIN/(0.09*TKEIN);FIINIT(FBP)=1.0
WHEN KWM,3
TEXT(Low-Re Menter k-w BF Step Y-X :T208
+ MESG(Low-Re Menter k-w model
+ TURMOD(KWMENTER-LOWRE);TLSC=OMEG
+ STORE(EP);EPSIN=EPSIN/(0.09*TKEIN)
+ STORE(CDWS) ! output cross-diffusion term
WHEN KWS,3
TEXT(Low-Re k-w SST BF StepY-X :T208
+ MESG(Menter 1992 Low-Re k-w SST model
+ TURMOD(KWSST-LOWRE);TLSC=OMEG
+ STORE(EP);EPSIN=EPSIN/(0.09*TKEIN)
+ STORE(CDWS)! output cross-diffusion term
ENDCASE
STORE(YPLS)
GROUP 8. Terms (in differential equations) & devices
TERMS(TEM1,N,P,P,P,P,P)
GROUP 9. Properties of the medium (or media)
RHO1=1.0;ENUL=UIN*HEIGHT/REYNO;PRT(TEM1)=0.86
PRNDTL(TEM1)=1.0;CP1=1000.0
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)
STORE(PRPS)
PATCH(BLOCK,INIVAL,1,NXS,1,NYS,1,1,1,1)
INIT(BLOCK,PRPS,0.,198)
EGWF=T
FIINIT(U1)=0.1*UIN;FIINIT(V1)=0.0;FIINIT(P1)=1.3E-4
FIINIT(:TLSC:)=EPSIN;FIINIT(KE)=0.04*UIN*UIN;FIINIT(TEM1)=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,TEM1,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=1000. & set Qin for unit temperature rise
QIN=RHO1*CP1*2.*HEIGHT*UIN;AREAO=ZWLAST*(CLEN-SLEN)
PATCH(HEATIN,SOUTH,NXS+1,NX,1,1,1,NZ,1,1)
COVAL(HEATIN,TEM1,FIXFLU,QIN/AREAO)
GROUP 14. Downstream pressure for PARAB=.TRUE.
GROUP 15. Termination of sweeps
LSWEEP=1500;TSTSWP=-1
GROUP 16. Termination of iterations
LITER(U1)=2;LITER(V1)=2;LITER(KE)=5;LITER(:TLSC:)=5;LITER(TEM1)=10
GROUP 17. Under-relaxation devices
REAL(DTF);DTF=10.*CLEN/UIN/NX; RELAX(P1,LINRLX,1.0)
CASE :CTURB: OF
WHEN LB,2
** use heavier relaxation on KE
+ KELIN=3; RELAX(KE,LINRLX,0.3); RELAX(EP,LINRLX,0.3)
WHEN CK,2
** use heavier relaxation on KE
+ KELIN=3; RELAX(KE,LINRLX,0.3); RELAX(EP,LINRLX,0.5)
WHEN 2L,2
** incomplete convergence - stalled residuals
+ KELIN=3; RELAX(KE,LINRLX,0.3); RELAX(EP,LINRLX,0.3)
+ RELAX(KE,FALSDT,DTF/10.);RELAX(:TLSC:,FALSDT,DTF/10.)
+ RELAX(ENUT,LINRLX,0.3)
+ OUTPUT(REYN,Y,N,Y,N,Y,Y)
+ OUTPUT(FTWO,Y,N,Y,N,Y,Y)
+ OUTPUT(FMU,Y,N,Y,N,Y,Y)
+(STORED of EPCO is CORR(EP) with IMAT<100)
+(STORED of EPRE is RESI(EP) with IMAT<100)
+ OUTPUT(EPCO,Y,N,Y,N,Y,Y);OUTPUT(EPRE,Y,N,Y,N,Y,Y)
+ IXMON=18;IYMON=39
WHEN KW,2
+ DTF=DTF/10
+ RELAX(KE,FALSDT,DTF); RELAX(:TLSC:,FALSDT,DTF)
WHEN KWR,3
+ DTF=DTF/10;LSWEEP=2000
+ RELAX(KE,FALSDT,DTF);RELAX(:TLSC:,FALSDT,DTF)
WHEN KWM,3
+ DTF=DTF/10
+ RELAX(KE,FALSDT,DTF); RELAX(:TLSC:,FALSDT,DTF)
WHEN KWS,3
+ DTF=DTF/10
+ RELAX(KE,FALSDT,DTF); RELAX(:TLSC:,FALSDT,DTF)
ENDCASE
RELAX(U1,FALSDT,DTF); RELAX(V1,FALSDT,DTF); RELAX(TEM1,FALSDT,100.*DTF)
GROUP 18. Limits on variables or increments to them
GROUP 21. Print-out of variables
SPEDAT(SET,OUTPUT,NOFIELD,L,T)
GROUP 22. Monitor print-out
SPEDAT(SET,GXMONI,PLOTALL,L,T)
ITABL=3;IXMON=NXS+2;IPLTL=2000;IYMON=NYS-2;NPRMON=10000
DISTIL=T
CASE :CTURB: OF
WHEN LB,2
+ 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(TEM1)=4.132E+00
+ EX(PRPS)=9.125E-01; EX(EPKE)=2.147E+02
+ 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; EX(YPLS)=7.523E-02
WHEN CK,2
+ 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(TEM1)=5.842E+00
+ EX(PRPS)=9.125E-01; EX(EPKE)=2.307E+02
+ 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; EX(YPLS)=6.854E-02
WHEN 2L,2
+ 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(TEM1)=5.698E+00
+ EX(PRPS)=9.125E-01;
+ EX(REYN)=2.106E+03; EX(LTLS)=6.673E-04; EX(WDIS)=1.796E-02
+ EX(FTWO)=4.456E+00; EX(FMU )=8.314E-01; EX(LEN1)=2.855E-03
+ EX(ENUT)=2.038E-03; EX(EPKE)=5.982E+02; EX(YPLS)=6.262E-02
+ EX(EPKE)=4.910E+02
WHEN KW,2
+ EX(P1 )=2.313E+01;EX(U1 )=6.828E+00
+ EX(PRPS)=9.125E-01
+ EX(FTWO)=8.430E-01;EX(FMU )=7.824E-01
+ EX(LEN1)=2.812E-03;EX(V1 )=1.865E-01
+ EX(KE )=1.461E+00;EX(EP )=4.378E+02
+ EX(TEM1)=7.865E+00;EX(REYT)=1.647E+02
+ EX(OMEG)=1.317E+04;EX(ENUT)=1.763E-03
+ EX(YPLS)=5.269E-02
WHEN KWR,3
+ EX(P1 )=2.317E+01;EX(U1 )=6.836E+00
+ EX(V1 )=1.715E-01;EX(KE )=1.302E+00
+ EX(EP )=2.629E+02;EX(TEM1)=8.757E+00
+ EX(YPLS)=4.995E-02
+ EX(REYT)=9.124E-11;EX(DVDY)=1.434E+01
+ EX(DVDX)=6.857E+00;EX(DUDY)=1.018E+03
+ EX(DUDX)=1.416E+01;EX(GEN1)=1.356E+09
+ EX(FBP )=9.115E-01;EX(XWP )=2.779E-04
+ EX(OMEG)=1.407E+04;EX(PRPS)=9.125E-01
+ EX(LEN1)=2.602E-03;EX(ENUT)=1.553E-03
WHEN KWM,3
+ EX(P1 )=2.206E+01;EX(U1 )=6.865E+00
+ EX(V1 )=2.007E-01;EX(KE )=1.476E+00
+ EX(EP )=2.756E+02;EX(TEM1)=7.585E+00
+ EX(YPLS)=5.109E-02
+ EX(CDWS)=6.780E+04;EX(LTLS)=6.672E-04
+ EX(WDIS)=1.796E-02;EX(BF1 )=6.030E-01
+ EX(OMEG)=1.333E+04;EX(PRPS)=9.125E-01
+ EX(LEN1)=2.681E-03;EX(ENUT)=1.762E-03
WHEN KWS,3
+ EX(P1 )=2.243E+01;EX(U1 )=6.858E+00
+ EX(V1 )=1.858E-01;EX(KE )=1.399E+00
+ EX(EP )=2.645E+02;EX(TEM1)=8.092E+00
+ EX(YPLS)=4.956E-02
+ EX(CDWS)=6.546E+04;EX(LTLS)=6.672E-04
+ EX(WDIS)=1.796E-02;EX(GEN1)=2.738E+09
+ EX(BF2 )=8.100E-01;EX(BF1 )=5.905E-01
+ EX(OMEG)=1.332E+04;EX(PRPS)=9.125E-01
+ EX(LEN1)=2.613E-03;EX(ENUT)=1.604E-03
ENDCASE
STOP