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 Press RETURN to continue readvdu(nphi,int,nphi) 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 STEADY=F;GRDPWR(T,20,40.0,1.0) 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.. READVDU(NX,INT,NX) 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.. READVDU(NY,INT,NY) 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 READVDU(RAYLEIGH,REAL,RAYLEIGH) 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