```
photon use
p
parphi
5 5 1

view 1 1 1
msg the grid
msg press RETURN for temperature contours at successive times
gr z 1;pause
do kk=2,20,2
con temp z kk fi;0.01
enddo
msg press RETURN for velocity vectors at successive times
pause;con off;red
do kk=2,20,2
vec z kk sh
enddo
msg press RETURN for C1 contours at successive times
pause
do kk=2,20,2
con c1 z kk fi;0.01
enddo
enduse

GROUP 1. Run title and other preliminaries
TEXT(Transient Free Convection In A Box
TITLE
DISPLAY
A stationary fluid in a cavity is suddenly exposed to walls of
unequal temperature. Flow is generated by buoyancy forces
acting oppositely on hot and cold parts of the fluid.

----------------------------
/|        -------->         |/
/|                          |/
/|           --->           |/        |
Hot     /| ^                        |/ Cold   |
wall;   /| |    ^                 | |/ wall;  |g
T=1.0   /| |    |            |    | |/ T=-1.0 |
/| |         *--     v    | |/        |
/|              \         v |/        v
/|          <--- \          |/
/|                \         |/
/|        <--------\        |/
-------------------\--------
^                        \
y|                        P1=0.0 at PATCH REFP
|--->
x
The temperature conditions are applied only after the first time
step so that the correct pressure field can first be established.

Under-relaxation is applied to V1 for this time step only, so as
to facilitate quick establishment of the correct pressure, by
means of the patch called STEP1, with a large CO and with VAL set
as SAME.

The 10 * 10 grid is too coarse for accuracy, but may easily be
refined.

Other interesting variants include changes to:
* the aspect ratio of the cavity;
* the Prandtl number, here taken as 5.0; and
* the Rayleigh number, here taken as about 1.E6.

ENDDIS
GROUP 2. Transience; time-step specification

GROUP 3. X-direction grid specification
** Set a symmetrical grid, consisting of power-law
grids (varying as IX**2.0), which start from
each edge and meet in the middle.
NX=10
mesg(NX=:NX: If not OK, insert desired value..
GRDPWR(X,-NX,1.0,2.0)

GROUP 4. Y-direction grid specification
** Set a symmetrical grid as in GROUP 3
NY=10
mesg(NY=:NY: If not OK, insert desired value..
GRDPWR(Y,-NY,1.0,2.0)
GROUP 7. Variables stored, solved & named
SOLVE(P1,U1,V1,H1)

GROUP 8. Terms (in differential equations) & devices
TERMS(H1,N,Y,Y,Y,Y,Y)

GROUP 9. Properties of the medium (or media)
** Set the temperature as TMP1A+TMP1B*H1
ENUL=1.e-3;TMP1=LINH;TMP1A=0.0;TMP1B=1.0;cp1=1./tmp1b
mesga(The dimensionless time step DT*ENUL/YVLAST**2 equals
mesgb(          :tlast*enul/(yvlast**2*lstep):
** Set the density as RHO1A+RHO1B*Temperature
RHO1=LINTEMP;RHO1A=1.0;RHO1B=-0.01;PRNDTL(H1)=5.0
REAL(RAYLEIGH,RAYL1,AGRAV)
AGRAV=9.81
RAYLEIGH=AGRAV*(-RHO1B)*2.0*YVLAST*YVLAST*YVLAST*PRNDTL(H1)
RAYLEIGH=RAYLEIGH/(ENUL*ENUL)
RAYL1=RAYLEIGH
mesg(Rayleigh Number = :rayleigh:  If not OK, insert desired number
AGRAV=AGRAV*RAYLEIGH/RAYL1
RAYLEIGH
AGRAV
GROUP 11. Initialization of variable or porosity fields
FIINIT(P1)=0.0;FIINIT(U1)=0.0;FIINIT(V1)=0.0;FIINIT(H1)=0.0

GROUP 13. Boundary conditions and special sources
** Pressure relief
PATCH(REFP,CELL,NX/2,NX/2,NY/2,NY/2,1,1,1,LSTEP)
COVAL(REFP,P1,FIXP,0.0)
COVAL(REFP,U1,ONLYMS,0.0);COVAL(REFP,V1,ONLYMS,0.0)
** Set the heat flux to be the prevailing value of H1
in the cell
COVAL(REFP,H1,FIXVAL,SAME)
** Hot wall
WALL (HOTWALL,WEST,1,1,1,NY,1,1,2,LSTEP)
COVAL(HOTWALL,V1,1.0,0.0);COVAL(HOTWALL,H1,1.0,1.0)
** Buoyancy
PATCH(BUOYANCY,PHASEM,1,NX,1,NY,1,1,1,LSTEP)
COVAL(BUOYANCY,V1,FIXFLU,-AGRAV)
** Cold wall
WALL (COLDWALL,EAST,NX,NX,1,NY,1,1,2,LSTEP)
COVAL(COLDWALL,V1,1.0,0.0);COVAL(COLDWALL,H1,1.0,-1.0)
** Top wall
WALL (TOPWALL,NORTH,1,NX,NY,NY,1,1,2,LSTEP)
COVAL(TOPWALL,U1,1.0,0.0)
** Bottom wall
WALL (BOTMWALL,SOUTH,1,NX,1,1,1,1,2,LSTEP)
COVAL(BOTMWALL,U1,1.0,0.0)
** temporary under-relaxation for V1 during first time step
PATCH(STEP1,CELL,1,NX,1,NY,1,1,1,1)
COVAL(STEP1,V1,1.E5,SAME)

GROUP 15. Termination of sweeps
LSWEEP=50;SELREF=T;RESFAC=0.01
GROUP 21. Print-out of variables
OUTPUT(P1,Y,Y,Y,Y,Y,Y);OUTPUT(U1,Y,Y,Y,Y,Y,Y)
OUTPUT(V1,Y,Y,Y,Y,Y,Y);OUTPUT(H1,Y,Y,Y,Y,Y,Y)

GROUP 22. Spot-value print-out
IXMON=1;IYMON=9
SPEDAT(SET,GXMONI,TRANSIENT,L,F)
GROUP 23. Field print-out and plot control
NTPRIN=5;NXPRIN=3;NYPRIN=2
PATCH(CONT,CONTUR,1,NX,1,NY,1,1,1,LSTEP)
PLOT(CONT,H1,0.0,10.0)
PATCH(IYEQ5,PROFIL,1,NX,NY/2,NY/2,1,1,1,LSTEP)
PLOT(IYEQ5,V1,-0.1,0.1);PLOT(IYEQ5,H1,0.0,1.0)
PATCH(IYEQ9,PROFIL,1,1,NY-1,NY-1,1,1,1,LSTEP)
PLOT(IYEQ9,H1,0.0,1.0)
IDISPA=1;TSTSWP=-1

The following statements activate coding in Group 19 of
GREX3,HTM :

# the number of commands
SPEDAT(PRINT,NUMBER,I,5)

# giving h1 the necessary four-character name
NAME(H1)= TEMP

# storing density because it is needed for computation of totals
STORE(DEN1)
# the command to compute the total of temperature
SPEDAT(PRINT,COMMAND1,C,TOTAL_TEMP)

# storing c1 as a place to put
(temp - tempmin)/(tempmax - tempmin)**power
STORE(C1)

# the command to compute and store the above quantity
SPEDAT(PRINT,COMMAND2,C,POW_TEMP_C1)

# the power to be used
SPEDAT(PRINT,POWER,R,0.5)

# the command to print the minimum and maximum values of u1
SPEDAT(PRINT,COMMAND3,C,MINMAX_U1)

# the command to print the average value of c1
SPEDAT(PRINT,COMMAND4,C,AVE_C1)

# the commands to print the value of v1 prevailing at
ix=nx/2; iy=ny/2; iz=1
SPEDAT(PRINT,COMMAND5,C,VALUE_V1)
SPEDAT(PRINT,IXLOC,I,:NX/2:)
SPEDAT(PRINT,IYLOC,I,:NY/2:)
SPEDAT(PRINT,IZLOC,I,1)

LIBREF  =     340
STOP
```