c
c.... File Name ..... gxswfan.htm ... 070616
C ... Swirl functions for Fans
SUBROUTINE SWIRFAN
INCLUDE 'farray'
INCLUDE 'satear'
INCLUDE 'grdloc'
INCLUDE 'grdear'
PARAMETER (NPRP=50)
COMMON/GENI/NXNY,IGFIL1(8),NFM,IGF(21),IPRL,IBTAU,ILTLS,IGFIL(15)
1 ,ITEM1,ITEM2,ISPH1,ISPH2,ICON1,ICON2,IPRPS,IRADX,IRADY
1 ,IRADZ,IVFOL
COMMON/DRHODP/ITEMP,IDEN/DVMOD/IDVCGR
COMMON/NAMFN/NAMFUN,NAMSUB
COMMON/SWRFAN/RADU1,RADU2,GX0,GY0,GZ0,SWN,ICFAN,IFANN(NPRP)
1 ,INJDIR
COMMON/FANMS1/IFANT,VELARR(NPRP),NAMPCH(NPRP)
CHARACTER*6 NAMFUN,NAMSUB,NAMPCH*16
CHARACTER*14 COBJTP,OBNAME,CTEMP*8,INJDIR*1
LOGICAL MASKPT, L2DVRO
C
C***********************************************************************
NAMSUB='SWRFAN'
IXL=IABS(IXL)
C*****************************************************************
C--- GROUP 13. Boundary conditions and special sources
C Index for Coefficient - CO
C Index for Value - VAL
IF(IGR==13) THEN
IF(ISC==12.OR.ISC==17) THEN
C------------------- SECTION 12 ------------------- value = GRND
C
C USWIRL ----- swirl velocity & direction (+ve for clockwise
C into the solution domain)
C GXO,GY0,GZO ----- coordinates of burner axis
C INJDIR ----- direction of injection
C
ICFAN=0 ! initialise
IF(MASKPT(IPAT)) THEN
OBNAME=COBJTP(IPAT)
IF(OBNAME(1:3)=='FAN') THEN
C... Decomposition of Swirl Velocity Vector
CTEMP='^'//NPATCH
CALL GETSDI(CTEMP,'ICFAN',ICFAN)
CALL GETSDR(CTEMP,'XVEL1',XVEL1)
CALL GETSDR(CTEMP,'YVEL1',YVEL1)
CALL GETSDR(CTEMP,'ZVEL1',ZVEL1)
CALL GETSDR(CTEMP,'FVEL1',FVEL1)
CALL GETSDR(CTEMP,'DIAM1',DIAM1)
CALL GETSDR(CTEMP,'DIAM2',DIAM2)
CALL GETSDR(CTEMP,'GX0',GX0)
CALL GETSDR(CTEMP,'GY0',GY0)
CALL GETSDR(CTEMP,'GZ0',GZ0)
CALL GETSDR(CTEMP,'SWIRLN',SWN)
CALL GETSDC(CTEMP,'INJDIR',INJDIR)
RADU1=DIAM1/2.0
RADU2=DIAM2/2.0
USWIRL=ABS(FVEL1)*SWN
DO IC=1,IFANT
IF(NPATCH==NAMPCH(IC)) THEN
USWIRL=VELARR(IC)*SWN; FVEL1=VELARR(IC)
IFANN(IC)=ICFAN; EXIT
ENDIF
ENDDO
ENDIF
ENDIF
ENDIF
C ... If not swirl fan, return
IF(ICFAN/=2) RETURN
IF(ISC==12) THEN
C------------------- SECTION 12 ------------------- value = GRND
C
IF(MASKPT(IPAT)) THEN
OBNAME=COBJTP(IPAT)
IF(OBNAME(1:3)=='FAN') THEN
C
C... X-direction
C
IF(INJDIR=='X') THEN
C... Axial injection in x direction through y-z plane
IF(INDVAR==V1) THEN
L0VAL=L0F(VAL); L0ZGG=L0F(ZGNZ); L0YG=L0F(YG2D)
L0ZG=L0VAL
DO 13801 IX=IXF,IXL
IADD=NY*(IX-1)
DO 13801 IY=IYF,IYL
I=IY+IADD
F(L0ZG+I)=F(L0ZGG+IZ)
13801 CONTINUE
CALL GXSINA(L0VAL,L0YG,L0ZG,GY0,GZ0,NY)
IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; IPW=IPW+1
c ... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=F(L0VAL+I)*(-USWIRL)
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
ELSEIF(INDVAR==W1) THEN
L0VAL=L0F(VAL); L0ZGG=L0F(ZGNZ); L0YG=L0F(YG2D)
L0ZG=L0VAL
DO 13802 IX=IXF,IXL
IADD=NY*(IX-1)
DO 13802 IY=IYF,IYL
I=IY+IADD
F(L0ZG+I)=F(L0ZGG+IZ)
13802 CONTINUE
CALL GXCOSA(L0VAL,L0YG,L0ZG,GY0,GZ0,NY)
IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; IPW=IPW+1
c... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=F(L0VAL+I)*USWIRL
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
c... P1 fixflu
ELSEIF(INDVAR==P1) THEN
L0VAL=L0F(VAL); L0ZGG=L0F(ZGNZ); L0YG=L0F(YG2D)
L0DEN1=L0F(DEN1); IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; IPW=IPW+1
c... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=FVEL1*F(L0DEN1+I)
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
ELSEIF(INDVAR==U1) THEN
L0VAL=L0F(VAL); L0ZGG=L0F(ZGNZ); L0YG=L0F(YG2D); IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; ipw=ipw+1
c... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=XVEL1
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
ENDIF
C
C ... Y-direction
C
ELSEIF(INJDIR=='Y') THEN
C ... Axial injection in y direction through x-z plane
IF(INDVAR==U1) THEN
L0VAL=L0F(VAL); L0ZGG=L0F(ZGNZ); L0XG=L0F(XG2D)
L0ZG=L0VAL
DO 13803 IX=IXF,IXL
IADD=NY*(IX-1)
DO 13803 IY=IYF,IYL
I=IY+IADD
F(L0ZG+I)=F(L0ZGG+IZ)
13803 CONTINUE
CALL GXSINA(L0VAL,L0XG,L0ZG,GX0,GZ0,NY)
IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; IPW=IPW+1
c ... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=F(L0VAL+I)*USWIRL
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
ELSEIF(INDVAR==W1) THEN
L0VAL=L0F(VAL); L0ZGG=L0F(ZGNZ); L0XG=L0F(XG2D)
L0ZG=L0VAL
DO 13804 IX=IXF,IXL
IADD=NY*(IX-1)
DO 13804 IY=IYF,IYL
I=IY+IADD
F(L0ZG+I)=F(L0ZGG+IZ)
13804 CONTINUE
CALL GXCOSA(L0VAL,L0XG,L0ZG,GX0,GZ0,NY)
IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; IPW=IPW+1
c ... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=F(L0VAL+I)*(-USWIRL)
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
c ... P1 fixflu
ELSEIF(INDVAR==P1) THEN
L0VAL=L0F(VAL); L0ZGG=L0F(ZGNZ); L0XG=L0F(XG2D)
L0DEN1=L0F(DEN1); IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; IPW=IPW+1
c ... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=FVEL1*F(L0DEN1+I)
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
ELSEIF(INDVAR==V1) THEN
L0VAL=L0F(VAL); L0ZGG=L0F(ZGNZ); L0XG=L0F(XG2D); IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; IPW=IPW+1
c ... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=YVEL1
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
ENDIF
c
c... Z-direction
c
ELSEIF(INJDIR=='Z') THEN
C... Axial injection in z direction through y-x plane
IF(INDVAR==U1) THEN
L0VAL=L0F(VAL); L0YG=L0F(YG2D); L0XG=L0F(XG2D)
CALL GXSINA(L0VAL,L0XG,L0YG,GX0,GY0,NY)
IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; IPW=IPW+1
c ... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=F(L0VAL+I)*(-USWIRL)
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
ELSEIF(INDVAR==V1) THEN
L0VAL=L0F(VAL); L0YG=L0F(YG2D); L0XG=L0F(XG2D)
CALL GXCOSA(L0VAL,L0XG,L0YG,GX0,GY0,NY)
IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; IPW=IPW+1
c ... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=F(L0VAL+I)*USWIRL
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
c ... P1 fixflu
ELSEIF(INDVAR==P1) THEN
L0VAL=L0F(VAL); L0YG=L0F(YG2D); L0XG=L0F(XG2D)
L0DEN1=L0F(DEN1); IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; IPW=IPW+1
c ... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=FVEL1*F(L0DEN1+I)
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
ELSEIF(INDVAR==W1) THEN
L0VAL=L0F(VAL); L0YG=L0F(YG2D); L0XG=L0F(XG2D); IPW=0
DO IX=IXF,IXL
IADD=NY*(IX-1)
DO IY=IYF,IYL
I=IY+IADD; IPW=IPW+1
c ... Block centre area
IF(L2DVRO(IPW)) THEN
F(L0VAL+I)=ZVEL1
ELSE
F(L0VAL+I)=0.0
ENDIF
ENDDO
ENDDO
ENDIF
ENDIF
ENDIF
ENDIF
ENDIF
ENDIF
NAMSUB='swrfan'
END
c