PHOTON USE
  p
 
 
 
 
  msg the grid. press return for temperature contours
  gr z 1;pause; gr off;red;gr ou z 1
  msg temperature contours. press return for velocity vectors
  con temp z 1 fi;0.0002;pause; con off;red;se re ve 2;vec z 1 sh
  msg velocity vectors. press return for reduced-pressure contours
  pause;vec off;red;con p1 z 1 fi;0.0002
  ENDUSE
 
 
  GROUP 1. RUN TITLE AND OTHER PRELIMINARIES
TEXT(Laminar Free Convection In Cavity 
TITLE
  DISPLAY
  A two-dimensional square cavity is formed between two
  vertical walls, one of which is heated and the other cooled.
  the top and bottom of the cavity are bounded by walls at which
  there is friction but no heat transfer.
 
  The left wall is the hot wall of a temperature  Thot=1 K
  The right wall is the cold wall of a temperature Tcold=0 K
 
  ENDDIS
#pause
  DISPLAY
  Parameters:
  -----------
 
  1) The flow characteristics are set by way of
     Length of the cavity=CAVL
     Prandtl Number=Pr=enul/diffusity
     Rayleigh Number=Ra=Beta*g*(Th-Tc)/(enul*alfa)
 
     The Thermal Expansion Coefficient Beta as well as
     the Thermal diffusivity Alfa remain constant.
 
  2) This library case also enables the user to select
     the numerical solution procedure, namely:
    -the choice between two solvers:
     STONE or CNGR (Conjugate-residuals gradient)
    -the possibilty of the activation of EXPERT
     the Self-Adjusting-Relaxation Algorithm.
 
  ENDDIS
#pause
 
                   SPECIAL DATA AND INPUT
                   ======================
    RA        THE RAYLEIGH NUMBER
    DVO1DT    THE COEFFICIENT OF THERMAL EXPANSION      1/K
    AGRAV     GRAVITY                                   M/S^2
    TREF      REFERENCE TEMPERATURE                     C DEG.
    CAVL      THE LENGTH OF THE CAVITY                  M
 
REAL(TREF,AGRAV,CAVL,Alfa,Pr,RA)
REAL(ansca,anspr,ansra); CHAR(ansol,ansex)
    *Beta=0.01
DVO1DT = 0.1E-1
 
    *Alfa= 1.0
Alfa =1.0
 
    *tref
TREF=0.5;                              BUOYE= TREF
 
MESG(Please define the cavity length(default=1.)
READVDU(ansca,real,1.)
CAVL=ansca
 
 
MESG(Please define Rayleigh Number:(default=1.E3)
READVDU(ansra,real,1.E3)
RA=ansra
 
 
MESG(Please define Prandtl Number:(default=0.71)
READVDU(anspr,real,0.71)
Pr=anspr
Pr
 
    *kinematic viscosity
ENUL=Alfa*Pr
ENUL
 
    *Gravity
AGRAV=RA*ENUL*Alfa/(DVO1DT*CAVL);                BUOYB=-AGRAV
AGRAV
 
    *choice of the solver
MESG(Please choose a solver Stone or CNGR:(default=CNGR)
READVDU(ansol,char,c)
CASE :ansol: OF
 WHEN CNGR,1
 CSG3=CNGR
 MESG(CNGR solver activated
 ORELSE
 MESG(STONE  solver activated
ENDCASE
 
    *select EXPERT for relaxations
MESG(Select EXPERT ?:(default=y)
READVDU(ansex,char,y)
CASE :ansex: OF
 WHEN YES,1
 NOWIPE=T
    *load coefficients required for the authomatic relaxations
 L(030)
 EXPERT=T
 MESG(EXPERT Self-Adjustment Algorithm activated
 ORELSE
 MESG(EXPERT not activated
ENDCASE
 
 
    GROUP 2. Transience; time-step specification
    GROUP 3. X-direction grid specification
    GROUP 4. Y-DIRECTION GRID SPECIFICATION
 
NX=26; XULAST=CAVL; NY=26; YVLAST=CAVL
l(unigrid
 
    GROUP 6. Body-fitted coordinates or grid distortion
    GROUP 7. VARIABLES STORED, SOLVED & NAMED
 
   *** whole-field solver for p1 is activated.
SOLVE(P1,U1,V1,H1); SOLUTN(P1,Y,Y,Y,N,N,N); NAME(H1)=TEMP
SOLUTN(U1,Y,Y,Y,N,N,N); SOLUTN(V1,Y,Y,Y,N,N,N)
 
    GROUP 8. TERMS (IN DIFFERENTIAL EQUATIONS) & DEVICES
   *** deactivate the built-in source in temp equation.
TERMS(TEMP,N,Y,Y,Y,Y,Y)
 
    GROUP 9. PROPERTIES OF THE MEDIUM (OR MEDIA)
RHO1=1.207; PRNDTL(TEMP)=Pr
CP1=1.
 
    GROUP 13. BOUNDARY CONDITIONS AND SPECIAL SOURCES
 
   1. HOT WALL BOUNDARY: CONSTANT TEMPERATURE OF 1 DEG.
WALL (HOT,WEST,1,1,1,NY,1,NZ,1,1)
COVAL(HOT,V1,1.0,0.0)
COVAL(HOT,TEMP,1.0,1.0)
 
   2. COLD WALL BOUNDARY: CONSTANT TEMPERATURE OF 0 DEG.
WALL (COLD,EAST,NX,NX,1,NY,1,NZ,1,1)
COVAL(COLD,V1,1.0,0.0)
COVAL(COLD,TEMP,1.0,0.0)
 
   3. LOW WALL BOUNDARY: ADIABATIC BUT WITH FRICTION
WALL (LOWAL,SOUTH,1,NX,1,1,1,NZ,1,1)
COVAL(LOWAL,U1,1.0,0.0)
 
   4. HIGH WALL BOUNDARY: ADIABATIC BUT WITH FRICTION
WALL (HIWAL,NORTH,1,NX,NY,NY,1,NZ,1,1)
COVAL(HIWAL,U1,1.0,0.0)
 
   5. BUOYANCY FORCE
   *** buoyancy source for the y-direction momentum equation
       is set in grex3, group 13, sec.15, and gxbuoy;
       val=grnd3 activates the boussinesq approximation as follows:
       rho*volume*grav.*DVO1DT*(tref-t).
       the following data need to be set for this purpose:
       buoye=tref; buoyb=agrav., as stated in
       group 1 above.
 
PATCH(BUOY,PHASEM,1,NX,1,NY,1,NZ,1,1); COVAL(BUOY,V1,FIXFLU,GRND3)
 
   6. REFERENCE PRESSURE AT THE CENTRE OF THE CAVITY
 
PATCH(REFP,CELL,NX/2,NX/2,NY/2,NY/2,1,1,1,1)
COVAL(REFP,P1,FIXP,0.0); COVAL(REFP,U1,ONLYMS,0.0)
COVAL(REFP,V1,ONLYMS,0.0)
COVAL(REFP,TEMP,ONLYMS,SAME)
 
    GROUP 14. Downstream pressure for PARAB=.TRUE.
    GROUP 15. TERMINATION OF SWEEPS
LSWEEP=200; SELREF=T; RESFAC=0.001
ENDIT(P1)=1.E-6; ENDIT(V1)=1.E-6; ENDIT(U1)=1.E-6
ENDIT(TEMP)=1.E-6
 
    GROUP 16. Termination of iterations
LITER(U1)=2; LITER(V1)=2
 
    GROUP 17. Under-relaxation devices
RELAX(U1,FALSDT,1.E-03); RELAX(V1,FALSDT,1.E-03)
RELAX(TEMP,FALSDT,1.)
 
    GROUP 18. Limits on variables or increments to them
    GROUP 19. Data communicated by satellite to GROUND
    GROUP 20. Preliminary print-out
    GROUP 21. Print-out of variables
    GROUP 22. Spot-value print-out
IXMON=nx/4; IYMON=ny/4; NYPRIN=NY/5; NZPRIN=NZ/5; NPLT=1
TSTSWP=-1
 
    GROUP 23. Field print-out and plot control
  *** Temperature and velocity profiles
PATCH(PROF,PROFIL,1,NX,1,NY,NZ/2,NZ/2,1,1)
PLOT (PROF,W1,0.0,0.0); PLOT (PROF,TEMP,0.0,1.0)
 
  *** Temperature contours
PATCH(CONT,CONTUR,1,NX,1,NY,1,NZ,1,1)
PLOT (CONT,TEMP,0.0,10.0)