file ylink.q1 .... dbs .... november 1992 PHOTON USE p;;; vec z 1 msg - msg Pressto continue pause msg pressure contours con p1 fi z 1 1 10 1 5 0.01 msg - msg Press to continue pause con p1 fi z 1 1 10 8 12 0.01 msg - msg Press to continue pause con off red msg temperature contours con temp fi z 1 1 10 1 5 0.01 msg - msg Press to continue pause con temp fi z 1 1 10 8 12 0.01 vec z 1 msg - msg Press e to END enduse GROUP 1. Run title and other preliminaries TEXT(Y-Direction Link In XY Plane TITLE mesg(PC486/50 time last reported as 1.5 min DISPLAY y-direction link in xy plane 2-dimensional (x-y), Cartesian, steady, elliptic simulation This file, ylink.q1, tests the "link" built in to PHOENICS It concerns steady flow through a duct of shape formed by the lateral displacement of one half of a rectangular duct relative to the other. enddis REAL(XLENGTH,YLENGTH,ZLENGTH,FLORES,P1CO) INTEGER(NYNOM,IXSHFT,ISHIFT,IXPLUS,IXLFT) BOOLEAN(B2T) B2T=F FLORES=0.0 P1CO=1.0 XLENGTH=1.0;YLENGTH=1.0;ZLENGTH=1.0 NYNOM=10;NX=10;NZ=1 NPHI=20 NY=NYNOM+2 DO II=1,5 + MESG( ENDDO MESG(NX=:NX:; NYNOM=:NYNOM:; NZ=:NZ: MESGB(XLENGTH=:XLENGTH:; YLENGTH=:YLENGTH:; ZLENGTH=:ZLENGTH: DELAY(200) IF(NX.GT.1) THEN MESGm(What x-direction shift would you like ? (from 1-nx to nx-1 ) READVDU(IXSHFT,INT,0) ENDIF ISHIFT=2-IXSHFT*NY mesgm(ishift is the offset between linked points, =2 - ixshft*ny ISHIFT mesgm(Flow-resistance factor = :flores: If not OK, insert your value READVDU(FLORES,REAL,:FLORES:) FLORES MESGM(FLOW IS FROM BOTTOM TO TOP. OK? (Y/N) B2T=T ans=y READVDU(ANS,CHAR,Y) IF(:ANS:.EQ.N) THEN B2T=F B2T ENDIF IXLFT=NX/2+1 mesgm(inflow is from :ixlft: to nx. If not OK, insert lowest ix READVDU(IXLFT,INT,IXLFT) IXLFT GROUP 3. X-direction grid specification **Domain is XLENGTH m long in x-direction, with equal intervals GRDPWR(X,NX,XLENGTH,1.0) YVLAST=YLENGTH YFRAC(1)=-NYNOM/2;YFRAC(2)=1/NYNOM YFRAC(3)=2;YFRAC(4)=0.01*YFRAC(2) YFRAC(5)=-YFRAC(1);YFRAC(6)=YFRAC(2) GROUP 5. Z-direction grid specification **Domain is ZLENGTH m long in z-direction, with equal intervals GRDPWR(Z,NZ,ZLENGTH,1.0) GROUP 7. Variables stored, solved & named **Choose first-phase enthalpy (H1) as dependent variable and activate the whole-field elliptic solver SOLUTN(H1,Y,Y,Y,N,N,N);NAME(H1)=TEMP STORE(VPOR);SOLVE(P1,U1,V1);SOLUTN(P1,Y,Y,Y,N,N,N) GROUP 8. Terms (in differential equations) & devices **For pure conduction, cut out built-in source and convection terms TERMS(TEMP,N,N,Y,N,Y,Y) GROUP 9. Properties of the medium (or media) **Thermal conductivity will be ENUL*RHO1/PRNDTL(TEMP), so : ENUL=1.0E-3;RHO1=1.0;PRNDTL(TEMP)=1.0 GROUP 11. Initialization of variable or porosity fields INIADD=F IF(B2T) THEN FIINIT(V1)=0.5 ELSE FIINIT(V1)=-0.5 ENDIF FIINIT(VPOR)=1.0 IF(IXSHFT.GT.0) THEN CONPOR(BAR1,0.0,CELL,1,IXSHFT,NYNOM/2+1,NYNOM/2+1,1,NZ) CONPOR(BAR2,0.0,CELL,NX-IXSHFT+1,NX,NYNOM/2+2,NYNOM/2+2,1,NZ) ENDIF IF(IXSHFT.LT.0) THEN CONPOR(BAR1,0.0,CELL,NX+1+IXSHFT,NX,NYNOM/2+1,NYNOM/2+1,1,NZ) CONPOR(BAR2,0.0,CELL,1,-IXSHFT,NYNOM/2+2,NYNOM/2+2,1,NZ) ENDIF FIINIT(TEMP)=0.0 PATCH(LOWER,INIVAL,1,NX,1,NYNOM/2,1,1,1,1) COVAL(LOWER,P1,0.0,1.0) COVAL(LOWER,TEMP,0.0,-0.9) GROUP 13. Boundary conditions and special sources ** inflow boundary IF(B2T) THEN PATCH(COLD,SOUTH,IXLFT,NX,1,1,1,1,1,1) COVAL(COLD,V1,ONLYMS,1.0) ELSE PATCH(COLD,SOUTH,IXLFT,NX,NY,NY,1,1,1,1) COVAL(COLD,V1,ONLYMS,-1.0) ENDIF COVAL(COLD,TEMP,1.E5,-0.9) COVAL(COLD,P1,FIXFLU,1.0) ** hot boundary IF(B2T) THEN PATCH(HOT,CELL,1,NX,NY,NY,1,1,1,1) ELSE PATCH(HOT,CELL,1,NX,1,1,1,1,1,1) ENDIF COVAL(HOT,TEMP,1.E5,0.9) COVAL(HOT,P1,1.E-2,0.0) ixplus is introduced in an attempt to counrteract the effect of th being no u cells at the west and east boundaries IXPLUS=0 IF(IXSHFT.GT.0) THEN IXPLUS=-1 ENDIF IF(IXSHFT.LT.0) THEN IXPLUS=1 ENDIF IXPLUS=-IXPLUS IXPLUS=0 MESGM(IXPLUS = :IXPLUS: CORRECT? IF(IXSHFT.GE.0) THEN PATCH(+1,SOUTH,1+IXSHFT,NX,NYNOM/2+1,NYNOM/2+1,1,1,1,1) ELSE PATCH(+1,SOUTH,1,NX+IXSHFT,NYNOM/2+1,NYNOM/2+1,1,1,1,1) ENDIF COVAL(+1,TEMP,FIXVAL,ISHIFT) COVAL(+1,P1,FIXVAL,ISHIFT) COVAL(+1,U1,FIXVAL,ISHIFT+IXPLUS) COVAL(+1,V1,FIXVAL,ISHIFT) INTEGER(III);III=NYNOM/2+2 IF(IXSHFT.GE.0) THEN PATCH(+2, SOUTH,1,NX-IXSHFT,III,III,1,1,1,1) ELSE PATCH(+2,SOUTH, 1-IXSHFT,NX,III,III,1,1,1,1) ENDIF COVAL(+2,TEMP,FIXVAL,-ISHIFT) COVAL(+2,P1,FIXVAL,-ISHIFT) COVAL(+2,V1,FIXVAL,-ISHIFT) COVAL(+2,U1,FIXVAL,-ISHIFT-IXPLUS) PATCH(FLOWRES,VOLUME,1,NX,1,NY,1,NZ,1,1) COVAL(FLOWRES,U1,FLORES,0.0) COVAL(FLOWRES,V1,FLORES,0.0) PATCH(GP12DFN1,NORTH,1,NX,NYNOM/2,NYNOM/2,1,1,1,1) COVAL(GP12DFN1,TEMP,0.5,0.0) PATCH(GP12DFN2,NORTH,1,NX,NYNOM/2+2,NYNOM/2+2,1,1,1,1) COVAL(GP12DFN2,TEMP,0.5,0.0) GROUP 15. Termination of sweeps LSWEEP=20;RESREF(P1)=1.E-10 GROUP 16. Termination of iterations ** Set the frequencies of application of the one-dimensional correction features in the linear-equation solver to once per iteration for each direction. ISOLX=1;ISOLY=1;ISOLZ=1 LITER(TEMP)=10;LITER(P1)=20 GROUP 17. Under-relaxation devices RELAX(P1,LINRLX,0.5);RELAX(V1,FALSDT,0.05);RELAX(U1,FALSDT,0.05) GROUP 21. Print-out of variables **Print fields of temperature OUTPUT(TEMP,Y,Y,Y,Y,Y,Y) GROUP 22. Spot-value print-out IYMON=NY/2+1;IZMON=NZ/2+1;IXMON=NX-1;ITABL=1 GROUP 23. Field print-out and plot control IXPRF=NYNOM/2-1;IXPRL=NYNOM/2+3 NXPRIN=1;NYPRIN=1 **Plot a profile along the line Iy=ny/2 PATCH(YLINE,PROFIL,1,NX,NY/2,NY/2,1,1,1,1) PLOT(YLINE,TEMP,0.0,0.0);PLOT(YLINE,P1,0.0,0.0) PLOT(YLINE,U1,0.0,0.0) ** Plot contour diagrams for the plane PATCH(FIRST,CONTUR,1,NYNOM/2,1,NY,1,NZ,1,1) PLOT(FIRST,TEMP,0.0,20.0);PLOT(FIRST,P1,0.0,20.0) PATCH(SECOND,CONTUR,NYNOM/2+3,NX,1,NY,1,NZ,1,1) PLOT(SECOND,TEMP,0.0,20.0);PLOT(SECOND,P1,0.0,20.0) GROUP 24. Dumps for restarts UWATCH=T;NPLT=1;TSTSWP=-1 LSWEEP=200 SELREF=T;RESFAC=1.E-2 IXPRF=IXPRF+1;IXPRL=IXPRL+1 LSWEEP=5;NPRINT=1;TSTSWP=-1;IXPRF=1;IXPRL=NX IF(NY.EQ.1) THEN SOLUTN(V1,N,N,N,N,N,N);OUTPUT(V1,N,N,N,N,N,N) ENDIF NPRINT=10000 LSWEEP=200;TSTSWP=-1 NYPRIN=1;IYPRF=5;IYPRL=8;NXPRIN=1 IXPRF=4;IXPRL=8;IXPRF=1;IXPRL=NX INIFLD=T FIINIT(TEMP)=0.9 CARTES=F RINNER=YVLAST/10;NPRINT=LSWEEP lsg57=t