TEXT(Decay Of Turbulence. K-Epsilon Model  
TITLE
  DISPLAY
  Homogeneous turbulence decays behind a grid in a stream of uniform
  velocity. The transport equations of the k-epsilon model reduce to
  ordinary differential equations with analytical solutions of which
  that for k is:
 
    k=( k0**a - a*e0*(x-x0)/{u0*k0**(1-a)} ) ** (1/a)
 
  where k0 and e0 are k and e for x=0, and u0 is the flow velocity;
  a=1.0-C2E. The value of C2E that fits the experimental data is 1.9
 
  The analytical solution is computed in the Q1 file, by way of a DO
  loop, and stored in the variable ANKE. The small differences
  between the computed KE and ANKE result from the upwind
  differencing: ANKE is the cell-centre value, whereas KE is the
  convected-out value.
  ENDDIS
REAL(XLEN,UIN,REYNO,WIDTH,TKEIN,EPIN)
XLEN=1.0;REYNO=5.E5;UIN=10.;WIDTH=0.25
    GROUP 3. X-direction grid specification
GRDPWR(X,200,XLEN,1.0)
    GROUP 7. Variables stored, solved & named
STORE(U1,ENUT);TURMOD(KEMODL);NAME(ENUT)=VIST
    GROUP 8. Terms (in differential equations) & devices
DIFCUT=0.0;ENUL=UIN*WIDTH/REYNO
    GROUP 11. Initialization of variable or porosity fields
TKEIN=1.0;EPIN=7.5
FIINIT(U1)=UIN;FIINIT(KE)=TKEIN;FIINIT(EP)=EPIN
  ** Analytical solution placed in ANKE store
STORE(ANKE);REAL(AA,NN,BB,ANTKE,C2E,GX);INTEGER(JJM1)
C2E=1.92;NN=1./(1.0-C2E);AA=1./NN;BB=1.0-AA
REAL(CON1,CON2,CON3)
CON1 = TKEIN**AA ; CON2 = AA*EPIN/(UIN*TKEIN**BB)
CON3 = XULAST/2.0 ; CON2=CON2 * CON3
DO JJ=1,NX
+ PATCH(IN:JJ:,INIVAL,JJ,JJ,1,NY,1,NZ,1,1)
+ GX=XFRAC(JJ)
+ IF(JJ.NE.1) THEN
+   JJM1=JJ-1; GX=(XFRAC(JJ) + XFRAC(JJM1))
+ ENDIF
+ ANTKE=(CON1 - CON2*GX)**NN
+ INIT(IN:JJ:,ANKE,ZERO,ANTKE)
ENDDO
 
    GROUP 13. Boundary conditions and special sources
  ** Inlet Boundary
PATCH(INLET,WEST,1,1,1,1,1,1,1,1)
COVAL(INLET,KE,RHO1*UIN,TKEIN);COVAL(INLET,EP,RHO1*UIN,EPIN)
    GROUP 16. Termination of iterations
LSWEEP=1000;SELREF=T;RESFAC=0.0001
    GROUP 21. Print-out of variables
OUTPUT(U1,N,N,N,N,N,N);OUTPUT(VIST,N,N,N,N,N,N)
    GROUP 22. Monitor print-out
IXMON=NX
    GROUP 23. Field print-out and plot control
NXPRIN=nx/5;ORSIZ=0.4
PATCH(PROFILE,PROFIL,1,NX,1,1,1,1,1,1)
PLOT(PROFILE,KE,0.0,1.5*TKEIN);PLOT(PROFILE,EP,0.0,2.0*EPIN)
PLOT(PROFILE,VIST,0.0,0.0)
 **END OF LIBRARY CASE 230