PHOTON USE p phi msg EGWF TESTS IN 2D DUCT msg msg Pressure distribution: cl;red con p1 z max fi;.001 vec z max msg msg Pressto continue pause msg msg Temperature distribution: cl;red con tem1 z max fi;.001 vec z max msg msg Press to continue pause msg msg Concentration distribution: cl;red con c1 z max fi;.001 vec z max msg msg Press to continue pause msg msg Press e to END enduse CHAR(CUSEGW,CTURB,CBFC);BOOLEAN(TURBL) INTEGER(NXF,NXL,NYF,NYL,NZF,NZL,ISUM,N1) REAL(VIN,VFRAC,TIN,CIN) VIN=10.;TIN=100.;CIN=-1. MESG( MESG( TEST OF EGWF FOR FLOW AROUND A HEATED SOLUBLE CUBE MESG( MESG( The flow enters the W-S-L corner of the domain, MESG( impinges on the solid cube, which is heated and MESG( which has constant C1 concentration, and then MESG( leaves the domain through the E-N-H corner. MESG( MESG( To run the EGWF variant, enter T MESG( To run the SOWF variant, enter F MESG( READVDU(CUSEGW,CHAR,F) MESG( MESG( Enter the turbulence model option; default MESG( a LAMinar flow. MESG( The options are: MESG( LAM - laminar flow MESG( TURB - simple mixing-length turbulence MESG( LVEL - LVEL algebraic turbulence model MESG( KLMODL - the k-l model of turbulence MESG( KEMODL - the k-e model of turbulence MESG( MESG( For further info. please see the PHOENICS MESG( encyclopaedia MESG( READVDU(CTURB,CHAR,LAM) MESG( MESG( Enter the COefficient for the WALL-functions; MESG( the default is GRND2 (LOGLAW). MESG( The options are: MESG( 1/PRNDTL() - suitable for laminar flow MESG( GRND1 (BLASIUS) - suitable for laminar flow MESG( GRND2 (LOGLAW) - equilibrium log-law wall- MESG( function for turbulent flows MESG( GRND3 (GENLAW) - non-equilibrium log-law wall- MESG( function for turbulent flows MESG( (not implemented for EGWF=T) MESG( READVDU(WALLCO,REAL,GRND2) MESG( MESG( Flow is always in the Y direction: MESG( by default the domain has 7x7x7 cells MESG( Enter the new NX, NY & NZ; MESG( if(iqalib.ne.0) then + nx=7; ny=7; nz=7 endif READVDU(NX,INT,7) READVDU(NY,INT,7) READVDU(NZ,INT,7) MESG( MESG( Enter the frequency of block-corrections. MESG( The default is no block corrections. MESG( READVDU(ISOLBK,INT,0) MESG( MESG( Enter the thickness of the wall-blocks. MESG( The default is 1. MESG( READVDU(N1,INT,1) EGWF=:CUSEGW: IF (EGWF) THEN + CASE :CTURB: OF + WHEN TURB,4 + TEXT( TURBULENT TEST OF EGWF - EGWF VARIANT) + WHEN LAM,3 + TEXT( LAMINAR TEST OF EGWF - EGWF VARIANT) + ORELSE + TEXT( :CTURB: TEST OF EGWF - EGWF VARIANT) + TURBL = T + ENDCASE ELSE + CASE :CTURB: OF + WHEN TURB,4 + TEXT( TURBULENT TEST OF EGWF - SOWF VARIANT) + WHEN LAM,3 + TEXT( LAMINAR TEST OF EGWF - SOWF VARIANT) + ORELSE + TEXT( :CTURB: TEST OF EGWF - SOWF VARIANT) + TURBL = T + ENDCASE ENDIF IF (NX .LE. 0) THEN + NX=1 ENDIF IF (NY .LE. 0) THEN + NY=1 ENDIF IF (NZ .LE. 0) THEN + NZ=1 ENDIF IF (N1 .LE. 0) THEN + N1=1 ENDIF GROUP 2. Transience; time-step specification GROUP 3. X-direction grid specification GRDPWR(X,-NX,1.0,2.0) NXF=1+N1;NXL=NX-N1 GROUP 4. Y-direction grid specification GRDPWR(Y,NY,10.0,1.0) GROUP 5. Z-direction grid specification GRDPWR(Z,-NZ,1.0,2.0) NZF=1+N1;NZL=NZ-N1 GROUP 6. Body-fitted coordinates or grid distortion GROUP 7. Variables stored, solved & named SOLVE(P1);SOLUTN(P1,Y,Y,Y,P,P,Y) SOLVE(TEM1);SOLUTN(TEM1,Y,Y,Y,P,P,Y) SOLUTN(C1,Y,Y,Y,P,P,Y) ISUM=0 IF (NX .GT. 1) THEN + SOLVE(U1);ISUM=ISUM+1 + SOLUTN(U1,P,P,Y,P,P,Y) ENDIF IF (NY .GT. 1) THEN + SOLVE(V1);ISUM=ISUM+1 + SOLUTN(V1,P,P,Y,P,P,Y) ENDIF IF (NZ .GT. 1) THEN + SOLVE(W1);ISUM=ISUM+1 + SOLUTN(W1,P,P,Y,P,P,Y) ENDIF VFRAC=ISUM**0.5 STORE(PRPS) IF (ISOLBK .NE. 0) THEN + STORE(BLOK) + IVARBK=-1 ENDIF IF (TURBL) THEN + TURMOD(:CTURB:) + SOLUTN(KE,P,P,p,P,P,Y) + SOLUTN(EP,P,P,p,P,P,Y) ELSE IF (:CTURB: .EQ. TURB) THEN + STORE(LEN1,VIST) ENDIF GROUP 8. Terms (in differential equations) & devices *** Get rid of potentially confusing heat sources that are *** irrelevant to the tests being performed TERMS(TEM1,N,P,P,P,P,P) *** Use the upwind scheme to get more easily understandable *** behaviour DIFCUT=0.0 GROUP 9. Properties of the medium (or media) *** Set up the properties of the materials to be used *** NB. the viscosity is boosted to make the viscous *** effects more noticeable: note that the conductivity *** of the gas must also be increased or else the wall- *** function code will fail because of the unreasonably *** small laminar Prandtl No. CSG10='q1' MATFLG=T;NMAT=2 56 1.0 1.e-5 1000. 0.025 4.0e-3 155 8000. 1.0 500. 50. 0.0 *** Set up the turbulence models for the turbulent cases *** NB. at present limited to an artificially simple mixing *** length model IF (:CTURB: .EQ. TURB .OR. :CTURB: .EQ. KLMODL) THEN + IF (NX.GT.1) THEN + EL1=LINEARX + ELSE + EL1=LINEARY + ENDIF + EL1A=(XULAST**2+YVLAST**2+ZWLAST**2)**0.5;EL1B=0.0 ENDIF IF (:CTURB: .EQ. TURB) THEN + ENUT=PROPLEN;ENUTA=0.0;ENUTB=0.1 ENDIF GROUP 11. Initialization of variable or porosity fields FIINIT(U1)=0.0;FIINIT(V1)=VIN;FIINIT(W1)=0.0 FIINIT(PRPS)=56;PRNDTL(TEM1)=CONDFILE;ENUL=FILE; FIINIT(TEM1)=10.;FIINIT(C1)=0.0 IF (.NOT. EGWF) THEN *** An unsuccesful attempt to bamboozle CONPOR + IF (NX .GT. 1) THEN + CONPOR(BLCKACP,-1,CELL,1,-(NXF-1),1,NY,1,NZ) + CONPOR(BLCKBCP,-1,CELL,-(NXL+1),NX,1,NY,1,NZ) + ENDIF + IF (NZ .GT. 1) THEN + CONPOR(BLCKCCP,-1,CELL,1,NX,1,NY,1,-(NZF-1)) + CONPOR(BLCKDCP,-1,CELL,1,NX,1,NY,-(NZL+1),NZ) + ENDIF ENDIF CASE :CTURB: OF WHEN KEMODL,6 + FIINIT(KE) = 0.01*VIN*VIN + FIINIT(EP) = FIINIT(KE)**1.5 WHEN KLMODL,6 + FIINIT(KE) = 0.01*VIN*VIN ENDCASE INIADD=F IF (NX .GT. 1) THEN + PATCH(BLCKA ,CELL,1,NXF-1,1,NY,1,NZ,1,1) + PATCH(BLCKAI,INIVAL,1,NXF-1,1,NY,1,NZ,1,1) + PATCH(BLCKB ,CELL,NXL+1,NX,1,NY,1,NZ,1,1) + PATCH(BLCKBI,INIVAL,NXL+1,NX,1,NY,1,NZ,1,1) + INIT (BLCKAI,PRPS,0.0,155.);INIT (BLCKBI,PRPS,0.0,155.) ENDIF IF (NZ .GT. 1) THEN + PATCH(BLCKC ,CELL,1,NX,1,NY,1,NZF-1,1,1) + PATCH(BLCKCI,INIVAL,1,NX,1,NY,1,NZF-1,1,1) + PATCH(BLCKD ,CELL,1,NX,1,NY,NZL+1,NZ,1,1) + PATCH(BLCKDI,INIVAL,1,NX,1,NY,NZL+1,NZ,1,1) + INIT (BLCKCI,PRPS,0.0,155.);INIT (BLCKDI,PRPS,0.0,155.) ENDIF IF (ISOLBK .NE. 0) THEN + INTEGER(IBL);IBL=1 + FIINIT(BLOK)=IBL + IF (NX .GT. 1) THEN + IBL=IBL+1;INIT(BLCKAI,BLOK,0.0,IBL) + IBL=IBL+1;INIT(BLCKBI,BLOK,0.0,IBL) + ENDIF + IF (NZ .GT. 1) THEN + IBL=IBL+1;INIT(BLCKCI,BLOK,0.0,IBL) + IBL=IBL+1;INIT(BLCKDI,BLOK,0.0,IBL) + ENDIF ENDIF GROUP 12. Patchwise adjustment of terms GROUP 12. Convection and diffusion adjustments GROUP 13. Boundary conditions and special sources *** Set the C1 and TEM1 sbources in the blocks: fixed *** flux for heat(??), and fixed values of C1 (heated *** rock-salt dissolving into water??) IF (NX .GT. 1) THEN + COVAL(BLCKA,C1,FIXVAL,-CIN) + COVAL(BLCKB,C1,FIXVAL,-CIN) + COVAL(BLCKA,TEM1,FIXFLU,20.) + COVAL(BLCKB,TEM1,FIXFLU,20.) ENDIF IF (NZ .GT. 1) THEN + COVAL(BLCKC,C1,FIXVAL,-CIN) + COVAL(BLCKD,C1,FIXVAL,-CIN) + COVAL(BLCKC,TEM1,FIXFLU,20.) + COVAL(BLCKD,TEM1,FIXFLU,20.) ENDIF *** Setup inlet conditions IF (NY .GT. 1) THEN + INLET(YINL,SOUTH,NXF,NXL,1,1,nzf,nzl,1,LSTEP) + VALUE(YINL,P1,RHO1*VIN) + VALUE(YINL,V1,VIN) + VALUE(YINL,TEM1,TIN);VALUE(YINL,C1,CIN) + IF (FIINIT(KE) .GT. 1.E-10) THEN + VALUE(YINL,KE,FIINIT(KE)) + ENDIF + IF (FIINIT(EP) .GT. 1.E-10) THEN + VALUE(YINL,EP,FIINIT(EP)) + ENDIF + OUTLET(YOUTL,NORTH,NXF,NXL,NY,NY,nzf,nzl,1,LSTEP) + VALUE(YOUTL,P1,0.0) ENDIF *** For the SOWF method set up the WALL PATCHes IF (.NOT. EGWF) THEN + IF (NX .GT. 1) THEN WALL(WALLA,WEST,NXF,NXF,1,NY,1,NZ,1,LSTEP) COVAL(WALLA,V1,WALLCO,OPPVAL) COVAL(WALLA,TEM1,WALLCO,OPPVAL);COVAL(WALLA,C1,WALLCO,OPPVAL) WALL(WALLB,EAST,NXL,NXL,1,NY,1,NZ,1,LSTEP) COVAL(WALLB,V1,WALLCO,OPPVAL) COVAL(WALLB,TEM1,WALLCO,OPPVAL);COVAL(WALLB,C1,WALLCO,OPPVAL) INIT(BLCKA,EPOR,0.0,1.0);INIT(BLCKB,EPOR,0.0,1.0) + ENDIF + IF (NZ .GT. 1) THEN WALL(WALLC,LOW ,1,NX,1,NY,NZF,NZF,1,LSTEP) COVAL(WALLC,V1,WALLCO,OPPVAL) COVAL(WALLC,TEM1,WALLCO,OPPVAL);COVAL(WALLC,C1,WALLCO,OPPVAL) WALL(WALLD,HIGH,1,NX,1,NY,NZL,NZL,1,LSTEP) COVAL(WALLD,V1,WALLCO,OPPVAL) COVAL(WALLD,TEM1,WALLCO,OPPVAL);COVAL(WALLD,C1,WALLCO,OPPVAL) INIT(BLCKC,EPOR,0.0,1.0);INIT(BLCKD,EPOR,0.0,1.0) + ENDIF ENDIF GROUP 15. Termination of sweeps LSWEEP=50 GROUP 16. Termination of iterations LITER(TEM1)=20;LITER(C1)=20;LITER(U1)=10;LITER(V1)=10 GROUP 17. Under-relaxation devices RELAX(U1,FALSDT,5.*XULAST/(NX*VIN)) RELAX(V1,FALSDT,5.*YVLAST/(NY*VIN)) RELAX(W1,FALSDT,5.*ZWLAST/(NZ*VIN)) RELAX(KE,FALSDT,5.*YVLAST/(NY*VIN)) RELAX(EP,FALSDT,5.*YVLAST/(NY*VIN)) RELAX(TEM1,FALSDT,1000.) RELAX(C1, FALSDT,1000.) GROUP 21. Print-out of variables OUTPUT(P1,Y,P,P,Y,Y,Y);OUTPUT(U1,Y,P,P,Y,Y,Y) OUTPUT(V1,Y,P,P,Y,Y,Y);OUTPUT(W1,Y,P,P,Y,Y,Y) OUTPUT(TEM1,Y,P,P,Y,Y,Y);OUTPUT(C1,Y,P,P,Y,Y,Y) IF (:CTURB: .EQ. KEMODL .OR. :CTURB: .EQ. KLMODL) THEN + OUTPUT(KE,Y,P,P,Y,Y,Y) ENDIF GROUP 22. Spot-value print-out IXMON=NXF;IYMON=NY/2+1;IZMON=NZ/2+1 TSTSWP=-1 GROUP 23. Field print-out and plot control NXPRIN=1;NYPRIN=1;NZPRIN=1;IYPRF=NY-2;IYPRL=NY IF (NX .GT. 1 .AND. NY .GT. 1) THEN + PATCH(CONZ,CONTUR,1,NX,1,NY,NZ/2+1,NZ/2+1,1,1) + PLOT (CONZ,P1,0.0,20.) + PLOT (CONZ,TEM1,0.0,20.) + PLOT (CONZ,C1,0.0,20.) ENDIF IF (NY .GT. 1 .AND. NZ .GT. 1) THEN + PATCH(CONX,CONTUR,NX/2+1,NX/2+1,1,NY,1,NZ,1,1) + PLOT (CONX,P1,0.0,20.) + PLOT (CONX,TEM1,0.0,20.) + PLOT (CONX,C1,0.0,20.) ENDIF IF (NZ .GT. 1 .AND. NX.GT. 1) THEN + PATCH(CONY,CONTUR,1,NX,NY/2+1,NY/2+1,1,NZ,1,1) + PLOT (CONY,P1,0.0,20.) + PLOT (CONY,TEM1,0.0,20.) + PLOT (CONY,C1,0.0,20.) ENDIF IF (NX .GT. 1) THEN + PATCH(PRFX,PROFIL,1,NX,NY/2+1,NY/2+1,NZ/2+1,NZ/2+1,1,1) + PLOT (PRFX,TEM1,100.0,1000.) + PLOT (PRFX,C1,-1.,1.) ENDIF IF (NY .GT. 1) THEN + PATCH(PRFY,PROFIL,NX/2+1,NX/2+1,1,NY,NZ/2+1,NZ/2+1,1,1) + PLOT (PRFY,TEM1,100.0,1000.) + PLOT (PRFY,C1,-1.,1.) ENDIF IF (NZ .GT. 1) THEN + PATCH(PRFZ,PROFIL,NX/2+1,NX/2+1,NY/2+1,NY/2+1,1,NZ,1,1) + PLOT (PRFZ,TEM1,100.0,1000.) + PLOT (PRFZ,C1,-1.,1.) ENDIF GROUP 24. Dumps for restarts