Ability to read multiple secondary energy distributions on a single reaction.

This commit is contained in:
Paul Romano 2011-11-11 10:33:38 -05:00
parent ad2983918b
commit 6f73c0aaa0
3 changed files with 106 additions and 70 deletions

View file

@ -1,7 +1,8 @@
module cross_section
use constants
use cross_section_header, only: Nuclide, Reaction, SAB_Table, XsListing
use cross_section_header, only: Nuclide, Reaction, SAB_Table, XsListing, &
DistEnergy
use datatypes, only: dict_create, dict_add_key, dict_get_key, &
dict_has_key, dict_delete, dict_keys
use datatypes_header, only: DictionaryCI, ListKeyValueCI
@ -841,21 +842,11 @@ contains
type(Nuclide), pointer :: nuc
integer :: LED ! location of energy distribution locators
integer :: LDIS ! location of all energy distributions
integer :: LOCC ! location of energy distributions for given MT
integer :: LNW ! location of next energy distribution if multiple
integer :: LAW ! secondary energy distribution law
integer :: NR ! number of interpolation regions
integer :: NE ! number of incoming energies
integer :: IDAT ! location of first energy distribution for given MT
integer :: loc ! locator
integer :: length ! length of data to allocate
integer :: length_interp_data ! length of interpolation data
integer :: i ! loop index
type(Reaction), pointer :: rxn => null()
integer :: i ! loop index
type(Reaction), pointer :: rxn => null()
LED = JXS(10)
LDIS = JXS(11)
! Loop over all reactions
do i = 1, NXS(5)
@ -864,55 +855,94 @@ contains
! find location of energy distribution data
LOCC = XSS(LED + i - 1)
LNW = XSS(LDIS + LOCC - 1)
LAW = XSS(LDIS + LOCC)
IDAT = XSS(LDIS + LOCC + 1)
NR = XSS(LDIS + LOCC + 2)
rxn % edist % law = LAW
rxn % edist % n_interp = NR
! allocate space for ENDF interpolation parameters
if (NR > 0) then
allocate(rxn % edist % nbt(NR))
allocate(rxn % edist % int(NR))
end if
! read ENDF interpolation parameters
XSS_index = LDIS + LOCC + 3
if (NR > 0) then
rxn % edist % nbt = get_real(NR)
rxn % edist % int = get_real(NR)
end if
! allocate energy distribution
allocate(rxn % edist)
! allocate space for law validity data
NE = XSS(LDIS + LOCC + 3 + 2*NR)
allocate(rxn % edist % energy(NE))
allocate(rxn % edist % pvalid(NE))
length_interp_data = 5 + 2*(NR + NE)
! read law validity data
XSS_index = LDIS + LOCC + 4 + 2*NR
rxn % edist % energy = get_real(NE)
rxn % edist % pvalid = get_real(NE)
! Set index to beginning of IDAT array
loc = LDIS + IDAT - 2
! determine length of energy distribution
length = length_energy_dist(loc, LAW, LOCC, length_interp_data)
! allocate secondary energy distribution array
allocate(rxn % edist % data(length))
! read secondary energy distribution
XSS_index = loc + 1
rxn % edist % data = get_real(length)
! read data for energy distribution
call get_energy_dist(rxn % edist, LOCC)
end do
end subroutine read_energy_dist
!===============================================================================
! GET_ENERGY_DIST reads in data for a single law for an energy distribution and
! calls itself recursively if there are multiple energy distributions for a
! single reaction
!===============================================================================
recursive subroutine get_energy_dist(edist, loc_law)
type(DistEnergy), pointer :: edist ! energy distribution
integer, intent(in) :: loc_law ! locator for data
integer :: LDIS ! location of all energy distributions
integer :: LNW ! location of next energy distribution if multiple
integer :: LAW ! secondary energy distribution law
integer :: NR ! number of interpolation regions
integer :: NE ! number of incoming energies
integer :: IDAT ! location of first energy distribution for given MT
integer :: loc ! locator
integer :: length ! length of data to allocate
integer :: length_interp_data ! length of interpolation data
LDIS = JXS(11)
! locator for next law and information on this law
LNW = XSS(LDIS + loc_law - 1)
LAW = XSS(LDIS + loc_law)
IDAT = XSS(LDIS + loc_law + 1)
NR = XSS(LDIS + loc_law + 2)
edist % law = LAW
edist % n_interp = NR
! allocate space for ENDF interpolation parameters
if (NR > 0) then
allocate(edist % nbt(NR))
allocate(edist % int(NR))
end if
! read ENDF interpolation parameters
XSS_index = LDIS + loc_law + 3
if (NR > 0) then
edist % nbt = get_real(NR)
edist % int = get_real(NR)
end if
! allocate space for law validity data
NE = XSS(LDIS + loc_law + 3 + 2*NR)
edist % n_energy = NE
allocate(edist % energy(NE))
allocate(edist % pvalid(NE))
length_interp_data = 5 + 2*(NR + NE)
! read law validity data
XSS_index = LDIS + loc_law + 4 + 2*NR
edist % energy = get_real(NE)
edist % pvalid = get_real(NE)
! Set index to beginning of IDAT array
loc = LDIS + IDAT - 2
! determine length of energy distribution
length = length_energy_dist(loc, LAW, loc_law, length_interp_data)
! allocate secondary energy distribution array
allocate(edist % data(length))
! read secondary energy distribution
XSS_index = loc + 1
edist % data = get_real(length)
! read next energy distribution if present
if (LNW > 0) then
allocate(edist % next)
call get_energy_dist(edist % next, LNW)
end if
end subroutine get_energy_dist
!===============================================================================
! LENGTH_ENERGY_DIST determines how many values are contained in an LDAT energy
! distribution array based on the secondary energy law and location in XSS

View file

@ -31,6 +31,10 @@ module cross_section_header
real(8), allocatable :: energy(:) ! energy grid for law validity
real(8), allocatable :: pvalid(:) ! probability of law validity
real(8), allocatable :: data(:) ! energy distribution data
! For reactions that may have multiple energy distributions such as (n.2n),
! this pointer allows multiple laws to be stored
type(DistEnergy), pointer :: next => null()
end type DistEnergy
!===============================================================================
@ -39,15 +43,15 @@ module cross_section_header
!===============================================================================
type Reaction
integer :: MT ! ENDF MT value
real(8) :: Q_value ! Reaction Q value
integer :: TY ! Number of neutrons released
integer :: IE ! Starting energy grid index
real(8), allocatable :: sigma(:) ! Cross section values
logical :: has_angle_dist ! Angle distribution present?
logical :: has_energy_dist ! Energy distribution present?
type(DistAngle) :: adist ! Secondary angular distribution
type(DistEnergy) :: edist ! Secondary energy distribution
integer :: MT ! ENDF MT value
real(8) :: Q_value ! Reaction Q value
integer :: TY ! Number of neutrons released
integer :: IE ! Starting energy grid index
real(8), allocatable :: sigma(:) ! Cross section values
logical :: has_angle_dist ! Angle distribution present?
logical :: has_energy_dist ! Energy distribution present?
type(DistAngle) :: adist ! Secondary angular distribution
type(DistEnergy), pointer :: edist ! Secondary energy distribution
end type Reaction
!===============================================================================
@ -103,7 +107,7 @@ module cross_section_header
integer :: n_precursor
real(8), allocatable :: nu_d_data(:)
real(8), allocatable :: nu_d_precursor_data(:)
type(DistEnergy), allocatable :: nu_d_edist(:)
type(DistEnergy), pointer :: nu_d_edist(:) => null()
! Unresolved resonance data
logical :: urr_present

View file

@ -1025,7 +1025,8 @@ contains
real(8) :: yield ! delayed neutron precursor yield
real(8) :: prob ! cumulative probability
logical :: actual_event ! did fission actually occur? (no survival biasing)
type(Nuclide), pointer :: nuc
type(Nuclide), pointer :: nuc
type(DistEnergy), pointer :: edist => null()
! Get pointer to nuclide
nuc => nuclides(index_nuclide)
@ -1130,11 +1131,12 @@ contains
! sample from energy distribution for group j
law = nuc % nu_d_edist(j) % law
edist => nuc % nu_d_edist(j)
do
if (law == 44 .or. law == 61) then
call sample_energy(nuc%nu_d_edist(j), E, E_out, mu)
call sample_energy(edist, E, E_out, mu)
else
call sample_energy(nuc%nu_d_edist(j), E, E_out)
call sample_energy(edist, E, E_out)
end if
! resample if energy is >= 20 MeV
if (E_out < 20) exit
@ -1425,7 +1427,7 @@ contains
subroutine sample_energy(edist, E_in, E_out, mu_out, A, Q)
type(DistEnergy), intent(inout) :: edist
type(DistEnergy), pointer :: edist
real(8), intent(in) :: E_in
real(8), intent(out) :: E_out
real(8), intent(inout), optional :: mu_out