photon use
  p
 
 
 
  msg Temperature contours
  con temp z 1 fi;0.1
  msg
  msg Press return for velocity vectors
  pause;con off; red
  vec z 1 sh
  msg Press return for pressure contours
  pause;vec off; red
  msg
  msg Pressure contours
  con p1 z 1 fi;0.1
  msg Press return for streamlines
  pause;con off;red;stre 2d z 1;int 10
  msg Press e to END
  enduse
 
    GROUP 1. Run title and other preliminaries
TEXT(Free Convection In A Porous Medium
TITLE
  DISPLAY
  ******************************************************************
  *  A horizontal cold pipe, eg containing liquefied natural gas,  *
  *  is insulated by a cylindrical blanket of glass wool, through  *
  *  which air can circulate slowly by natural convection, as a    *
  *  result of which the heat transfer rate is considerably larger *
  *  than if conduction were the only mechanism.                   *
  *                                                                *
  *  The governing dimensionless parameters are:                   *
  *1.The Prandtl number, in which it is the effective conductivity *
  *  of the glass wool which should enter,                         *
  *2.the Rayleigh number, ie                                       *
  *             grav * DVO1DT * deltaT * diam**3/vis**2,           *
  *3.and the ratio of average air-passage width in the glass wool  *
  *  to the pipe diameter.                                         *
  *                                                                *
  *  The values chosen here are Pr = 0.5, Ra = 1.e8, Wi = 1.e-4    *
  *                                                                *
  ******************************************************************
  ENDDIS
INTEGER(IDUM)
NX=20;NY=20
REAL(PR,RAYL,WI,RESCON,DUMMY,POWERX,POWERY)
PR=0.5;RAYL=1.E8;WI=1.E-4
 
    GROUP 3. X-direction grid specification
CARTES=F;GRDPWR(X,-NX,3.14159,POWERX)
    GROUP 4. Y-direction grid specification
RINNER=1.0;GRDPWR(Y,-NY,1.0,POWERY)
    GROUP 5. Z-direction grid specification
ZWLAST=1.1054
    GROUP 7. Variables stored, solved & named
SOLUTN(P1,Y,Y,N,N,N,Y);SOLUTN(U1,Y,Y,N,Y,N,Y)
SOLUTN(V1,Y,Y,N,Y,N,Y);SOLUTN(H1,Y,Y,N,N,N,Y)
NAME(H1)=TEMP
    GROUP 8. Terms (in differential equations) & devices
    Suppress enthalpy source term
TERMS(TEMP,N,Y,Y,N,Y,N)
    GROUP 9. Properties of the medium (or media)
ENUL=1.0;PR=5.0;PRNDTL(TEMP)=PR
    GROUP 11. Initialization of variable or porosity fields
FIINIT(TEMP)=0.4;FIINIT(U1)=0.0;FIINIT(V1)=0.0;FIINIT(P1)=0.0
    GROUP 13. Boundary conditions and special sources
    ...south boundary
WALL (INNERWAL,SOUTH,1,NX,1,1,1,1,1,1)
COVAL(INNERWAL,U1,1.0,0.0);COVAL(INNERWAL,TEMP,1.0,0.0)
    ...north boundary
WALL (OUTERWAL,NORTH,1,NX,NY,NY,1,1,1,1)
COVAL(OUTERWAL,U1,1.0,0.0);COVAL(OUTERWAL,TEMP,1.0,1.0)
    ...pressure relief
PATCH(PRELIEF,CELL,NX/2,NX/2,1,1,1,1,1,1)
COVAL(PRELIEF,P1,FIXP,0.0);COVAL(PRELIEF,TEMP,ONLYMS,SAME)
 
mesga( Rayleigh number is :RAYL:. OK? If not, insert new value
readvdu(RAYL,real,rayl)
 
do ii=1,5
+ mesg(
enddo
mesga( Average porous medium pore size/pipe diam = :WI:. OK?
mesga(  If not, insert new value.
readvdu(WI,real,WI)
 
    ...buoyancy
PATCH(BUOYANCY,VOLUME,1,NX,1,NY,1,1,1,1)
BUOYB=-1.0;BUOYE=0.0;BUOYD=RAYL/PR
COVAL(BUOYANCY,U1,FIXFLU,BOUSS);COVAL(BUOYANCY,V1,FIXFLU,BOUSS)
    ...resistance exerted by porous medium
RESCON=ENUL/(RINNER*2.0*WI)
PATCH(MEDIUM,VOLUME,1,NX,1,NY,1,1,1,1)
COVAL(MEDIUM,U1,RESCON,0.0);COVAL(MEDIUM,V1,RESCON,0.0)
    GROUP 15. Termination of sweeps
LSWEEP=200;SELREF=T;RESFAC=0.1
    GROUP 16. Termination of iterations
LITER(P1)=100;LITER(TEMP)=10;liter(u1)=10;liter(v1)=10
    GROUP 17. Under-relaxation devices
DUMMY=10*(RAYL**0.5+1.0)
RELAX(U1,FALSDT,1.0/DUMMY);RELAX(V1,FALSDT,1.0/DUMMY)
    GROUP 21. Print-out of variables
OUTPUT(P1,N,Y,Y,Y,Y,Y);OUTPUT(U1,N,Y,Y,Y,Y,Y)
OUTPUT(V1,N,Y,Y,Y,Y,Y);OUTPUT(TEMP,N,Y,Y,Y,Y,Y)
    GROUP 22. Spot-value print-out
TSTSWP=-1;IXMON=10;IYMON=5;
    GROUP 23. Field print-out and plot control
NXPRIN=5;NYPRIN=3;IPLTL=LSWEEP;IPROF=1
PATCH(TX1,PROFIL,1,1,1,NY,1,1,1,1);PLOT(TX1,TEMP,0.0,0.0)
PATCH(TNX,PROFIL,NX,NX,1,NY,1,1,1,1);PLOT(TNX,TEMP,0.0,0.0)
PATCH(TNX2,PROFIL,NX/2,NX/2,1,NY,1,1,1,1)
PLOT(TNX2,TEMP,0.0,0.0);PLOT(TNX2,U1,0.0,0.0)
PATCH(TCON,CONTUR,1,30,1,30,1,1,1,1);PLOT(TCON,TEMP,0.0,20.0)
    GROUP 24. Dumps for restarts