Encyclopaedia Index
TALK=f;RUN( 1, 1);VDU=VGAMOUSE

TEXT( Library case Y623: CONFINED JET FLOW: 17 FLUID MFM

  >>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
    PLANT information :
     * Data input groups used: 13, 19
     * Ground groups planted : 13, 19-6
     * Headings used  : SORC??, SC06??
     * Functions used : SUM
     * Commands used  : IF, REGION
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<

REAL(HIN,GMIXL,CLEN,WIDTH,WIN1,WIN2,REYNO,WD2)
REAL(TKEIN1,EPIN1,TKEIN2,EPIN2)
INTEGER(IYJ);IYJ=3
REYNO=1.E6;WIDTH=0.3;HIN=1.
WIN1=10.;WIN2=2.0
    GROUP 3. X-direction grid specification
CARTES=F;XULAST=0.1
    GROUP 4. Y-direction grid specification
NY=15;WD2=0.5*WIDTH;GRDPWR(Y,NY,WD2,1.0)
    GROUP 5. Z-direction grid specification
NZ=20;CLEN=20.*WD2;GRDPWR(Z,NZ,CLEN,1.0)
    GROUP 7. Variables stored, solved & named
  ** H1 - single fluid concentration variable;
      G - square of concentration fluctuation;
SOLVE(P1,W1,V1,H1,G)
  ** GENG - production for the square of concentration fluctuation
STORE(ENUT,LEN1,GEN1,EPKE,GENG)
SOLUTN(P1,Y,Y,Y,N,N,N)
TURMOD(KEMODL)
    GROUP 8. Terms (in differential equations) & devices
TERMS(H1,N,Y,Y,Y,Y,Y)
TERMS( G,N,Y,Y,Y,Y,Y)
    GROUP 9. Properties of the medium (or media)
RHO1=1.0;ENUL=WIN1*WIDTH/REYNO
PRT(H1)= 0.86;PRNDTL(H1)= 0.71
PRT(G) = 0.7 ;PRNDTL(G) = 0.7
    GROUP 11. Initialization of variable or porosity fields
FIINIT(W1)=0.5*(WIN1+WIN2);FIINIT(H1)=HIN;FIINIT(LEN1)=0.1*YVLAST
FIINIT(ENUT)=0.01*WIN1*YVLAST
  ** TKEIN = 0.25*WIN1*WIN1*FRIC where FRIC=0.018 AT REYNO=1.E5
TKEIN1=0.25*WIN1*WIN1*0.018
TKEIN2=0.25*WIN2*WIN2*0.018
FIINIT(KE)=0.5*(TKEIN1+TKEIN2)
  ** EPIN = 0.1643*KIN**1.5/LMIX where LMIX=0.045*WIDTH
GMIXL=0.011*WD2
EPIN2=TKEIN2**1.5/GMIXL*0.1643
EPIN1=TKEIN1**1.5/GMIXL*0.1643
FIINIT(EP)=0.5*(EPIN1+EPIN2)
FIINIT(P1)=1.3E-4
    GROUP 13. Boundary conditions and special sources
  ** Inlet Boundaries
INLET(IN1,LOW,1,1,1,IYJ,1,1,1,1)
VALUE(IN1,P1 , WIN1)
VALUE(IN1,W1 , WIN1)
VALUE(IN1,H1 , 1.0)
VALUE(IN1,KE , TKEIN1)
VALUE(IN1,EP , EPIN1)
VALUE(IN1,G  , 0.0)
 INLET(IN2,LOW,1,1,IYJ+1,NY,1,1,1,1)
 VALUE(IN2,P1, WIN2)
 VALUE(IN2,W1, WIN2)
 VALUE(IN2,H1, 0.0)
 VALUE(IN2,KE, TKEIN2)
 VALUE(IN2,EP, EPIN2)
 VALUE(IN2,G ,  0.0)
  **Outlet boundary
PATCH(OUTLET,HIGH,1,NX,1,NY,NZ,NZ,1,1)
COVAL(OUTLET,P1,1.0e05,0.0)
COVAL(OUTLET,W1,ONLYMS,0.0);COVAL(OUTLET,V1,ONLYMS,0.0)
COVAL(OUTLET,KE,ONLYMS,0.0);COVAL(OUTLET,EP,ONLYMS,0.0)
  **North-Wall boundary
WALL (WFNN,NORTH,1,NX,NY,NY,1,NZ,1,1)
  **Source term for g
PATCH(SORG,VOLUME,1,NX,1,NY,1,NZ,1,1)
   CO=2.0*:RHO1:*EPKE
   VAL=GENG/(2.0*:RHO1:*EPKE+TINY)
COVAL(SORG,G  , GRND  ,GRND )
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
    Above statements contain the formulae for combined source/sink term
    for the production and dissipation of the concentration
    fluctuations.
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
PATCH(WG,VOLUME,1,NX,NY,NY,1,NZ,1,1)
   VAL=GENG/(2.0*:RHO1:*EPKE+TINY)
COVAL(WG,G  , FIXVAL  ,GRND )
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
    Production of G is made equal to its dissipation at North-Wall
    boundary by FIXVALing its value to the production rate divided by
    twice product of density and EPKE. The latter is built-in variable
    standing for EP/KE.
  <<<<<<<<<<<<<<<<<<<<<<<<< Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
    GROUP 15. Termination of sweeps
LSWEEP=250
RESFAC=0.01
  ** Provide the re-calculation of reference residuals for G
   RES=SUM(VOL*(GENG-2.*:RHO1:*EPKE*G)/(NY*NZ))
   RESREF(G)=RES
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
    By above two statements the reference residuals for G is calculated
    at the and of z-slab as a sum of its generation rate per cell.
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
    GROUP 16. Termination of iterations
LITHYD=10
    GROUP 17. Under-relaxation devices
KELIN=3
RELAX(P1,LINRLX,0.25)
RELAX(V1,FALSDT,0.025);RELAX(W1,FALSDT,0.025)
RELAX(KE,FALSDT,0.025);RELAX(EP,FALSDT,0.025)
RELAX(G,FALSDT ,0.025)
    GROUP 19. Data communicated by SATELLITE to GROUND
   **Calculation of GENG
     * Auxiliary variables
STORE(DFZ,DFY,DFZH,DFYN)
FIINIT(DFZ) =0.0;FIINIT(DFY) =0.0
FIINIT(DFZH)=0.0;FIINIT(DFYN)=0.0
   DFZ=((H1[,,+1]-H1)/DZGNZ)**2
   REGION(1,NX,1,NY,1,NZ-1,1,1)
   DFY=((H1[,+1,]-H1)/DYG2D)**2
   REGION(1,1,1,NY-1,1,NZ,1,1)
   DFZH=((H1-H1[,,-1])/DZGNZ[,,-1])**2
   REGION(1,NX,1,NY,NZ,NZ,1,1)
   DFYN=((H1-H1[,-1,])/DYG2D[,-1,])**2
   REGION(1,NX,NY,NY,1,NZ,1,1)
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
   The above statements calculate the square concentration derivatives
   separately for internal  and  near domain bounadary regions.
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
   GENG=2.8*:RHO1:*ENUT*(DFZ+DFY+DFZH+DFYN)
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
   The sum of the radial and longitudinal derivatives is multiplied by
   density and turbulent viscosity times 2.8 to get the generation term.
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
  **Output calculations
    GG  - concentration fluctuation;
    GGF - concentration fluctuation normalised by local concentration
          of central jet fluid.
STORE(GG,GGF)
FIINIT(GGF)=0.0
   GG=SQRT(G)
  IF(ISWEEP.EQ.LSWEEP)
   GGF=GG/(H1+TINY)
  IF(ISWEEP.EQ.LSWEEP)
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
    At the end of z-slab for the last sweep the concentration
    fluctuation is calculated and normalized by the local average
    concentration.
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
    GROUP 21. Print-out of variables
WALPRN=T;OUTPUT(KE,Y,Y,Y,Y,Y,Y);OUTPUT(H1,Y,Y,Y,Y,Y,Y)
    GROUP 22. Monitor print-out
IZMON=NZ-1;IYMON=NY-1;UWATCH=T
    GROUP 23. Field print-out and plot control
NPLT=1;NZPRIN=1
NYPRIN=1;IYPRF=1;IYPRL=30
TSTSWP=-1
NAMSAT=MOSG
     ********************     MFTM section    ***************
     ** Number of fluids in population
INTEGER(NFLUIDS)
NFLUIDS=17
     ** Micro-mixing constant
REAL(MMC)
MMC = 5.0 ; RG(1) = MMC
     ** Solve for fluid mass fractions F1, F2, ..., F17
DO II=1,NFLUIDS
 SOLVE(F:II:)
 TERMS(F:II:,N,Y,y,y,y,y)
 PRT(F:II:)= 0.86;PRNDTL(F:II:)= 0.71
 RELAX(F:II:,linrlx,0.15)
 VARMIN(F:II:)=0.0;VARMAX(F:II:)=1.0
 PATCH(PROF:II:,PROFIL,1,1,1,1,1,20,1,1)
 PLOT(PROF:II:,F:II:,0.000E+00, 0.000E+00)
ENDDO
ABSIZ=0.5; ORSIZ=0.2
     ** Fluid population boundary conditions
DO II=1,NFLUIDS
 VALUE(IN1,F:II:,0.0)
 VALUE(IN2,F:II:,0.0)
ENDDO
VALUE(IN1,F1 , 1.0)
VALUE(IN2,F:NFLUIDS:, 1.0)
     ** Coupling/splitting rates
PATCH(MIX,PHASEM,1,NX,1,NY,1,NZ,1,1)
      *  Fluid 1
   CO=RG(1)*EPKE*(F3+F5+F7+F9+F11+F13+F15+F17)
COVAL(MIX,F1   , GRND  ,0.0 )
      *  Fluid 2
   VAL=2.*RG(1)*EPKE*(F1*F3)-$
           RG(1)*EPKE*(F4+F6+F8+F10+F12+F14+F16)*F2
COVAL(MIX,F2, FIXFLU,GRND)
       *  Fluid 3
   VAL=2.*RG(1)*EPKE*(F2*F4+F1*F5)-$
           RG(1)*EPKE*(F1+F17+F5+F7+F9+F11+F13+F15)*F3
COVAL(MIX,F3, FIXFLU,GRND)
       *  Fluid 4
   VAL=2.*RG(1)*EPKE*(F3*F5+F2*F6+F1*F7)-$
           RG(1)*EPKE*(F2+F6+F8+F10+F12+F14+F16)*F4
COVAL(MIX,F4, FIXFLU,GRND)
       *  Fluid 5
   VAL=2.*RG(1)*EPKE*(F4*F6+F3*F7+F2*F8+F1*F9)-$
           RG(1)*EPKE*(F1+F3+F17+F7+F9+F11+F13+F15)*F5
COVAL(MIX,F5, FIXFLU,GRND)
       *  Fluid 6
   VAL=2.*RG(1)*EPKE*(F5*F7+F4*F8+F3*F9+$
                                 F2*F10+F1*F11)-$
           RG(1)*EPKE*(F2+F4+F8+F10+F12+F14+F16)*F6
COVAL(MIX,F6, FIXFLU,GRND)
        * Fluid 7
   VAL=2.*RG(1)*EPKE*(F6*F8+F5*F9+F4*F10+$
                           F3*F11+F2*F12+F1*F13)-$
           RG(1)*EPKE*(F1+F3+F5+F17+F9+F11+F13+F15)*F7
COVAL(MIX,F7, FIXFLU,GRND)
         * Fluid 8
   VAL=2.*RG(1)*EPKE*(F7*F9+F6*F10+F5*F11+$
                     F4*F12+F3*F13+F2*F14+F1*F15)-$
           RG(1)*EPKE*(F2+F4+F6+F10+F12+F14+F16)*F8
COVAL(MIX,F8, FIXFLU,GRND)
         * Fluid 9
   VAL=2.*RG(1)*EPKE*(F8*F10+F7*F11+F6*F12+$
               F5*F13+F4*F14+F3*F15+F2*F16+F1*F17)-$
           RG(1)*EPKE*(F1+F3+F5+F7+F17+F11+F13+F15)*F9
COVAL(MIX,F9, FIXFLU,GRND)
         * Fluid 10
   VAL=2.*RG(1)*EPKE*(F9*F11+F8*F12+F7*F13+$
                      F6*F14+F5*F15+F4*F16+F3*F17)-$
           RG(1)*EPKE*(F2+F4+F6+F8+F12+F14+F16)*F10
COVAL(MIX,F10, FIXFLU,GRND)
         * Fluid 11
   VAL=2.*RG(1)*EPKE*(F10*F12+F9*F13+F8*F14+$
                              F7*F15+F6*F16+F5*F17)-$
           RG(1)*EPKE*(F1+F3+F5+F7+F9+F17+F13+F15)*F11
COVAL(MIX,F11, FIXFLU,GRND)
         * Fluid 12
   VAL=2.*RG(1)*EPKE*(F11*F13+F10*F14+F9*F15+$
                                      F8*F16+F7*F17)-$
           RG(1)*EPKE*(F2+F4+F6+F8+F10+F14+F16)*F12
COVAL(MIX,F12, FIXFLU,GRND)
         * Fluid 13
   VAL=2.*RG(1)*EPKE*(F12*F14+F11*F15+$
                              F10*F16+F9*F17)-$
           RG(1)*EPKE*(F1+F3+F5+F7+F9+F11+F15+F17)*F13
COVAL(MIX,F13, FIXFLU,GRND)
         * Fluid 14
   VAL=2.*RG(1)*EPKE*(F13*F15+F12*F16+F11*F17)-$
           RG(1)*EPKE*(F2+F4+F6+F8+F10+F12+F16)*F14
COVAL(MIX,F14, FIXFLU,GRND)
         * Fluid 15
   VAL=2.*RG(1)*EPKE*(F14*F16+F13*F17)-$
           RG(1)*EPKE*(F1+F3+F5+F7+F9+F11+F13+F17)*F15
COVAL(MIX,F15, FIXFLU,GRND)
         * Fluid 16
   VAL=2.*RG(1)*EPKE*(F15*F17)-$
           RG(1)*EPKE*(F2+F4+F6+F8+F10+F12+F14)*F16
COVAL(MIX,F16, FIXFLU,GRND)
         * Fluid 17
   CO=RG(1)*EPKE*(F1+F3+F5+F7+F9+F11+F13+F15)
COVAL(MIX,F17  , GRND  ,0.0 )
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
    The above   source/sink   terms   in   the  fluid-mass-fraction
    equations are shared according to a  coupling/splitting  scheme
    derived from Spalding concept.

    The scheme hypotheses is that the coupling may only take place
    between those parent fluids which would produce the appropriate
    offsprings inheriting the ATTRIBUTES of either parent in EQUAL
    proportion.
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
  ** Output calculations
STORE(CAV,MAS,GAV,GF)
FIINIT(GF)=0.0
   CAV=16./16.*F1 + 15./16.*F2 + 14./16.*F3 +$
               13./16.*F4 + 12./16.*F5 + 11./16.*F6 +$
               10./16.*F7 +  9./16.*F8 +  8./16.*F9 +$
                7./16.*F10+  6./16.*F11+  5./16.*F12+$
                4./16.*F13+  3./16.*F14+  2./16.*F15+$
                1./16.*F16+  0./16.*F17
  IF(ISWEEP.EQ.LSWEEP)
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
   At the end of the iz-slab for the last sweep, CAV, averaged
    concentration of central jet fluid, is calculated from the
    individual fluid mass-fractionsand their arrributes;
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
   MAS=F1+F2+F3+F4+F5+F6+F7+F8+F9+F10+F11+F12+$
               F13+F14+F15+F16+F17
  IF(ISWEEP.EQ.LSWEEP)
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
   At the end of the iz-slab for the last sweep, MAS, sum of fluid mass
   fractions, is calculated to check its equality to unity;
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
   GAV=ABS(CAV-16./16)*F1 + ABS(CAV-15./16.)*F2 +$
               ABS(CAV-14./16)*F3 + ABS(CAV-13./16.)*F4 +$
               ABS(CAV-12./16)*F5 + ABS(CAV-11./16.)*F6 +$
               ABS(CAV-10./16)*F7 + ABS(CAV- 9./16.)*F8 +$
               ABS(CAV- 8./16)*F9 + ABS(CAV- 7./16.)*F10+$
               ABS(CAV- 6./16)*F11+ ABS(CAV- 5./16.)*F12+$
               ABS(CAV- 4./16)*F13+ ABS(CAV- 3./16.)*F14+$
               ABS(CAV- 2./16)*F15+ ABS(CAV- 1./16.)*F16+$
               ABS(CAV- 0./16)*F17
  IF(ISWEEP.EQ.LSWEEP)
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
   At the end of the iz-slab for the last sweep, GAV, averaged
   concentration fluctuation, is calculated as the sum of local
   deviations of averaged concentrations from the individual
   concentration attributes;
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
   GF =GAV/(CAV+TINY)
  IF(ISWEEP.EQ.LSWEEP)
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
   At the end of the iz-slab for the last sweep, GF, averaged
   concentration fluctuation normalised by local averaged concentration
   of central jet fluid, is calculated.
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
    ** Output data processing for plotting PDF
     * Specify the cell in question: IY=IG(1), IZ=IG(2)
IG(1)=4; IG(2)=4
STORE(FPD);FIINIT(FPD)=0.0
   FPD=F1[1,IG(1),IG(2)]*AMAX1(ABS(F1[1,IG(1),IG(2)]$
                     -YV2D)/(F1[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,1,1) /ISWEEP.EQ.LSWEEP
   FPD=F2[1,IG(1),IG(2)]*AMAX1(ABS(F2[1,IG(1),IG(2)]$
                     -YV2D)/(F2[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,2,2) /ISWEEP.EQ.LSWEEP
   FPD=F3[1,IG(1),IG(2)]*AMAX1(ABS(F3[1,IG(1),IG(2)]$
                     -YV2D)/(F3[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,3,3) /ISWEEP.EQ.LSWEEP
   FPD=F4[1,IG(1),IG(2)]*AMAX1(ABS(F4[1,IG(1),IG(2)]$
                     -YV2D)/(F4[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,4,4) /ISWEEP.EQ.LSWEEP
   FPD=F5[1,IG(1),IG(2)]*AMAX1(ABS(F5[1,IG(1),IG(2)]$
                     -YV2D)/(F5[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,5,5) /ISWEEP.EQ.LSWEEP
   FPD=F6[1,IG(1),IG(2)]*AMAX1(ABS(F6[1,IG(1),IG(2)]$
                     -YV2D)/(F6[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,6,6) /ISWEEP.EQ.LSWEEP
   FPD=F7[1,IG(1),IG(2)]*AMAX1(ABS(F7[1,IG(1),IG(2)]$
                     -YV2D)/(F7[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,7,7) /ISWEEP.EQ.LSWEEP
   FPD=F8[1,IG(1),IG(2)]*AMAX1(ABS(F8[1,IG(1),IG(2)]$
                     -YV2D)/(F8[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,8,8) /ISWEEP.EQ.LSWEEP
   FPD=F9[1,IG(1),IG(2)]*AMAX1(ABS(F9[1,IG(1),IG(2)]$
                     -YV2D)/(F9[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,9,9) /ISWEEP.EQ.LSWEEP
   FPD=F10[1,IG(1),IG(2)]*AMAX1(ABS(F10[1,IG(1),IG(2)]$
                      -YV2D)/(F10[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,10,10) /ISWEEP.EQ.LSWEEP
   FPD=F11[1,IG(1),IG(2)]*AMAX1(ABS(F11[1,IG(1),IG(2)]$
                      -YV2D)/(F11[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,11,11) /ISWEEP.EQ.LSWEEP
   FPD=F12[1,IG(1),IG(2)]*AMAX1(ABS(F12[1,IG(1),IG(2)]$
                      -YV2D)/(F12[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,12,12) /ISWEEP.EQ.LSWEEP
   FPD=F13[1,IG(1),IG(2)]*AMAX1(ABS(F13[1,IG(1),IG(2)]$
                      -YV2D)/(F13[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,13,13) /ISWEEP.EQ.LSWEEP
   FPD=F14[1,IG(1),IG(2)]*AMAX1(ABS(F14[1,IG(1),IG(2)]$
                      -YV2D)/(F14[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,14,14) /ISWEEP.EQ.LSWEEP
   FPD=F15[1,IG(1),IG(2)]*AMAX1(ABS(F15[1,IG(1),IG(2)]$
                      -YV2D)/(F15[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,15,15) /ISWEEP.EQ.LSWEEP
   FPD=F16[1,IG(1),IG(2)]*AMAX1(ABS(F16[1,IG(1),IG(2)]$
                      -YV2D)/(F16[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,16,16) /ISWEEP.EQ.LSWEEP
   FPD=F17[1,IG(1),IG(2)]*AMAX1(ABS(F17[1,IG(1),IG(2)]$
                      -YV2D)/(F17[1,IG(1),IG(2)]+0.-YV2D) ,0.0)
   REGION(1,1,1,NY,17,17) /ISWEEP.EQ.LSWEEP
  >>>>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>
   The above operations are made at the end of the iz-slab for the last
   sweep to fill each IY-column of the domain by Fi value.
  <<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<

    PHOTON USE
    P



    msg(  Single fluid central-jet concentration contours
    con h1 x 1 fil;.001
    pause
    con cl; red
    msg(  Averaged 17-fluid contours
    con cav x 1 fil;.001
    pause
    con cl; red
    msg( Concentration fluctuations by transport equation
    con gg  x 1 fil;.001
    pause
    con cl; red
    msg( Averaged 17-fluid concentration fluctuation
    con gav x 1 fil;.001
    msg
    msg Hit Enter for FPD hystogram
    msg
    pause
     p

     20 1

     con fpd x 1 fil;10.
     msg Hit Enter to continue
    ENDUSE
STOP