From e0f0e44c7a98bcf9e90836096682c97ef4f62ce0 Mon Sep 17 00:00:00 2001 From: ChasingNeutrons Date: Tue, 24 Feb 2026 19:54:38 +0000 Subject: [PATCH 1/2] Added fission * energy release as a cross section --- .../Reactions/reactionMG/fissionMG_class.f90 | 21 ++++++++++++++++- .../fissionCE_class.f90 | 17 ++++++++++++++ .../aceDatabase/aceNeutronDatabase_class.f90 | 3 +++ .../aceDatabase/aceNeutronNuclide_class.f90 | 23 +++++++++++++++++-- .../ceNeutronData/ceNeutronDatabase_inter.f90 | 6 +++++ .../baseMgNeutronMaterial_class.f90 | 2 ++ .../testNeutronDatabase_class.f90 | 11 ++++++++- .../xsPackages/neutronXsPackages_class.f90 | 12 ++++++++-- SharedModules/endfConstants.f90 | 8 ++++--- SharedModules/universalVariables.f90 | 3 +-- .../Tests/macroResponse_test.f90 | 13 +++++++++-- .../Tests/microResponse_test.f90 | 14 +++++++++-- .../TallyResponses/macroResponse_class.f90 | 3 +++ .../TallyResponses/microResponse_class.f90 | 3 +++ docs/Input Manual.rst | 7 +++--- 15 files changed, 128 insertions(+), 18 deletions(-) diff --git a/NuclearData/Reactions/reactionMG/fissionMG_class.f90 b/NuclearData/Reactions/reactionMG/fissionMG_class.f90 index 3d51ee384..47d41d89f 100644 --- a/NuclearData/Reactions/reactionMG/fissionMG_class.f90 +++ b/NuclearData/Reactions/reactionMG/fissionMG_class.f90 @@ -26,14 +26,17 @@ module fissionMG_class !! Fission spectrum is independent of the energy group !! !! Public members: - !! data -> data for fissionMG [energyGroup, reaction (nu or chi)] + !! data -> data for fissionMG [energyGroup, reaction (nu or chi)] + !! kappa -> Energy per fission. 200 MeV by default. !! !! Interface: !! reactionMG interface !! buildFromDict -> builds fissionMG from a SCONE dictionary + !! getKappa -> obtains the energy per fission in MeV !! type, public, extends(reactionMG) :: fissionMG real(defReal),dimension(:,:),allocatable :: data + real(defReal) :: kappa = 200.0_defReal contains ! Superclass procedures procedure :: init @@ -46,6 +49,7 @@ module fissionMG_class ! Local procedures procedure :: buildFromDict + procedure :: getKappa end type fissionMG @@ -216,6 +220,9 @@ subroutine buildFromDict(self, dict) ! Get number of groups call dict % get(nG, 'numberOfGroups') + + ! Get energy per fission + call dict % getOrDefault(self % kappa, 'kappa', 200.0_defReal) ! Allocate space allocate(self % data(nG, 2)) @@ -246,6 +253,18 @@ subroutine buildFromDict(self, dict) end subroutine buildFromDict + !! + !! Get the energy release from fission + !! + pure function getKappa(self) result(kappa) + class(fissionMG), intent(in) :: self + real(defReal) :: kappa + + kappa = self % kappa + + end function getKappa + + !! !! Cast reactionHandle pointer to fissionMG pointer !! diff --git a/NuclearData/Reactions/uncorrelatedReactionCE/fissionCE_class.f90 b/NuclearData/Reactions/uncorrelatedReactionCE/fissionCE_class.f90 index bb64e01a2..c782f9672 100644 --- a/NuclearData/Reactions/uncorrelatedReactionCE/fissionCE_class.f90 +++ b/NuclearData/Reactions/uncorrelatedReactionCE/fissionCE_class.f90 @@ -54,6 +54,7 @@ module fissionCE_class !! eLawPrompt -> Energy Law for the prompt fission neutrons !! nuBarDelayed -> Delayed release table for incindent energy [MeV] !! delayed -> Information about Delayed emission precursors + !! Q -> Q value for the reaction [MeV] !! !! Interface: !! uncorrelatedReactionCE interface @@ -65,6 +66,7 @@ module fissionCE_class class(energyLawENDF), allocatable :: eLawPrompt class(releaseLawENDF),allocatable :: nuBarDelayed type(precursor),dimension(:),allocatable :: delayed + real(defReal) :: Q = ZERO contains ! Superclass procedures @@ -79,6 +81,7 @@ module fissionCE_class ! Type specific procedures procedure :: buildFromACE + procedure :: getQ end type fissionCE contains @@ -148,8 +151,21 @@ elemental subroutine kill(self) deallocate(self % delayed) end if + self % Q = ZERO + end subroutine kill + !! + !! Returns the Q-value + !! + pure function getQ(self) result(Q) + class(fissionCE), intent(in) :: self + real(defReal) :: Q + + Q = self % Q + + end function getQ + !! !! Returns true if reaction is in Centre-Of-Mass frame !! @@ -349,6 +365,7 @@ subroutine buildFromACE(self, ACE, MT) ! Read basic data call new_totalNU(self % nuBarTotal, ACE) call new_energyLawENDF(self % eLawPrompt, ACE, MT) + self % Q = ACE % QforMT(MT) ! Read Delayed Data if (withDelayed) then diff --git a/NuclearData/ceNeutronData/aceDatabase/aceNeutronDatabase_class.f90 b/NuclearData/ceNeutronData/aceDatabase/aceNeutronDatabase_class.f90 index 6899dd076..f6d7c8528 100644 --- a/NuclearData/ceNeutronData/aceDatabase/aceNeutronDatabase_class.f90 +++ b/NuclearData/ceNeutronData/aceDatabase/aceNeutronDatabase_class.f90 @@ -854,6 +854,9 @@ subroutine init(self, dict, ptr, silent ) end if + ! Obtain the energy per fission with which to scale fission heating + call dict % getOrDefault(self % H235, 'energyPerFission', 202.27_defReal) + ! Get path to ACE library call dict % get(aceLibPath,'aceLibrary') diff --git a/NuclearData/ceNeutronData/aceDatabase/aceNeutronNuclide_class.f90 b/NuclearData/ceNeutronData/aceDatabase/aceNeutronNuclide_class.f90 index 2f1979edf..c92cbee1f 100644 --- a/NuclearData/ceNeutronData/aceDatabase/aceNeutronNuclide_class.f90 +++ b/NuclearData/ceNeutronData/aceDatabase/aceNeutronNuclide_class.f90 @@ -47,6 +47,7 @@ module aceNeutronNuclide_class integer(shortInt), parameter :: CAPTURE_XS = 4 integer(shortInt), parameter :: FISSION_XS = 5 integer(shortInt), parameter :: NU_FISSION = 6 + integer(shortInt), parameter :: KAPPA_XS = 7 !! @@ -436,9 +437,11 @@ elemental subroutine microXSs(self, xss, idx, f) if (self % isFissile()) then xss % fission = data(FISSION_XS, 2) * f + (ONE-f) * data(FISSION_XS, 1) xss % nuFission = data(NU_FISSION, 2) * f + (ONE-f) * data(NU_FISSION, 1) + xss % kappaXS = data(KAPPA_XS, 2) * f + (ONE-f) * data(KAPPA_XS, 1) else xss % fission = ZERO xss % nuFission = ZERO + xss % kappaXS = ZERO end if end associate @@ -483,9 +486,11 @@ subroutine getThXSs(self, xss, idx, f, E, kT, rand) if (self % isFissile()) then xss % fission = data(FISSION_XS, 2) * f + (ONE-f) * data(FISSION_XS, 1) xss % nuFission = data(NU_FISSION, 2) * f + (ONE-f) * data(NU_FISSION, 1) + xss % kappaXS = data(KAPPA_XS, 2) * f + (ONE-f) * data(KAPPA_XS, 1) else xss % fission = ZERO xss % nuFission = ZERO + xss % kappaXS = ZERO end if ! Read S(a,b) tables for elastic scatter: return zero if elastic scatter is off. @@ -555,9 +560,11 @@ subroutine getUrrXSs(self, xss, idx, f, E, xi) if (self % isFissile()) then xss % fission = data(FISSION_XS, 2) * f + (ONE-f) * data(FISSION_XS, 1) xss % nuFission = data(NU_FISSION, 2) * f + (ONE-f) * data(NU_FISSION, 1) + xss % kappaXS = data(KAPPA_XS, 2) * f + (ONE-f) * data(KAPPA_XS, 1) else xss % fission = ZERO xss % nuFission = ZERO + xss % kappaXS = ZERO end if ! Check if flag for multiplication factor (IFF) is true, and apply it to elastic scattering, @@ -587,6 +594,7 @@ subroutine getUrrXSs(self, xss, idx, f, E, xi) if (self % isFissile()) then xss % nuFission = xss % nuFission/xss % fission * val(3) + xss % kappaXS = xss % kappaXS / xss % fission * val(3) xss % fission = val(3) end if @@ -725,6 +733,7 @@ subroutine init(self, ACE, nucIdx, database) top, firstIdxMT4 real(defReal), dimension(:), allocatable :: xsMT4 type(stackInt) :: scatterMT, absMT + real(defReal) :: H_Q character(100), parameter :: Here = "init (aceNeutronNuclide_class.f90)" ! Reset nuclide just in case @@ -745,7 +754,7 @@ subroutine init(self, ACE, nucIdx, database) ! Allocate space for main XSs if (self % isFissile()) then - N = 6 + N = 7 else N = 4 end if @@ -807,10 +816,20 @@ subroutine init(self, ACE, nucIdx, database) call self % fission % init(ACE, N_f) end if - ! Calculate nuFission + ! Obtain Heating/Q scaling ratio + ! Check if database is associated in order to satisfy tests where it might not be! + if (associated(database)) then + H_Q = database % H235 / database % Q235 + else + H_Q = ONE + end if + + ! Calculate nuFission and kappaXS do i = bottom, Ngrid self % mainData(NU_FISSION,i) = self % mainData(FISSION_XS,i) * & self % fission % release(self % eGrid(i)) + self % mainData(KAPPA_XS,i) = self % mainData(FISSION_XS,i) * & + self % fission % getQ() * H_Q end do end if diff --git a/NuclearData/ceNeutronData/ceNeutronDatabase_inter.f90 b/NuclearData/ceNeutronData/ceNeutronDatabase_inter.f90 index c044c707e..ef368c55b 100644 --- a/NuclearData/ceNeutronData/ceNeutronDatabase_inter.f90 +++ b/NuclearData/ceNeutronData/ceNeutronDatabase_inter.f90 @@ -39,6 +39,8 @@ module ceNeutronDatabase_inter !! !! Public Members: !! mapDBRCnuc -> map to link indexes of DBRC nuclides with their corresponding 0K + !! H235 -> heating fission of U235, used for scaling fission energy + !! Q235 -> Q-value of U235 fission, used for scaling fission energy !! !! Interface: !! nuclearDatabase Interface @@ -53,6 +55,10 @@ module ceNeutronDatabase_inter !! type, public, abstract, extends(nuclearDatabase) :: ceNeutronDatabase type(intMap) :: mapDBRCnuc + + ! Default fission scaling [MeV] + real(defReal) :: H235 = 202.27_defReal + real(defReal) :: Q235 = 193.406_defReal contains diff --git a/NuclearData/mgNeutronData/baseMgNeutron/baseMgNeutronMaterial_class.f90 b/NuclearData/mgNeutronData/baseMgNeutron/baseMgNeutronMaterial_class.f90 index dc98bdc35..068a16549 100644 --- a/NuclearData/mgNeutronData/baseMgNeutron/baseMgNeutronMaterial_class.f90 +++ b/NuclearData/mgNeutronData/baseMgNeutron/baseMgNeutronMaterial_class.f90 @@ -129,9 +129,11 @@ subroutine getMacroXSs_byG(self, xss, G, rand) if(self % isFissile()) then xss % fission = self % data(FISSION_XS, G) xss % nuFission = self % data(NU_FISSION, G) + xss % kappaXS = xss % fission * self % fission % getKappa() else xss % fission = ZERO xss % nuFission = ZERO + xss % kappaXS = ZERO end if end subroutine getMacroXSs_byG diff --git a/NuclearData/testNeutronData/testNeutronDatabase_class.f90 b/NuclearData/testNeutronData/testNeutronDatabase_class.f90 index e6032e73b..dbe12118f 100644 --- a/NuclearData/testNeutronData/testNeutronDatabase_class.f90 +++ b/NuclearData/testNeutronData/testNeutronDatabase_class.f90 @@ -70,11 +70,12 @@ module testNeutronDatabase_class !! captureXS [in] -> Optional. Value of Capture XS !! fissionXS [in] -> Optional. Value of Fission XS !! nuFissionXS [in] -> Optional Value of nuFission + !! kappaXS [in] -> Optional Value of kappa XS !! !! Errors: !! None !! - subroutine build(self, xsVal, eScatterXS, ieScatterXS ,captureXS, fissionXS, nuFissionXS) + subroutine build(self, xsVal, eScatterXS, ieScatterXS ,captureXS, fissionXS, nuFissionXS, kappaXS) class(testNeutroNDatabase), intent(inout) :: self real(defReal), intent(in) :: xsVal real(defReal), intent(in),optional :: eScatterXS @@ -82,6 +83,7 @@ subroutine build(self, xsVal, eScatterXS, ieScatterXS ,captureXS, fissionXS, nuF real(defReal), intent(in),optional :: captureXS real(defReal), intent(in),optional :: fissionXS real(defReal), intent(in),optional :: nuFissionXS + real(defReal), intent(in),optional :: kappaXS self % xsVal = xsVal @@ -125,6 +127,13 @@ subroutine build(self, xsVal, eScatterXS, ieScatterXS ,captureXS, fissionXS, nuF else self % mat % xss % nuFission = xsVal end if + + ! kappa * fission + if(present(kappaXS)) then + self % mat % xss % kappaXS = kappaXS + else + self % mat % xss % kappaXS = xsVal + end if end subroutine build diff --git a/NuclearData/xsPackages/neutronXsPackages_class.f90 b/NuclearData/xsPackages/neutronXsPackages_class.f90 index efb3cf159..6ad3e9771 100644 --- a/NuclearData/xsPackages/neutronXsPackages_class.f90 +++ b/NuclearData/xsPackages/neutronXsPackages_class.f90 @@ -2,7 +2,6 @@ !! This module brakes standard rules !! It contains a library of XS Packages for Neutron particle type !! -!! module neutronXsPackages_class use numPrecision @@ -22,6 +21,7 @@ module neutronXsPackages_class !! capture -> sum of all reactions without secendary neutrons excluding fission [1/cm] !! fission -> total Fission MT=18 Cross-section [1/cm] !! nuFission -> total average neutron production Cross-section [1/cm] + !! kappaXS -> heating cross-section [MeV/cm] !! !! Interface: !! clean -> Set all XSs to 0.0 @@ -35,6 +35,7 @@ module neutronXsPackages_class real(defReal) :: capture = ZERO real(defReal) :: fission = ZERO real(defReal) :: nuFission = ZERO + real(defReal) :: kappaXS = ZERO contains procedure :: clean => clean_neutronMacroXSs procedure :: add => add_neutronMacroXSs @@ -54,6 +55,7 @@ module neutronXsPackages_class !! capture -> all reactions without secendary neutrons excluding fission [barn] !! fission -> total Fission MT=18 Cross-section [barn] !! nuFission -> total average neutron production Cross-section [barn] + !! kappaXS -> heating Cross-section [MeV*barn] !! type, public :: neutronMicroXSs real(defReal) :: total = ZERO @@ -62,6 +64,7 @@ module neutronXsPackages_class real(defReal) :: capture = ZERO real(defReal) :: fission = ZERO real(defReal) :: nuFission = ZERO + real(defReal) :: kappaXS = ZERO contains procedure :: invert => invert_microXSs end type neutronMicroXSs @@ -88,6 +91,7 @@ elemental subroutine clean_neutronMacroXSs(self) self % capture = ZERO self % fission = ZERO self % nuFission = ZERO + self % kappaXS = ZERO end subroutine clean_neutronMacroXSs @@ -114,6 +118,7 @@ elemental subroutine add_neutronMacroXSs(self, micro, dens) self % capture = self % capture + dens * micro % capture self % fission = self % fission + dens * micro % fission self % nuFission = self % nuFission + dens * micro % nuFission + self % kappaXS = self % kappaXS + dens * micro % kappaXS end subroutine add_neutronMacroXSs @@ -159,6 +164,9 @@ elemental function get(self, MT) result(xs) case(macroNuFission) xs = self % nuFission + case(macroKappaFission) + xs = self % kappaXS + case(macroAbsorbtion) xs = self % fission + self % capture @@ -233,7 +241,7 @@ end function invert_macroXSs !! !! Use a real r in <0;1> to sample reaction from Microscopic XSs !! - !! This function involves a bit of code so is written for conviniance + !! This function involves a bit of code so is written for convenience !! !! Args: !! r [in] -> Real number in <0;1> diff --git a/SharedModules/endfConstants.f90 b/SharedModules/endfConstants.f90 index e5ed867a2..97d68b000 100644 --- a/SharedModules/endfConstants.f90 +++ b/SharedModules/endfConstants.f90 @@ -87,6 +87,7 @@ module endfConstants N_Xt = 205 ,& N_X3He = 206 ,& N_Xa = 207 ,& + N_KAPPA = 301 ,& nubar_tot = 452 ,& nubar_del = 455 ,& nubar_prompt = 456 ,& @@ -114,7 +115,8 @@ module endfConstants macroEscatter = -3 ,& macroIEscatter = -4 ,& macroFission = -6 ,& - macroNuFission = -7 + macroNuFission = -7 ,& + macroKappaFission = -80 ! List of Macro MT numbers for macroscopic XSs. Unique to SCONE (not from Serpent) integer(shortInt), parameter :: macroAllScatter = -20 ,& @@ -129,10 +131,10 @@ module endfConstants 36, 37, 38, 41, 42, 44, 45, 91, 101, 102, 103, & 104, 105, 106, 107, 108, 109, 111, 112, 113, & 114, 115, 116, 117, 203, 204, 205, 206, 207, & - 891, N_Nl, N_2Nl] + 301, 891, N_Nl, N_2Nl] integer(shortInt), dimension(*), parameter :: availableMacroMTs = & - [-1, -2, -3, -4, -6, -7, -20, -21, -22] + [-1, -2, -3, -4, -6, -7, -20, -21, -22, -80] diff --git a/SharedModules/universalVariables.f90 b/SharedModules/universalVariables.f90 index 790fff3e5..14a6d15b8 100644 --- a/SharedModules/universalVariables.f90 +++ b/SharedModules/universalVariables.f90 @@ -80,8 +80,7 @@ module universalVariables ! Neutron mass and speed of light in vacuum from from https://physics.nist.gov/cuu/Constants/index.html real(defReal), parameter :: neutronMass = 939.56542194_defReal, & ! Neutron mass in MeV (m*c^2) lightSpeed = 2.99792458e10_defReal, & ! Light speed in cm/s - kBoltzmann = 1.380649e-23_defReal, & ! Bolztmann constant in J/K - energyPerFission = 200.0_defReal ! MeV + kBoltzmann = 1.380649e-23_defReal ! Bolztmann constant in J/K ! Unit conversion real(defReal), parameter :: joulesPerMeV = 1.60218e-13_defReal ,& ! Convert MeV to J diff --git a/Tallies/TallyResponses/Tests/macroResponse_test.f90 b/Tallies/TallyResponses/Tests/macroResponse_test.f90 index bfa87292e..103c61007 100644 --- a/Tallies/TallyResponses/Tests/macroResponse_test.f90 +++ b/Tallies/TallyResponses/Tests/macroResponse_test.f90 @@ -18,6 +18,7 @@ module macroResponse_test type(macroResponse) :: response_fission type(macroResponse) :: response_nuFission type(macroResponse) :: response_absorbtion + type(macroResponse) :: response_kappa type(macroResponse) :: response_MT type(testNeutronDatabase) :: xsData contains @@ -36,8 +37,8 @@ subroutine setUp(this) type(dictionary) :: tempDict ! Allocate and initialise test nuclearData - ! Cross-sections: Total eScattering IeScatter Capture Fission nuFission - call this % xsData % build(6.0_defReal, 3.0_defReal, ZERO, 2.0_defReal, 1.0_defReal, 1.5_defReal) + ! Cross-sections: Total eScattering IeScatter Capture Fission nuFission kappFission + call this % xsData % build(6.0_defReal, 3.0_defReal, ZERO, 2.0_defReal, 1.0_defReal, 1.5_defReal, 9.0_defReal) ! Set up responses ! Total @@ -75,6 +76,13 @@ subroutine setUp(this) call this % response_absorbtion % init(tempDict) call tempDict % kill() + ! Energy deposition + call tempDict % init(2) + call tempDict % store('type','macroResponse') + call tempDict % store('MT', macroKappaFission) + call this % response_kappa % init(tempDict) + call tempDict % kill() + ! MT number call tempDict % init(2) call tempDict % store('type','macroResponse') @@ -117,6 +125,7 @@ subroutine testGettingResponse(this) @assertEqual(1.0_defReal, this % response_fission % get(p, this % xsData), TOL) @assertEqual(1.5_defReal, this % response_nuFission % get(p, this % xsData), TOL) @assertEqual(3.0_defReal, this % response_absorbtion % get(p, this % xsData), TOL) + @assertEqual(9.0_defReal, this % response_kappa % get(p, this % xsData), TOL) @assertEqual(105.0_defReal, this % response_MT % get(p, this % xsData), TOL) p % isMG = .true. diff --git a/Tallies/TallyResponses/Tests/microResponse_test.f90 b/Tallies/TallyResponses/Tests/microResponse_test.f90 index ec1e810f9..8f96e8adf 100644 --- a/Tallies/TallyResponses/Tests/microResponse_test.f90 +++ b/Tallies/TallyResponses/Tests/microResponse_test.f90 @@ -19,6 +19,7 @@ module microResponse_test type(microResponse) :: response_capture type(microResponse) :: response_fission type(microResponse) :: response_absorbtion + type(microResponse) :: response_kappa type(microResponse) :: response_MT type(testNeutronDatabase) :: xsData contains @@ -38,8 +39,8 @@ subroutine setUp(this) ! Allocate and initialise test nuclearData - ! Cross-sections: Total eScattering IeScatter Capture Fission nuFission - call this % xsData % build(6.0_defReal, 3.0_defReal, ZERO, 2.0_defReal, 1.0_defReal, 1.5_defReal) + ! Cross-sections: Total eScattering IeScatter Capture Fission nuFission kappFission + call this % xsData % build(6.0_defReal, 3.0_defReal, ZERO, 2.0_defReal, 1.0_defReal, 1.5_defReal, 9.0_defReal) ! Set dictionaries to initialise material call dictMat1 % init(1) @@ -97,6 +98,14 @@ subroutine setUp(this) call this % response_absorbtion % init(tempDict) call tempDict % kill() + ! Energy deposition + call tempDict % init(3) + call tempDict % store('type', 'microResponse') + call tempDict % store('MT', N_KAPPA) + call tempDict % store('material', 'Xenon') + call this % response_kappa % init(tempDict) + call tempDict % kill() + ! MT number call tempDict % init(2) call tempDict % store('type','microResponse') @@ -140,6 +149,7 @@ subroutine testGettingResponse(this) @assertEqual(0.5_defReal, this % response_fission % get(p, this % xsData), TOL) @assertEqual(1.5_defReal, this % response_eScatter % get(p, this % xsData), TOL) @assertEqual(1.5_defReal, this % response_absorbtion % get(p, this % xsData), TOL) + @assertEqual(4.5_defReal, this % response_kappa % get(p, this % xsData), TOL) @assertEqual(8.0_defReal, this % response_MT % get(p, this % xsData), TOL) p % isMG = .true. diff --git a/Tallies/TallyResponses/macroResponse_class.f90 b/Tallies/TallyResponses/macroResponse_class.f90 index 44ba0c9f8..d5d488d84 100644 --- a/Tallies/TallyResponses/macroResponse_class.f90 +++ b/Tallies/TallyResponses/macroResponse_class.f90 @@ -115,6 +115,9 @@ subroutine build(self, MT) case(N_ABSORPTION) self % MT = macroAbsorbtion + case(N_KAPPA) + self % MT = macroKappaFission + case default self % mainData = .false. self % MT = MT diff --git a/Tallies/TallyResponses/microResponse_class.f90 b/Tallies/TallyResponses/microResponse_class.f90 index fd7c63808..9bd5bb72f 100644 --- a/Tallies/TallyResponses/microResponse_class.f90 +++ b/Tallies/TallyResponses/microResponse_class.f90 @@ -143,6 +143,9 @@ subroutine build(self, MT) case(N_ABSORPTION) self % MT = macroAbsorbtion + case(N_KAPPA) + self % MT = macroKappaFission + case default self % MT = MT self % mainData = .false. diff --git a/docs/Input Manual.rst b/docs/Input Manual.rst index cc4ceb1bd..874d54cb2 100644 --- a/docs/Input Manual.rst +++ b/docs/Input Manual.rst @@ -1370,8 +1370,9 @@ Example: :: * macroResponse: used to score macroscopic reaction rates - MT: MT number of the desired reaction. The options are: -1 (total), -2 (disappearance), - -3 (elastic scattering), -4 (total inelastic scattering), -6 (fission), -7 nu*fission), - -20 (total scattering), -21 (absorption), -22 (total non elastic, i.e., absorption + inelastic). + -3 (elastic scattering), -4 (total inelastic scattering), -6 (fission), -7 (nu*fission), + -20 (total scattering), -21 (absorption), -22 (total non elastic, i.e., absorption + inelastic), + -80 (kappa*fission). Additionally, all the MT numbers allowed by microResponse can be used here. Example: :: @@ -1385,7 +1386,7 @@ Example: :: * microResponse: used to score microscopic reaction rates - MT: MT number of the desired reaction. The options are: 1, 2, 3, 4, 5, 11, 16-25, 27-30, - 32-38, 41, 42, 44, 45, 51-90, 91, 101-109, 111-117, 203-207, 875-890. These MT numbers + 32-38, 41, 42, 44, 45, 51-90, 91, 101-109, 111-117, 203-207, 301, 875-890. These MT numbers are defined in the conventional way, i.e., following the ENDF standard - material: material name where to score the reaction. The material must be defined to include only one nuclide; its density could be anything, it doesn't From c3683540a978082e3098b3a62c09d3201a0ea2b1 Mon Sep 17 00:00:00 2001 From: ChasingNeutrons Date: Tue, 3 Mar 2026 20:05:31 +0000 Subject: [PATCH 2/2] Fixes --- .../Reactions/reactionMG/fissionMG_class.f90 | 6 ++- SharedModules/endfConstants.f90 | 2 +- Tallies/TallyClerks/Tests/mgXsClerk_test.f90 | 18 ++++--- Tallies/TallyClerks/mgXsClerk_class.f90 | 49 +++++++++++++------ .../Tests/macroResponse_test.f90 | 2 +- .../Tests/microResponse_test.f90 | 2 +- 6 files changed, 50 insertions(+), 29 deletions(-) diff --git a/NuclearData/Reactions/reactionMG/fissionMG_class.f90 b/NuclearData/Reactions/reactionMG/fissionMG_class.f90 index 47d41d89f..ed7bd9d9c 100644 --- a/NuclearData/Reactions/reactionMG/fissionMG_class.f90 +++ b/NuclearData/Reactions/reactionMG/fissionMG_class.f90 @@ -18,6 +18,8 @@ module fissionMG_class !! public :: fissionMG_TptrCast + real(defReal), private, parameter :: KAPPA_DEFAULT = 202.27 ! [MeV] + !! !! A special type of MG reaction that contains data related to fission !! @@ -36,7 +38,7 @@ module fissionMG_class !! type, public, extends(reactionMG) :: fissionMG real(defReal),dimension(:,:),allocatable :: data - real(defReal) :: kappa = 200.0_defReal + real(defReal) :: kappa = KAPPA_DEFAULT contains ! Superclass procedures procedure :: init @@ -222,7 +224,7 @@ subroutine buildFromDict(self, dict) call dict % get(nG, 'numberOfGroups') ! Get energy per fission - call dict % getOrDefault(self % kappa, 'kappa', 200.0_defReal) + call dict % getOrDefault(self % kappa, 'kappa', KAPPA_DEFAULT) ! Allocate space allocate(self % data(nG, 2)) diff --git a/SharedModules/endfConstants.f90 b/SharedModules/endfConstants.f90 index ff3b1ad34..de8ee3dba 100644 --- a/SharedModules/endfConstants.f90 +++ b/SharedModules/endfConstants.f90 @@ -120,7 +120,7 @@ module endfConstants macroFission = -6 ,& macroNuFission = -7 ,& macroPromptNuFission = -8 ,& - macroDelayedNuFission = -9 + macroDelayedNuFission = -9 ,& macroKappaFission = -80 ! List of Macro MT numbers for macroscopic XSs. Unique to SCONE (not from Serpent) diff --git a/Tallies/TallyClerks/Tests/mgXsClerk_test.f90 b/Tallies/TallyClerks/Tests/mgXsClerk_test.f90 index d14e284c4..d193cdb13 100644 --- a/Tallies/TallyClerks/Tests/mgXsClerk_test.f90 +++ b/Tallies/TallyClerks/Tests/mgXsClerk_test.f90 @@ -74,7 +74,7 @@ subroutine setUp(this) ! Build test neutronDatabase call this % nucData % build(ONE, captureXS = 2.0_defReal, & - fissionXS = 1.5_defReal, nuFissionXS = 3.0_defReal) + fissionXS = 1.5_defReal, nuFissionXS = 3.0_defReal, kappaXS = 300.0_defReal) end subroutine setUp @@ -105,8 +105,8 @@ subroutine testScoring_clerk1(this) type(particleState) :: pFiss type(outputFile) :: out real(defReal), dimension(:,:), allocatable :: fiss, capt, transFL, transOS, & - nu, chi, P0, P1, P2, P3, P4, & - P5, P6, P7, prod + nu, chi, kappa, P0, P1, P2, P3, & + P4, P5, P6, P7, prod real(defReal), parameter :: TOL = 1.0E-9 ! Configure memory @@ -143,7 +143,7 @@ subroutine testScoring_clerk1(this) call mem % closeCycle(ONE) ! Process and get results - call this % clerk_test1 % processRes(mem, capt, fiss, transFL, transOS, nu, chi, P0, P1, prod) + call this % clerk_test1 % processRes(mem, capt, fiss, transFL, transOS, nu, chi, kappa, P0, P1, prod) call this % clerk_test1 % processPN(mem, P2, P3, P4, P5, P6, P7) ! Verify results of scoring @@ -151,6 +151,7 @@ subroutine testScoring_clerk1(this) @assertEqual([ZERO, ZERO, 1.5_defReal, ZERO], fiss(1,:), TOL, 'Fission XS' ) @assertEqual([ZERO, ZERO, TWO, ZERO], nu(1,:), TOL, 'NuFission XS' ) @assertEqual([ZERO, ZERO, HALF, HALF], chi(1,:), TOL, 'Chi' ) + @assertEqual([ZERO, ZERO, 300.0_defReal, ZERO], kappa(1,:), TOL, 'Kappa XS' ) @assertEqual([ZERO, ZERO, 4.0_defReal, ZERO], transOS(1,:), TOL, 'Transport XS O.S.' ) @assertEqual([ZERO, ZERO, 5.5_defReal, ZERO], transFL(1,:), TOL, 'Transport XS F.L.' ) @assertEqual([ZERO, ZERO, ZERO, ZERO, ZERO, TWO, ZERO, ZERO], P0(1,:), TOL, 'P0' ) @@ -165,7 +166,7 @@ subroutine testScoring_clerk1(this) @assertEqual([ZERO, ZERO, ZERO, ZERO, ZERO, -0.0683670044_defReal, ZERO, ZERO], P7(1,:), TOL, 'P7' ) ! Test getting size - @assertEqual(100, this % clerk_test1 % getSize(),'Test getSize():') + @assertEqual(104, this % clerk_test1 % getSize(),'Test getSize():') ! Test correctness of output calls call out % init('dummyPrinter', fatalErrors = .false.) @@ -185,7 +186,7 @@ subroutine testScoring_clerk2(this) type(particleState) :: pFiss type(outputFile) :: out real(defReal), dimension(:,:), allocatable :: fiss, capt, transFL, transOS, & - nu, chi, P0, P1, prod + nu, chi, kappa, P0, P1, prod real(defReal), parameter :: TOL = 1.0E-9 ! Configure memory @@ -221,13 +222,14 @@ subroutine testScoring_clerk2(this) call mem % closeCycle(ONE) ! Process and get results - call this % clerk_test2 % processRes(mem, capt, fiss, transFL, transOS, nu, chi, P0, P1, prod) + call this % clerk_test2 % processRes(mem, capt, fiss, transFL, transOS, nu, chi, kappa, P0, P1, prod) ! Verify results of scoring @assertEqual([ZERO, ZERO, TWO], capt(1,:), TOL, 'Capture XS' ) @assertEqual([ZERO, ZERO, 1.5_defReal], fiss(1,:), TOL, 'Fission XS' ) @assertEqual([ZERO, ZERO, TWO], nu(1,:), TOL, 'NuFission XS' ) @assertEqual([HALF, ZERO, HALF], chi(1,:), TOL, 'Chi' ) + @assertEqual([ZERO, ZERO, 300.0_defReal], kappa(1,:), TOL, 'Kappa XS' ) @assertEqual([ZERO, ZERO, 4.0_defReal], transOS(1,:), TOL, 'Transport XS O.S.' ) @assertEqual([ZERO, ZERO, 5.5_defReal], transFL(1,:), TOL, 'Transport XS F.L.' ) @assertEqual([ZERO, ZERO, ZERO, ZERO, ZERO, ZERO, ZERO, TWO, ZERO], P0(1,:), TOL, 'P0' ) @@ -235,7 +237,7 @@ subroutine testScoring_clerk2(this) @assertEqual([ONE, ONE, ONE, ONE, ONE, ONE, ONE, TWO, ONE], prod(1,:), TOL, 'prod' ) ! Test getting size - @assertEqual(48, this % clerk_test2 % getSize(),'Test getSize():') + @assertEqual(51, this % clerk_test2 % getSize(),'Test getSize():') ! Test correctness of output calls call out % init('dummyPrinter', fatalErrors = .false.) diff --git a/Tallies/TallyClerks/mgXsClerk_class.f90 b/Tallies/TallyClerks/mgXsClerk_class.f90 index a65d4aaf0..567cbc1e6 100644 --- a/Tallies/TallyClerks/mgXsClerk_class.f90 +++ b/Tallies/TallyClerks/mgXsClerk_class.f90 @@ -29,7 +29,7 @@ module mgXsClerk_class private !! Size of clerk memory - integer(shortInt), parameter :: ARRAY_SCORE_SIZE = 7 ,& ! Size of data to store as 1D arrays + integer(shortInt), parameter :: ARRAY_SCORE_SIZE = 8 ,& ! Size of data to store as 1D arrays MATRIX_SCORE_SMALL = 3 ,& ! Size of data to store as 2D arrays when scoring up to P1 MATRIX_SCORE_FULL = 9 ! Size of data to store as 2D arrays when scoring up to P7 @@ -40,14 +40,16 @@ module mgXsClerk_class FISS_idx = 4 ,& ! Fission macroscopic reaction rate NUBAR_idx = 5 ,& ! NuBar CHI_idx = 6 ,& ! Fission neutron spectrum - SCATT_EV_idx = 7 ! Analog: number of scattering events + KAPPAFISS_idx = 7 ,& ! Fission heating macroscopic reaction rate + SCATT_EV_idx = 8 ! Analog: number of scattering events !! !! Multi-group macroscopic cross section calculation !! !! It prints out: - !! capture xs, fission xs, transport xs, nu, chi, the P0 and P1 scattering matrices, - !! the P0 scattering production matrix. On request, also the P2 -> P7 scattering matrices. + !! capture xs, fission xs, transport xs, nu, chi, kappa * fission xs, the P0 and P1 scattering + !! matrices, the P0 scattering production matrix. On request, also the P2 -> P7 scattering + !! matrices. !! !! NOTE: !! - the cross sections are tallied with a collision estimator; @@ -230,7 +232,7 @@ subroutine reportInColl(self, p, xsData, mem, virtual) type(particleState) :: state type(neutronMacroXSs) :: xss class(neutronMaterial), pointer :: mat - real(defReal) :: nuFissXS, captXS, fissXS, scattXS, flux + real(defReal) :: nuFissXS, captXS, fissXS, scattXS, kappaXS, flux integer(shortInt) :: enIdx, locIdx, binIdx integer(longInt) :: addr character(100), parameter :: Here =' reportInColl (mgXsClerk_class.f90)' @@ -288,6 +290,7 @@ subroutine reportInColl(self, p, xsData, mem, virtual) captXS = xss % capture * flux fissXS = xss % fission * flux scattXS = (xss % elasticScatter + xss % inelasticScatter) * flux + kappaXS = xss % kappaXS * flux else @@ -296,6 +299,7 @@ subroutine reportInColl(self, p, xsData, mem, virtual) captXS = ZERO fissXS = ZERO scattXS = ZERO + kappaXS = ZERO end if @@ -305,6 +309,7 @@ subroutine reportInColl(self, p, xsData, mem, virtual) call mem % score(captXS, addr + CAPT_idx) call mem % score(fissXS, addr + FISS_idx) call mem % score(scattXS, addr + SCATT_idx) + call mem % score(kappaXS, addr + KAPPAFISS_idx) end subroutine reportInColl @@ -500,7 +505,7 @@ end subroutine reportSpawn !! none !! pure subroutine processRes(self, mem, capt_res, fiss_res, transFL_res, transOS_res, & - nu_res, chi_res, P0_res, P1_res, prod_res) + nu_res, chi_res, kappa_res, P0_res, P1_res, prod_res) class(mgXsClerk), intent(in) :: self type(scoreMemory), intent(in) :: mem real(defReal), dimension(:,:), allocatable, intent(out) :: capt_res @@ -509,6 +514,7 @@ pure subroutine processRes(self, mem, capt_res, fiss_res, transFL_res, transOS_r real(defReal), dimension(:,:), allocatable, intent(out) :: transOS_res real(defReal), dimension(:,:), allocatable, intent(out) :: nu_res real(defReal), dimension(:,:), allocatable, intent(out) :: chi_res + real(defReal), dimension(:,:), allocatable, intent(out) :: kappa_res real(defReal), dimension(:,:), allocatable, intent(out) :: P0_res real(defReal), dimension(:,:), allocatable, intent(out) :: P1_res real(defReal), dimension(:,:), allocatable, intent(out) :: prod_res @@ -518,7 +524,7 @@ pure subroutine processRes(self, mem, capt_res, fiss_res, transFL_res, transOS_r integer(shortInt) :: N, M, i, j, k, g1, gEnd, idx real(defReal) :: capt, fiss, scatt, nu, chi, P0, P1, prod, sumChi, flux, scattProb, & captStd, fissStd, scattStd, nuStd, chiStd, P0std, P1std, prodStd, & - fluxStd, scattProbStd, scattXS, scattXSstd + fluxStd, scattProbStd, scattXS, scattXSstd, kappa, kappaStd ! Get number of bins N = self % energyN @@ -526,9 +532,9 @@ pure subroutine processRes(self, mem, capt_res, fiss_res, transFL_res, transOS_r ! Allocate arrays for MG xss allocate( capt_res(2,N*M), fiss_res(2,N*M), transFL_res(2,N*M), transOS_res(2,N*M), & - nu_res(2,N*M), chi_res(2,N*M), P0_res(2,N*N*M), P1_res(2,N*N*M), & - prod_res(2,N*N*M), tot(N*M), fluxG(N*M), delta(M, N), totStd(N*M), & - fluxGstd(N*M), deltaStd(M, N) ) + nu_res(2,N*M), chi_res(2,N*M), kappa_res(2,N*M), P0_res(2,N*N*M), & + P1_res(2,N*N*M), prod_res(2,N*N*M), tot(N*M), fluxG(N*M), delta(M, N), & + totStd(N*M), fluxGstd(N*M), deltaStd(M, N) ) ! Initialise values sumChi = 0 ! to normalise chi @@ -548,6 +554,7 @@ pure subroutine processRes(self, mem, capt_res, fiss_res, transFL_res, transOS_r call mem % getResult(nu, nuStd, addr + NUBAR_idx) call mem % getResult(chi, chiStd, addr + CHI_idx) call mem % getResult(scattProb, scattProbStd, addr + SCATT_EV_idx) + call mem % getResult(kappa, kappaStd, addr + KAPPAFISS_idx) ! Calculate MG constants, being careful to avoid division by zero ! If flux is zero all cross sections must be set to zero @@ -564,13 +571,16 @@ pure subroutine processRes(self, mem, capt_res, fiss_res, transFL_res, transOS_r ! Calculate fission production term and uncertainties if (fiss == ZERO) then - fiss_res(1:2,i) = ZERO - nu_res(1:2,i) = ZERO + fiss_res(1:2,i) = ZERO + nu_res(1:2,i) = ZERO + kappa_res(1:2,i) = ZERO else fiss_res(1,i) = fiss/flux fiss_res(2,i) = fiss_res(1,i) * sqrt((fissStd/fiss)**2 + (fluxStd/flux)**2) nu_res(1,i) = nu/fiss nu_res(2,i) = nu_res(1,i) * sqrt((nuStd/nu)**2 + (fissStd/fiss)**2) + kappa_res(1,i) = kappa/flux + kappa_res(2,i) = kappa_res(1,i) * sqrt((kappaStd/kappa)**2 + (fluxStd/flux)**2) end if ! Store fission spectrum @@ -758,8 +768,8 @@ subroutine print(self, outFile, mem) type(scoreMemory), intent(in) :: mem integer(shortInt),dimension(:),allocatable :: resArrayShape real(defReal), dimension(:,:), allocatable :: fiss, capt, transFL, transOS, & - nu, chi, P0, P1, P2, P3, P4, & - P5, P6, P7, prod + nu, chi, kappa, P0, P1, P2, P3, & + P4, P5, P6, P7, prod character(nameLen) :: name integer(shortInt) :: i @@ -784,7 +794,7 @@ subroutine print(self, outFile, mem) end if ! Process and get results - call self % processRes(mem, capt, fiss, transFL, transOS, nu, chi, P0, P1, prod) + call self % processRes(mem, capt, fiss, transFL, transOS, nu, chi, kappa, P0, P1, prod) ! Print results name = 'capture' @@ -801,6 +811,13 @@ subroutine print(self, outFile, mem) end do call outFile % endArray() + name = 'kappaFission' + call outFile % startArray(name, resArrayShape) + do i=1,product(resArrayShape) + call outFile % addResult(kappa(1,i),kappa(2,i)) + end do + call outFile % endArray() + name = 'transportFluxLimited' call outFile % startArray(name, resArrayShape) do i=1,product(resArrayShape) @@ -853,7 +870,7 @@ subroutine print(self, outFile, mem) call outFile % endArray() ! Deallocate to limit memory consumption when writing to the output file - deallocate(capt, fiss, transFL, transOS, nu, chi, P0, P1, prod) + deallocate(capt, fiss, transFL, transOS, nu, chi, kappa, P0, P1, prod) ! If high order scattering is requested, print the other matrices if (self % PN) then diff --git a/Tallies/TallyResponses/Tests/macroResponse_test.f90 b/Tallies/TallyResponses/Tests/macroResponse_test.f90 index 103c61007..69c3815f0 100644 --- a/Tallies/TallyResponses/Tests/macroResponse_test.f90 +++ b/Tallies/TallyResponses/Tests/macroResponse_test.f90 @@ -37,7 +37,7 @@ subroutine setUp(this) type(dictionary) :: tempDict ! Allocate and initialise test nuclearData - ! Cross-sections: Total eScattering IeScatter Capture Fission nuFission kappFission + ! Cross-sections: Total eScattering IeScatter Capture Fission nuFission kappaFission call this % xsData % build(6.0_defReal, 3.0_defReal, ZERO, 2.0_defReal, 1.0_defReal, 1.5_defReal, 9.0_defReal) ! Set up responses diff --git a/Tallies/TallyResponses/Tests/microResponse_test.f90 b/Tallies/TallyResponses/Tests/microResponse_test.f90 index 8f96e8adf..462fac181 100644 --- a/Tallies/TallyResponses/Tests/microResponse_test.f90 +++ b/Tallies/TallyResponses/Tests/microResponse_test.f90 @@ -39,7 +39,7 @@ subroutine setUp(this) ! Allocate and initialise test nuclearData - ! Cross-sections: Total eScattering IeScatter Capture Fission nuFission kappFission + ! Cross-sections: Total eScattering IeScatter Capture Fission nuFission kappaFission call this % xsData % build(6.0_defReal, 3.0_defReal, ZERO, 2.0_defReal, 1.0_defReal, 1.5_defReal, 9.0_defReal) ! Set dictionaries to initialise material