LINK: The PHOENICS grid-restructuring feature
Fitting the structured grid of PHOENICS to practically-arising problems can be accomplished in a crude way by "blocking off", by way of zero porosity settings, spaces covered by the grid which are not in fact accessible to fluid or heat. This practice may be excessively wasteful of both storage and computer time.
Grid-restructuring (also called domain-decomposition) enables the the waste to be avoided by (in effect) removal of cells from regions where they are not needed to places where they are needed.
During the removal process, the cells can be reduced in size sufficiently to be placed inside other cells. This enables patches to be covered with small cells, without the necessity (usually entailed by a fully structured grid) of extending the refinement into places where it is not needed.
The following working Q1 file illustrates the use of the link scheme for the simulation of a flow in a duct having a shape of the following type:
------------------------
|
| ---->
-------------------|
---> ----------------------- | | -------------------| ^ | y | x -------------->
The entries in the Q1 are in upper-case type. Explanatory comments are in lower case. The case is one which features in the "active-demo" set, under the heading "grid-restructuring".
XLINK.Q1 ..... DBS ...... 01.11.92
TALK=F;RUN( 1, 1);VDU=VGACURSR
GROUP 1. RUN TITLE AND OTHER PRELIMINARIES
TEXT(X-DIRECTION LINK IN XY PLANE
MESGM(X-DIRECTION LINK IN XY PLANE
REAL(XLENGTH,YLENGTH,ZLENGTH,FLORES,P1CO)
INTEGER(NXNOM,IYSHFT,ISHIFT,IYPLUS,IYBOT)
CHAR(ANS)
BOOLEAN(L2R)
L2R=F
FLORES=0.0
P1CO=2.0E10
XLENGTH=1.0;YLENGTH=1.0;ZLENGTH=1.0
NXNOM=10;NY=10;NZ=1
NPHI=20
NX=NXNOM+2
MESGA(THIS FILE, XLINKXY.Q1, TESTS THE "LINK" BUILT INTO PHOENICS
MESG(IT CONCERNS STEADY FLOW THROUGH A DUCT OF SHAPE FORMED BY THE
MESG(LATERAL DISPLACEMENT OF ONE HALF OF A RECTANGULAR DUCT RELATIVE
MESG(TO THE OTHER.
MESG(NXNOM=:NXNOM:; NY=:NY:; NZ=:NZ:
MESGB(XLENGTH=:XLENGTH:; YLENGTH=:YLENGTH:; ZLENGTH=:ZLENGTH:
DELAY(500)
IF(NY.GT.1) THEN
MESGM(WHAT Y-DIRECTION SHIFT WOULD YOU LIKE ? (FROM 1-NY TO NY-1 )
READVDU(IYSHFT,INT,0)
ENDIF
ISHIFT=2*NY-IYSHFT
MESGM(ISHIFT IS THE OFFSET BETWEEN LINKED POINTS, =2*NY - IYSHFT
ISHIFT
MESGM(FLOW-RESISTANCE FACTOR = :FLORES: IF NOT OK, INSERT YOUR VALUE
READVDU(FLORES,REAL,:FLORES:)
FLORES
MESGM(PRESSURE COEFFICIENT = :P1CO: IF NOT OK, INSERT YOUR VALUE
READVDU(P1CO,REAL,:P1CO:)
P1CO
MESGM(FLOW IS FROM LEFT TO RIGHT. OK? (Y/N)
L2R=T
READVDU(ANS,CHAR,Y)
IF(:ANS:.EQ.N) THEN
L2R=F
L2R
ENDIF
IYBOT=NY/2+1
MESGM(INFLOW IS FROM :IYBOT: TO NY. IF NOT OK, INSERT BOTTOM IY
READVDU(IYBOT,INT,IYBOT)
IYBOT
GROUP 3. X-DIRECTION GRID SPECIFICATION
**DOMAIN IS XLENGTH M LONG IN X-DIRECTION, WITH EQUAL INTERVALS
GRDPWR(X,NX,XLENGTH,1.0)
XULAST=XLENGTH
XFRAC(1)=-NXNOM/2
XFRAC(2)=1/NXNOM
XFRAC(3)=2
XFRAC(4)=0.01*XFRAC(2)
XFRAC(5)=-XFRAC(1)
XFRAC(6)=XFRAC(2)
GROUP 4. Y-DIRECTION GRID SPECIFICATION
**DOMAIN IS YLENGTH M LONG IN Y-DIRECTION, WITH EQUAL INTERVALS
GRDPWR(Y,NY,YLENGTH,1.0)
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)
STORE(IMB1,PCOR)
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(L2R) THEN
FIINIT(U1)=0.5
ELSE
FIINIT(U1)=-0.5
ENDIF
FIINIT(VPOR)=1.0
IF(IYSHFT.GT.0) THEN
CONPOR(BAR1,0.0,CELL,NXNOM/2+1,NXNOM/2+1,1,IYSHFT,1,NZ)
CONPOR(BAR2,0.0,CELL,NXNOM/2+2,NXNOM/2+2, NY-IYSHFT+1,NY,1,NZ)
ENDIF
IF(IYSHFT.LT.0) THEN
CONPOR(BAR1,0.0,CELL,NXNOM/2+1,NXNOM/2+1,NY+1+IYSHFT,NY,1,NZ)
CONPOR(BAR2,0.0,CELL,NXNOM/2+2,NXNOM/2+2, 1,-IYSHFT,1,NZ)
ENDIF
GROUP 13. BOUNDARY CONDITIONS AND SPECIAL SOURCES
** COLD BOUNDARY ON THE LEFT
IF(L2R) THEN
PATCH(COLD,WEST,1,1,IYBOT,NY,1,1,1,1)
COVAL(COLD,U1,ONLYMS,1.0)
ELSE
PATCH(COLD,WEST,NX,NX,IYBOT,NY,1,1,1,1)
COVAL(COLD,U1,ONLYMS,-1.0)
ENDIF
COVAL(COLD,TEMP,1.E5,-0.9)
COVAL(COLD,P1,FIXFLU,1.0)
** HOT BOUNDARY ON THE RIGHT
IF(L2R) THEN
PATCH(HOT,CELL,NX,NX,1,NY,NZ,NZ,1,1)
ELSE
PATCH(HOT,CELL,1,1,IYBOT,NY,1,1,1,1)
ENDIF
COVAL(HOT,TEMP,1.E5,0.9)
COVAL(HOT,P1,1.E-2,0.0)
IYPLUS IS INTRODUCED IN AN ATTEMPT TO COUNTERACT THE EFFECT OF
THERE BEING NO V CELLS AT THE SOUTH AND NORTH BOUNDARIES
IYPLUS=0
IF(IYSHFT.GT.0) THEN
IYPLUS=-1
ENDIF
IF(IYSHFT.LT.0) THEN
IYPLUS=1
ENDIF
IYPLUS=-IYPLUS
IYPLUS=0
MESGM(IYPLUS = :IYPLUS: CORRECT?
IF(IYSHFT.GE.0) THEN
PATCH(+1,WEST,NXNOM/2+1,NXNOM/2+1,1+IYSHFT,NY,1,1,1,1)
PATCH(+1U,WEST,NXNOM/2,NXNOM/2,1+IYSHFT,NY,1,1,1,1)
ELSE
PATCH(+1,WEST,NXNOM/2+1,NXNOM/2+1,1,NY+IYSHFT,1,1,1,1)
PATCH(+1U,WEST,NXNOM/2,NXNOM/2,1,NY+IYSHFT,1,1,1,1)
ENDIF
COVAL(+1,TEMP,1.E5,ISHIFT)
COVAL(+1,P1,P1CO,ISHIFT)
COVAL(+1,U1,FIXVAL,ISHIFT)
COVAL(+1,V1,FIXVAL,ISHIFT+IYPLUS)
COVAL(+1U,U1,FIXVAL,ISHIFT)
INTEGER(III)
III=NXNOM/2+2
IF(IYSHFT.GE.0) THEN
PATCH(+2, WEST,III, III, 1,NY-IYSHFT,1,1,1,1)
PATCH(+2U,WEST,III-1,III-1,1,NY-IYSHFT,1,1,1,1)
ELSE
PATCH(+2,WEST,III,III, 1-IYSHFT,NY,1,1,1,1)
PATCH(+2U,WEST,III-1,III-1,1-IYSHFT,NY,1,1,1,1)
ENDIF
COVAL(+2,TEMP,1.E5,-ISHIFT)
COVAL(+2,P1,P1CO,-ISHIFT)
COVAL(+2,U1,FIXVAL,-ISHIFT)
COVAL(+2U,U1,FIXVAL,-ISHIFT)
COVAL(+2,V1,FIXVAL,-ISHIFT-IYPLUS)
PATCH(FLOWRES,VOLUME,1,NX,1,NY,1,NZ,1,1)
COVAL(FLOWRES,U1,FLORES,0.0)
COVAL(FLOWRES,V1,FLORES,0.0)
PATCH(GP12DFE1,EAST,NXNOM/2,NXNOM/2,1,NY,1,1,1,1)
COVAL(GP12DFE1,TEMP,0.5,0.0)
PATCH(GP12DFE2,EAST,NXNOM/2+2,NXNOM/2+2,1,NY,1,1,1,1)
COVAL(GP12DFE2,TEMP,0.5,0.0)
GROUP 15. TERMINATION OF SWEEPS
LSWEEP=20;RESREF(P1)=1.E-10
RESREF(U1)=1.E-30;RESREF(V1)=1.E-30
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(V1,FALSDT,0.05)
RELAX(U1,FALSDT,0.05)
GROUP 18. LIMITS ON VARIABLES OR INCREMENTS TO THEM
GROUP 19. DATA COMMUNICATED BY SATELLITE TO GROUND
GROUP 20. PRELIMINARY PRINT-OUT
ECHO=T
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=NXNOM/2-1;IXPRL=NXNOM/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,NXNOM/2,1,NY,1,NZ,1,1)
PLOT(FIRST,TEMP,0.0,20.0);PLOT(FIRST,P1,0.0,20.0)
PATCH(SECOND,CONTUR,NXNOM/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
SOLUTN(TEMP,N,N,N,N,N,N);OUTPUT(TEMP,N,N,N,N,N,N)
OUTPUT(VPOR,Y,Y,Y,Y,Y,Y)
NPRINT=10000;LSWEEP=200
STOP
See also the papaer by:..... in the 1992 PHOENICS User Conference.
wbs