diff options
Diffstat (limited to 'Dragon/src/DOORS_MOD.f90')
| -rw-r--r-- | Dragon/src/DOORS_MOD.f90 | 162 |
1 files changed, 162 insertions, 0 deletions
diff --git a/Dragon/src/DOORS_MOD.f90 b/Dragon/src/DOORS_MOD.f90 new file mode 100644 index 0000000..9538892 --- /dev/null +++ b/Dragon/src/DOORS_MOD.f90 @@ -0,0 +1,162 @@ +MODULE DOORS_MOD + USE GANLIB +CONTAINS + SUBROUTINE DOORS(CDOOR,IPTRK,NMAT,NANIS,NUN,SIGG,SUNKNO,FUNKNO) + ! + !--------------------------------------------------------------------- + ! + !Purpose: + ! compute the product of a cross section times a flux unknow vector. + ! + !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 + ! CDOOR name of the geometry/solution operator. + ! IPTRK pointer to the tracking (L_TRACK signature). + ! NMAT number of mixtures in the macrolib. + ! NANIS number of Legendre components in the macrolib (=0: isotropic). + ! NUN total number of unknowns in vectors SUNKNO and FUNKNO. + ! SIGG cross section. + ! FUNKNO optional unknown vector. If not present, a flat flux + ! approximation is assumed. + ! + !Parameters: output + ! SUNKNO source vector. Volumes are included with BIVAC and TRIVAC + ! trackings. + ! + !--------------------------------------------------------------------- + ! + !---- + ! SUBROUTINE ARGUMENTS + !---- + CHARACTER(LEN=12), INTENT(IN) :: CDOOR + TYPE(C_PTR), INTENT(IN) :: IPTRK + INTEGER, INTENT(IN) :: NMAT,NANIS,NUN + REAL, DIMENSION(0:NMAT,NANIS+1), INTENT(IN) :: SIGG + REAL, DIMENSION(NUN), INTENT(IN), OPTIONAL :: FUNKNO + REAL, DIMENSION(NUN), INTENT(INOUT) :: SUNKNO + !---- + ! LOCAL VARIABLES + !---- + INTEGER, PARAMETER :: NSTATE=40 + INTEGER, DIMENSION(NSTATE) :: ISTATE + !---- + ! ALLOCATABLE ARRAYS + !---- + INTEGER, ALLOCATABLE, DIMENSION(:) :: MATCOD + INTEGER, ALLOCATABLE, DIMENSION(:,:,:) :: KEYFLX + REAL, ALLOCATABLE, DIMENSION(:) :: VOL + !---- + ! RECOVER TRACKING PARAMETERS + ! NFUNL: number of spherical harmonics components used to expand the + ! flux and the sources. + ! NANIS_TRK: number of components in the angular expansion of the flux + !---- + CALL LCMGET(IPTRK,'STATE-VECTOR',ISTATE) + NREG=ISTATE(1) + IF(ISTATE(2).GT.NUN) CALL XABORT('DOORS: WRONG NUN.') + IF(ISTATE(4).GT.NMAT) CALL XABORT('DOORS: WRONG NMAT.') + NDIM=0 + NLIN=1 + NFUNL=1 + NANIS_TRK=1 + IF(CDOOR.EQ.'MCCG') THEN + NANIS_TRK=ISTATE(6) + NDIM=ISTATE(16) + CALL LCMGET(IPTRK,'MCCG-STATE',ISTATE) + NFUNL=ISTATE(19) + NLIN=ISTATE(20) + ELSE IF(CDOOR.EQ.'SN') THEN + NFUNL=ISTATE(7) + NLIN=ISTATE(8) + NDIM=ISTATE(9) + NLIN=NLIN**NDIM + NLIN=NLIN*ISTATE(35) + NANIS_TRK=ISTATE(16) + ELSE IF(CDOOR.EQ.'BIVAC') THEN + NLIN=ABS(ISTATE(8)) ! order of finite elements + NFUNL=MAX(1,ISTATE(14)) + NANIS_TRK=ABS(ISTATE(16)) + ELSE IF(CDOOR.EQ.'TRIVAC') THEN + NLIN=ABS(ISTATE(9)) ! order of finite elements + NLFUNL=MAX(1,ISTATE(30)) + NANIS_TRK=ABS(ISTATE(32)) + ENDIF + ALLOCATE(MATCOD(NREG),VOL(NREG),KEYFLX(NREG,NLIN,NFUNL)) + KEYFLX(:NREG,:NLIN,:NFUNL)=0 + CALL LCMLEN(IPTRK,'MATCOD',ILNLCM,ITYLCM) + IF(ILNLCM.NE.NREG) CALL XABORT('DOORS: INCOMPATIBLE NUMBER OF REGIONS.') + CALL LCMGET(IPTRK,'MATCOD',MATCOD) + CALL LCMGET(IPTRK,'VOLUME',VOL) + IF((CDOOR.EQ.'MCCG').OR.(CDOOR.EQ.'SN')) THEN + CALL LCMGET(IPTRK,'KEYFLX$ANIS',KEYFLX) + ELSE + CALL LCMGET(IPTRK,'KEYFLX',KEYFLX) + ENDIF + !---- + ! PERFORM SIGG*FUNKNO MULTIPLICATION + !---- + IF(CDOOR.EQ.'SN') THEN + CALL DOORS_SN(IPTRK,NANIS,NREG,NMAT,NUN,MATCOD,SIGG,SUNKNO,FUNKNO) + ELSE IF(CDOOR.EQ.'BIVAC') THEN + CALL DOORS_BIV(IPTRK,NANIS,NREG,NMAT,NUN,MATCOD,VOL,SIGG,SUNKNO,FUNKNO) + ELSE IF(CDOOR.EQ.'TRIVAC') THEN + CALL DOORS_TRI(IPTRK,NANIS,NREG,NMAT,NUN,MATCOD,VOL,SIGG,SUNKNO,FUNKNO) + ELSE IF(PRESENT(FUNKNO)) THEN + ! general case + IF((NANIS.EQ.0).OR.(NFUNL.EQ.1).OR.(NANIS_TRK.EQ.1)) THEN + ! LAB isotropy or transport correction + DO IR=1,NREG + IBM=MATCOD(IR) + IF(IBM.LE.0) CYCLE + DO IE=1,NLIN + IND=KEYFLX(IR,IE,1) + SUNKNO(IND)=SUNKNO(IND)+SIGG(IBM,1)*FUNKNO(IND) + ENDDO + ENDDO ! IR + ELSE + ! spherical harmonics expansion of the flux and source + DO IR=1,NREG + IBM=MATCOD(IR) + IF(IBM.LE.0) CYCLE + DO IAL=0,MIN(NFUNL-1,NANIS,NANIS_TRK-1) + FACT=REAL(2*IAL+1) + DO IAM=0,MIN(NFUNL-1,NANIS,NANIS_TRK-1) + DO IE=1,NLIN + IND=0 + IF(NDIM.EQ.3) THEN + IND=KEYFLX(IR,IE,1+IAL*NANIS_TRK+IAM) + ELSE IF((NDIM.EQ.2).AND.(IAM.LE.IAL)) THEN + IND=KEYFLX(IR,IE,1+IAL*(IAL+1)/2+IAM) + ELSE IF(IAM.EQ.IAL) THEN + IND=KEYFLX(IR,IE,1+IAL) + ENDIF + IF(IND.EQ.0) THEN + CYCLE + ELSE IF(IND.GT.NUN) THEN + CALL XABORT('DOORS: NUN OVERFLOW.') + ENDIF + SUNKNO(IND)=SUNKNO(IND)+FACT*SIGG(IBM,IAL+1)*FUNKNO(IND) + ENDDO ! IE + ENDDO ! IAM + ENDDO ! IAL + ENDDO ! IR + ENDIF + ELSE + ! general case (flat flux) + DO IR=1,NREG + IND=KEYFLX(IR,1,1) + SUNKNO(IND)=SUNKNO(IND)+SIGG(MATCOD(IR),1) + ENDDO + ENDIF + DEALLOCATE(KEYFLX,VOL,MATCOD) + END SUBROUTINE DOORS +END MODULE DOORS_MOD + |