summaryrefslogtreecommitdiff
path: root/Dragon/src/DOORS_TRI.f90
diff options
context:
space:
mode:
Diffstat (limited to 'Dragon/src/DOORS_TRI.f90')
-rw-r--r--Dragon/src/DOORS_TRI.f90252
1 files changed, 252 insertions, 0 deletions
diff --git a/Dragon/src/DOORS_TRI.f90 b/Dragon/src/DOORS_TRI.f90
new file mode 100644
index 0000000..4c17ebd
--- /dev/null
+++ b/Dragon/src/DOORS_TRI.f90
@@ -0,0 +1,252 @@
+SUBROUTINE DOORS_TRI(IPTRK,NANIS,NREG,NMAT,NUN,MATCOD,VOL,SIGG,SUNKNO,FUNKNO)
+ !
+ !-----------------------------------------------------------------------
+ !
+ !Purpose:
+ ! Compute the source for the solution of diffusion or PN equations.
+ ! Use a TRIVAC tracking.
+ !
+ !Copyright:
+ ! Copyright (C) 2025 Ecole Polytechnique de Montreal
+ ! This library is free software; you can redistribute it and/or
+ ! modify it under the terms of the GNU Lesser General Public
+ ! License as published by the Free Software Foundation; either
+ ! version 2.1 of the License, or (at your option) any later version
+ !
+ !Author(s): A. Hebert
+ !
+ !Parameters: input
+ ! IPTRK pointer to the tracking LCM object.
+ ! NANIS maximum cross section Legendre order (=0: isotropic).
+ ! NREG number of regions.
+ ! NMAT number of mixtures.
+ ! NUN number of unknowns per energy group including net current.
+ ! MATCOD mixture indices.
+ ! VOL volumes. Volumes are included in SUNKNO.
+ ! SIGG cross section.
+ ! FUNKNO optional unknown vector. If not present, a flat flux
+ ! approximation is assumed.
+ !
+ !Parameters: input/output
+ ! SUNKNO integrated sources.
+ !
+ !-----------------------------------------------------------------------
+ !
+ USE GANLIB
+ !----
+ ! SUBROUTINE ARGUMENTS
+ !----
+ TYPE(C_PTR) IPTRK
+ INTEGER NANIS,NREG,NMAT,NUN,MATCOD(NREG)
+ REAL VOL(NREG),SIGG(0:NMAT,NANIS+1),SUNKNO(NUN)
+ REAL, OPTIONAL :: FUNKNO(NUN)
+ !----
+ ! LOCAL VARIABLES
+ !----
+ PARAMETER(NSTATE=40)
+ INTEGER JPAR(NSTATE)
+ !----
+ ! RECOVER BIVAC SPECIFIC PARAMETERS.
+ !----
+ CALL LCMGET(IPTRK,'STATE-VECTOR',JPAR)
+ IF(JPAR(1).NE.NREG) CALL XABORT('DOORS_TRI: INCONSISTENT NREG.')
+ IF(JPAR(2).NE.NUN) CALL XABORT('DOORS_TRI: INCONSISTENT NUN.')
+ IF(NANIS.NE.0) CALL XABORT('DOORS_TRI: SPN NOT IMPLEMENTED.')
+ ITYPE=JPAR(6)
+ IF((ITYPE.EQ.2).OR.(ITYPE.EQ.5).OR.(ITYPE.EQ.7)) THEN
+ ! Cartesian 1D, 2D or 3D geometry
+ CALL DOORS_TRIGSO(IPTRK,NREG,NMAT,NUN,MATCOD,VOL,SIGG,SUNKNO,FUNKNO)
+ ELSE IF((ITYPE.EQ.8).OR.(ITYPE.EQ.9)) THEN
+ ! Hexagonal 2D or 3D geometry
+ CALL DOORS_TRIGSR(IPTRK,NREG,NMAT,NUN,MATCOD,SIGG,SUNKNO,FUNKNO)
+ ELSE
+ CALL XABORT('DOORS_TRI: GEOMETRY TYPE NOT AVAILABLE.')
+ ENDIF
+ RETURN
+CONTAINS
+ SUBROUTINE DOORS_TRIGSO(IPTRK,NREG,NMAT,NUN,MATCOD,VOL,SIGG,SUNKNO,FUNKNO)
+ !
+ !-----------------------------------------------------------------------
+ !
+ !Purpose:
+ ! Source term calculation for a mixed-dual formulation of the finite
+ ! element technique in a 3-D Cartesian geometry.
+ !
+ !-----------------------------------------------------------------------
+ !
+ USE GANLIB
+ !----
+ ! SUBROUTINE ARGUMENTS
+ !----
+ TYPE(C_PTR) IPTRK
+ INTEGER NREG,NMAT,NUN,MATCOD(NREG)
+ REAL VOL(NREG),SIGG(0:NMAT),SUNKNO(NUN)
+ REAL, OPTIONAL :: FUNKNO(NUN)
+ !----
+ ! LOCAL VARIABLES
+ !----
+ PARAMETER(NSTATE=40)
+ INTEGER IPAR(NSTATE)
+ CHARACTER HSMG*131
+ !----
+ ! ALLOCATABLE ARRAYS
+ !----
+ INTEGER, ALLOCATABLE, DIMENSION(:) :: KN
+ !----
+ ! RECOVER TRIVAC SPECIFIC PARAMETERS.
+ !----
+ CALL LCMGET(IPTRK,'STATE-VECTOR',IPAR)
+ ITYPE=IPAR(6)
+ IELEM=IPAR(9)
+ ICOL=IPAR(10)
+ LX=IPAR(14)
+ LY=IPAR(15)
+ LZ=IPAR(16)
+ ISCAT=IPAR(32)
+ IF((ITYPE.NE.2).AND.(ITYPE.NE.5).AND.(ITYPE.NE.7)) THEN
+ CALL XABORT('DOORS_TRIGSO: INVALID CARTESIAN GEOMETRY.')
+ ELSE IF((IELEM.LT.0).OR.(ICOL.GT.3)) THEN
+ CALL XABORT('DOORS_TRIGSO: RAVIART-THOMAS METHOD EXPECTED.')
+ ELSE IF(ISCAT.GT.1) THEN
+ WRITE(HSMG,'(56HDOORS_TRIGSO: MACRO-CALCULATION WITH ANISOTROPIC SCATTER, &
+ & 66HING CURRENTLY NOT IMPLEMENTED; USE SCAT 1 KEYWORD IN TRIVAT: DATA.)')
+ CALL XABORT(HSMG)
+ ELSE IF(LX*LY*LZ.NE.NREG) THEN
+ CALL XABORT('DOORS_TRIGSO: INVALID NREG.')
+ ENDIF
+ CALL LCMLEN(IPTRK,'KN',MAXKN,ITYLCM)
+ ALLOCATE(KN(MAXKN))
+ CALL LCMGET(IPTRK,'KN',KN)
+ IF(PRESENT(FUNKNO)) THEN
+ NUM1=0
+ DO K=1,NREG
+ L=MATCOD(K)
+ IF(L.LE.0) CYCLE
+ IF(VOL(K).EQ.0.0) GO TO 10
+ GARS=SIGG(L)
+ DO I0=1,IELEM**3
+ IND1=KN(NUM1+1)+I0-1
+ SUNKNO(IND1)=SUNKNO(IND1)+FUNKNO(IND1)*VOL(K)*GARS
+ ENDDO ! I0
+ 10 NUM1=NUM1+1+6*IELEM**2
+ ENDDO ! K
+ ELSE
+ ! a flat flux is assumed
+ NUM1=0
+ DO K=1,NREG
+ L=MATCOD(K)
+ IF(L.LE.0) CYCLE
+ IF(VOL(K).EQ.0.0) GO TO 20
+ IND1=KN(NUM1+1)
+ SUNKNO(IND1)=SUNKNO(IND1)+VOL(K)*SIGG(L)
+ 20 NUM1=NUM1+1+6*IELEM**2
+ ENDDO ! K
+ ENDIF
+ DEALLOCATE(KN)
+ RETURN
+ END SUBROUTINE DOORS_TRIGSO
+ !
+ SUBROUTINE DOORS_TRIGSR(IPTRK,NREG,NMAT,NUN,MATCOD,SIGG,SUNKNO,FUNKNO)
+ !
+ !-----------------------------------------------------------------------
+ !
+ !Purpose:
+ ! Source term calculation for a Thomas-Raviart-Schneider formulation of
+ ! the finite element technique in a 3-D hexagonal geometry.
+ !
+ !-----------------------------------------------------------------------
+ !
+ USE GANLIB
+ !----
+ ! SUBROUTINE ARGUMENTS
+ !----
+ TYPE(C_PTR) IPTRK
+ INTEGER NREG,NMAT,NUN,MATCOD(3,NREG/3)
+ REAL SIGG(0:NMAT),SUNKNO(NUN)
+ REAL, OPTIONAL :: FUNKNO(NUN)
+ !----
+ ! LOCAL VARIABLES
+ !----
+ PARAMETER(NSTATE=40)
+ INTEGER IPAR(NSTATE)
+ CHARACTER HSMG*131
+ !----
+ ! ALLOCATABLE ARRAYS
+ !----
+ INTEGER, ALLOCATABLE, DIMENSION(:) :: IPERT
+ INTEGER, ALLOCATABLE, DIMENSION(:,:) :: KN
+ REAL, ALLOCATABLE, DIMENSION(:) :: FRZ
+ REAL, ALLOCATABLE, DIMENSION(:,:) :: ZZ
+ !----
+ ! RECOVER TRIVAC SPECIFIC PARAMETERS.
+ !----
+ CALL LCMGET(IPTRK,'STATE-VECTOR',IPAR)
+ ITYPE=IPAR(6)
+ IELEM=IPAR(9)
+ ICOL=IPAR(10)
+ LX=IPAR(14)
+ LZ=IPAR(16)
+ ISCAT=IPAR(32)
+ NBLOS=LX*LZ/3
+ IF((ITYPE.NE.8).AND.(ITYPE.NE.9)) THEN
+ CALL XABORT('DOORS_TRIGSR: INVALID HEXAGONAL GEOMETRY.')
+ ELSE IF((IELEM.LT.0).OR.(ICOL.GT.3)) THEN
+ CALL XABORT('DOORS_TRIGSR: RAVIART-THOMAS METHOD EXPECTED.')
+ ELSE IF(ISCAT.GT.1) THEN
+ WRITE(HSMG,'(56HDOORS_TRIGSR: MACRO-CALCULATION WITH ANISOTROPIC SCATTER, &
+ & 66HING CURRENTLY NOT IMPLEMENTED; USE SCAT 1 KEYWORD IN TRIVAT: DATA.)')
+ CALL XABORT(HSMG)
+ ELSE IF(3*NBLOS.NE.NREG) THEN
+ CALL XABORT('DOORS_TRIGSR: INVALID NREG.')
+ ENDIF
+ CALL LCMLEN(IPTRK,'KN',MAXKN,ITYLCM)
+ ALLOCATE(ZZ(3,NBLOS),KN(NBLOS,MAXKN/NBLOS),IPERT(NBLOS),FRZ(NBLOS))
+ CALL LCMGET(IPTRK,'SIDE',SIDE)
+ CALL LCMGET(IPTRK,'ZZ',ZZ)
+ CALL LCMGET(IPTRK,'KN',KN)
+ CALL LCMGET(IPTRK,'IPERT',IPERT)
+ CALL LCMGET(IPTRK,'FRZ',FRZ)
+ NELEM=IELEM*(IELEM+1)
+ TTTT=0.5*SQRT(3.0)*SIDE*SIDE
+ IF(PRESENT(FUNKNO)) THEN
+ NUM=0
+ DO KEL=1,NBLOS
+ IF(IPERT(KEL).EQ.0) CYCLE
+ NUM=NUM+1
+ L=MATCOD(1,IPERT(KEL))
+ IF(L.EQ.0) CYCLE
+ VOL0=TTTT*ZZ(1,IPERT(KEL))*FRZ(KEL)
+ GARS=SIGG(L)
+ DO K3=0,IELEM-1
+ DO K2=0,IELEM-1
+ DO K1=0,IELEM-1
+ JND1=(NUM-1)*IELEM**3+K3*IELEM**2+K2*IELEM+K1+1
+ JND2=(KN(NUM,1)-1)*IELEM**3+K3*IELEM**2+K2*IELEM+K1+1
+ JND3=(KN(NUM,2)-1)*IELEM**3+K3*IELEM**2+K2*IELEM+K1+1
+ SUNKNO(JND1)=SUNKNO(JND1)+FUNKNO(JND1)*VOL0*GARS
+ SUNKNO(JND2)=SUNKNO(JND2)+FUNKNO(JND2)*VOL0*GARS
+ SUNKNO(JND3)=SUNKNO(JND3)+FUNKNO(JND3)*VOL0*GARS
+ ENDDO ! K1
+ ENDDO ! K2
+ ENDDO ! K3
+ ENDDO ! KEL
+ ELSE
+ ! a flat flux is assumed
+ NUM=0
+ DO KEL=1,NBLOS
+ IF(IPERT(KEL).EQ.0) CYCLE
+ NUM=NUM+1
+ L=MATCOD(1,IPERT(KEL))
+ IF(L.EQ.0) CYCLE
+ VOL0=TTTT*ZZ(1,IPERT(KEL))*FRZ(KEL)
+ JND1=(NUM-1)*IELEM**3+1
+ JND2=(KN(NUM,1)-1)*IELEM**3+1
+ JND3=(KN(NUM,2)-1)*IELEM**3+1
+ SUNKNO(JND1)=SUNKNO(JND1)+VOL0*SIGG(L)
+ SUNKNO(JND2)=SUNKNO(JND2)+VOL0*SIGG(L)
+ SUNKNO(JND3)=SUNKNO(JND3)+VOL0*SIGG(L)
+ ENDDO ! KEL
+ ENDIF
+ DEALLOCATE(FRZ,IPERT,KN,ZZ)
+ END SUBROUTINE DOORS_TRIGSR
+END SUBROUTINE DOORS_TRI