GROUP 1. Run title and other preliminaries
TEXT(LAM-BRE KE_1D PLANE COUETTE FLOW  :T205
TITLE
mesg(PC486/50 time last reported as appx. 40.sec
  DISPLAY
  The problem considered is plane turbulent couette flow in a
  channel at Re=1.E5, as described in detail for library case T100.
 
  The turbulence is simulated by use of the Lam-Bremhorst low-Re
  k-e model. The calculation integrates down to the wall and the
  solution is performed by use of the single-slab solver. A non-
  uniform grid is employed so as to concentrate cells very close to
  the walls. The calculation may also be performed with the Chen-
  Kim low-Re k-e model and the Wilcox low-Re k-omega model.
  ENDDIS
  The Lam-Bremhorst and Chen-Kim models predict skin-friction
  coefficients Cf of 3.78E-3 and 3.50e-3 respectively, which are in
  reasonable agreement with the experimental value of 3.07E-3
  suggested by the data of El Telbany and Reynolds [1982]. The
  k-omega model produces much closer agreement, yielding
  Cf=2.931e-3, where Cf=2.*tauw/(rho*Uav**2) and Uav is the
  average velocity.
 
   The following AUTOPLOT use file produces three plots;
   the first is the axial velocity profile; the second
   is the turbulence energy profile; and the third is the
   dissipation rate profile.
 
   AUTOPLOT USE
   file
   phi 5
 
   da 1 w1;col9 1
   msg Velocity (W1) profile
   msg Press RETURN to continue
   pause
   clear
   da 1 ke;col9 1
   msg KE profile
   msg Press RETURN to continue
   pause
   clear
   da 1 ep;col9 1
   msg EP profile
   msg Press e to END
   ENDUSE
 
CHAR(CTURB,TLSC);BOOLEAN(VARLAM);VARLAM=T
REAL(HEIGHT,WTOP,REY,TKEIN,EPSIN,MIXL,DTF)
REAL(WAV,US,MASIN,DELT1,DELY,KFAC,AA)
INTEGER(NY2,JJM,JJJ)
HEIGHT=0.1;WTOP=1.0; REY=1.E5;WAV=0.5*WTOP
  ** US from data of El Telbany & Reynolds [1982]
US=WAV*0.196/LOG10(REY);TKEIN=US*US/.3
MIXL=0.045*HEIGHT;EPSIN=TKEIN**1.5/MIXL*0.1643
    GROUP 4. Y-direction grid specification
ENULA=WAV*HEIGHT/REY
  ** define first dely from wall and the grid-expansion
     factor Kfac which defines a constant ratio of lengths of
     two adjacent cells.
DELT1=0.5*ENULA/US;KFAC=1.08;DELY=DELT1/(0.5*HEIGHT)
  ** calculate NY from dely & Kfac
AA=(0.5/DELY)*(KFAC-1.0)+1.0;AA=LOG(AA)/LOG(KFAC)+1.0001
NY2=AA;NY=2*NY2
  ** define uniform grid initially
IREGY=1;GRDPWR(Y,NY,YVLAST,1.0)
  ** compute expanding grid from south boundary over one
     half of the channel width
YFRAC(1)=DELY
DO JJ=2,NY2
+ JJM=JJ-1
+ DELY=KFAC*DELY
+ YFRAC(JJ)=YFRAC(JJM)+DELY
ENDDO
YFRAC(NY2)=0.5
  ** create symmetrical grid in the second half of the channel
JJJ=0
DO JJ=NY-1,NY2+1,-1
+ JJJ=JJJ+1
+ YFRAC(JJ)=1.-YFRAC(JJJ)
ENDDO
YFRAC(NY)=1.0;YVLAST=HEIGHT
    GROUP 7. Variables stored, solved & named
SOLVE(W1);STORE(ENUT,LEN1);SOLUTN(W1,P,P,P,P,P,N)
MESG( Enter the required turbulence model:
MESG(  CK   -  Chen-Kim low-Re k-e model
MESG(  LB   -  Lam-Bemhorst low-Re k-e model (default)
MESG(  KO   -  Wilcox low-re k-omega model
MESG(
READVDU(CTURB,CHAR,LB)
CASE :CTURB: OF
WHEN CK,2
+ TEXT(CHEN-KIM KE_1D PLANE COUETTE FLOW :T205
+ MESG(Chen-Kim low-Re k-e model
+ TURMOD(KECHEN-LOWRE);KELIN=1;TLSC=EP
WHEN LB,2
+ MESG(Lam-Bremhorst low-Re k-e model
+ TURMOD(KEMODL-LOWRE);KELIN=1;TLSC=EP
WHEN KO,2
+ TEXT(K-OMEGA_1D PLANE COUETTE FLOW   :T205
+ MESG(k-omega low-Re model
+ TURMOD(KOMODL-LOWRE);TLSC=OMEG
+ STORE(EP);EPSIN=EPSIN/(0.09*TKEIN)
ENDCASE
    GROUP 8. Terms (in differential equations) & devices
  ** Deactivate convection for single-slab solution
TERMS(W1,N,N,P,P,P,P);TERMS(KE,Y,N,P,P,P,P)
TERMS(:TLSC:,Y,N,P,P,P,P)
    GROUP 9. Properties of the medium (or media)
ENUL=ENULA
   ** test for ground-set enul
IF(VARLAM) THEN
+ TMP1=CONST;TMP1A=0.0;ENUL=LINTEM
ENDIF
    GROUP 11. Initialization of variable or porosity fields
FIINIT(:TLSC:)=EPSIN;FIINIT(KE)=TKEIN
PATCH(ICOUF,LINVLY,1,1,1,NY,1,NZ,1,1)
INIT(ICOUF,W1,WTOP/HEIGHT,0.0)
    GROUP 13. Boundary conditions and special sources
  ** moving upper wall
WALL(WALLN,NORTH,1,1,NY,NY,1,NZ,1,1);COVAL(WALLN,W1,LOGLAW,WTOP)
  ** stationary bottom wall
WALL(WALLS,SOUTH,1,1,1,1,1,NZ,1,1)
    GROUP 15. Termination of sweeps
LSWEEP=70;TSTSWP=-1;LITHYD=6
    GROUP 16. Termination of iterations
MASIN=RHO1*WAV*HEIGHT; RESREF(W1)=1.E-12*MASIN*WAV
RESREF(KE)=RESREF(W1)*TKEIN; RESREF(:TLSC:)=RESREF(W1)*EPSIN
    GROUP 17. Under-relaxation devices
DTF=0.1*ZWLAST/WAV
RELAX(W1,FALSDT,DTF); RELAX(KE,FALSDT,DTF/4.)
RELAX(:TLSC:,FALSDT,DTF/4.)
    GROUP 18. Limits on variables or increments to them
VARMIN(W1)=1.E-10
    GROUP 22. Spot-value print-out
IYMON=2;NPLT=5;NZPRIN=1;NYPRIN=2;IYPRF=1
    GROUP 24. Dumps for restarts
STORE(FMU,REYT,REYN,FONE,FTWO);WALPRN=T