From 7dfcc480ba1e19bd3232349fc733caef94034292 Mon Sep 17 00:00:00 2001 From: stainer_t Date: Mon, 8 Sep 2025 13:48:49 +0200 Subject: Initial commit from Polytechnique Montreal --- Dragon/src/BREANM.f | 340 ++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 340 insertions(+) create mode 100644 Dragon/src/BREANM.f (limited to 'Dragon/src/BREANM.f') diff --git a/Dragon/src/BREANM.f b/Dragon/src/BREANM.f new file mode 100644 index 0000000..99287ef --- /dev/null +++ b/Dragon/src/BREANM.f @@ -0,0 +1,340 @@ +*DECK BREANM + SUBROUTINE BREANM(IPMAC1,NG,LX1,NMIX1,IMIX,ICODE,ISPH,ZKEFF,B2, + 1 ENER,XXX1,VOL1,FLX1,DC1,TOT1,CHI1,SIGF1,SCAT1,JXM,JXP,FHETXM, + 2 FHETXP,ADF1,NGET,ADFREF,IPRINT) +* +*----------------------------------------------------------------------- +* +*Purpose: +* Implement the 1D DF-ANM reflector model. +* +*Copyright: +* Copyright (C) 2022 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 +* IPMAC1 nodal macrolib. +* NG number of energy groups. +* LX1 number of nodes in the reflector model. +* NMIX1 number of mixtures in the nodal calculation. +* IMIX mix index of each node. +* ICODE physical albedo index on each side of the domain. +* ISPH SPH flag (=0: use discontinuity factors; =1: use SPH factors). +* ZKEFF effective multiplication factor. +* B2 buckling. +* ENER energy limits. +* XXX1 spatial mesh. +* VOL1 volumes. +* FLX1 averaged fluxes +* DC1 diffusion coefficients. +* TOT1 total cross sections. +* CHI1 fission spectra. +* SIGF1 nu*fission cross sections. +* SCAT1 scattering P0 cross sections. +* JXM left boundary currents. +* JXP right boundary currents. +* FHETXM left boundary fluxes. +* FHETXP right boundary fluxes. +* ADF1 assembly discontinuity factors from macrolib. +* NGET type of NGET normalization if discontinuity factors +* (=0: simple; =1: imposed ADF on fuel assembly; =2: recover +* fuel assembly ADF from input macrolib). +* ADFREF imposed ADF values on fuel assembly side. +* IPRINT edition flag. +* +*----------------------------------------------------------------------- +* + USE GANLIB +*---- +* SUBROUTINE ARGUMENTS +*---- + TYPE(C_PTR) IPMAC1 + INTEGER NG,LX1,NMIX1,IMIX(LX1),ICODE(2),ISPH,NGET,IPRINT + REAL ZKEFF,B2,ENER(NG+1),XXX1(LX1+1),VOL1(NMIX1),FLX1(NMIX1,NG), + 1 DC1(NMIX1,NG),TOT1(NMIX1,NG),CHI1(NMIX1,NG),SIGF1(NMIX1,NG), + 2 SCAT1(NMIX1,NG,NG),JXM(NMIX1,NG),JXP(NMIX1,NG),FHETXM(NMIX1,NG), + 3 FHETXP(NMIX1,NG),ADF1(NMIX1,NG),ADFREF(NG) +*---- +* LOCAL VARIABLES +*---- + PARAMETER (NSTATE=40) + INTEGER ISTATE(NSTATE) + CHARACTER HADF*8 + TYPE(C_PTR) JPMAC1,KPMAC1 +*---- +* ALLOCATABLE ARRAYS +*---- + INTEGER, ALLOCATABLE, DIMENSION(:) :: IJJ,NJJ,IPOS + REAL, ALLOCATABLE, DIMENSION(:) :: WORK,AFACTOR,BETA,WORK1,WORK2, + 1 WORK4,WORK5,VOLTOT + REAL, ALLOCATABLE, DIMENSION(:,:) :: FDXM,FDXP,FHOMM,FHOMP,WORK3 + REAL(KIND=8), ALLOCATABLE, DIMENSION(:,:,:) :: L,R +*---- +* SCRATCH STORAGE ALLOCATION +*---- + ALLOCATE(L(NG,2*NG,LX1),R(NG,2*NG,LX1)) + ALLOCATE(FHOMM(NMIX1,NG),FHOMP(NMIX1,NG),FDXM(NMIX1,NG), + 1 FDXP(NMIX1,NG),AFACTOR(NG),BETA(NG),WORK1(NG),WORK2(NG), + 2 WORK3(NG,NG),WORK4(NG),WORK5(NG),VOLTOT(NMIX1)) +*---- +* COMPUTE BOUNDARY FLUXES +*---- + FDXM(:NMIX1,:NG)=0.0 + FDXP(:NMIX1,:NG)=0.0 + FHOMM(:NMIX1,:NG)=0.0 + FHOMP(:NMIX1,:NG)=0.0 + VOLTOT(:NMIX1)=0.0 + J_FUEL=0 + DO I=1,LX1 + IBM=IMIX(I) + IF(IBM.EQ.0) CYCLE + WORK1(:NG)=DC1(IBM,:NG) + WORK3(:NG,:NG)=SCAT1(IBM,:NG,:NG) + WORK4(:NG)=CHI1(IBM,:NG) + WORK5(:NG)=SIGF1(IBM,:NG) + DO IGR=1,NG + IF(SIGF1(IBM,IGR).GT.0.0) J_FUEL=I + WORK2(IGR)=TOT1(IBM,IGR)+B2*DC1(IBM,IGR)-SCAT1(IBM,IGR,IGR) + ENDDO + VOL=XXX1(I+1)-XXX1(I) + CALL NSSLR1(ZKEFF,NG,VOL,WORK1,WORK2,WORK3,WORK4,WORK5, + 1 L(1,1,I),R(1,1,I)) + ! + VOLTOT(IBM)=VOLTOT(IBM)+VOL + FHOMM(IBM,:NG)=FHOMM(IBM,:NG)+REAL(MATMUL(L(:NG,:NG,I), + 1 FLX1(IBM,:NG))+MATMUL(L(:NG,NG+1:2*NG,I),JXM(IBM,:NG)),4)*VOL + FHOMP(IBM,:NG)=FHOMP(IBM,:NG)+REAL(MATMUL(R(:NG,:NG,I), + 1 FLX1(IBM,:NG))+MATMUL(R(:NG,NG+1:2*NG,I),JXP(IBM,:NG)),4)*VOL + ENDDO + IF(IPRINT.GT.0) THEN + WRITE(6,'(/23H BREANM: SURFACE FLUXES)') + DO I=1,LX1 + IBM=IMIX(I) + IF(IBM.EQ.0) CYCLE + WRITE(6,'(/8H REGION=,I5)') I + WRITE(6,20) 'fluxm',REAL(MATMUL(L(:NG,:NG,I), + 1 FLX1(IBM,:NG))+MATMUL(L(:NG,NG+1:2*NG,I),JXM(IBM,:NG)),4) + WRITE(6,20) 'fluxp',REAL(MATMUL(R(:NG,:NG,I), + 1 FLX1(IBM,:NG))+MATMUL(R(:NG,NG+1:2*NG,I),JXP(IBM,:NG)),4) + ENDDO + ENDIF + DO IBM=1,NMIX1 + DO IGR=1,NG + FDXM(IBM,IGR)=VOLTOT(IBM)*FHETXM(IBM,IGR)/FHOMM(IBM,IGR) + FDXP(IBM,IGR)=VOLTOT(IBM)*FHETXP(IBM,IGR)/FHOMP(IBM,IGR) + ENDDO + ENDDO + IF(IPRINT.GT.0) THEN + WRITE(6,'(/48H BREANM: DISCONTINUITY FACTORS BEFORE NORMALIZAT, + 1 3HION)') + DO IBM=1,NMIX1 + WRITE(6,'(/9H MIXTURE=,I5)') IBM + WRITE(6,20) 'FDXM',FDXM(IBM,:NG) + WRITE(6,20) 'FDXP',FDXP(IBM,:NG) + ENDDO + ENDIF +*---- +* COMPUTE ALBEDOS +*---- + IF(ICODE(2).NE.0) THEN + BETA(:)=0.0 + IBM=IMIX(LX1) + DO IGR=1,NG + IF(IBM.EQ.0) CYCLE + AFACTOR(IGR)=FDXP(IBM,IGR)*JXP(IBM,IGR)/FHETXP(IBM,IGR) + BETA(IGR)=(1.0-2.0*AFACTOR(IGR))/(1.0+2.0*AFACTOR(IGR)) + ENDDO + IF(IPRINT.GT.0) THEN + WRITE(6,'(/16H BREANM: ALBEDOS)') + WRITE(6,20) 'BETA',BETA(:NG) + ENDIF + ENDIF +*---- +* NGET NORMALIZATION OF THE DISCONTINUITY FACTORS +*---- + IF(J_FUEL.GT.0) THEN + IF(NGET.GT.0) THEN + IBM=IMIX(J_FUEL) + DO IGR=1,NG + ! impose the adf on the fuel assembly side + IF(IBM.EQ.0) CYCLE + IF(NGET.EQ.1) THEN + FNORM=ADFREF(IGR)/FDXP(IBM,IGR) + ELSE + FNORM=ADF1(IBM,IGR)/FDXP(IBM,IGR) + ENDIF + FDXP(IBM,IGR)=FDXP(IBM,IGR)*FNORM + IF(J_FUEL1) THEN + IBMM=IMIX(J-1) + IF(IBMM.GT.0) FDXP(IBMM,IGR)=FDXP(IBMM,IGR)*FDXP(IBM,IGR)/ + 1 FDXM(IBM,IGR) + ENDIF + FDXM(IBM,IGR)=FDXP(IBM,IGR) + ENDDO + ENDDO + ENDIF + DO J=J_FUEL+1,LX1 + IBM=IMIX(J) + IF(IBM.EQ.0) CYCLE + DO IGR=1,NG + IF(J