crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-26 13:15:39 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions

Merge pull request #955 from paulromano/xs-improvements

Browse source 
Small optimizations for cross section lookups
This commit is contained in:
Sterling Harper 2018-01-18 14:30:53 -05:00 • committed by GitHub
commit 72e6ecbc1d
No known key found for this signature in database
GPG key ID: 4AEE18F83AFDEB23
5 changed files with 106 additions and 63 deletions
Download patch file Download diff file Expand all files Collapse all files

4
CMakeLists.txt
View file

@ -116,9 +116,9 @@ if(CMAKE_Fortran_COMPILER_ID STREQUAL GNU)
endif()
# GCC compiler options
list(APPEND f90flags -cpp -std=f2008ts -fbacktrace -O2)
list(APPEND f90flags -cpp -std=f2008ts -fbacktrace -O2 -fstack-arrays)
if(debug)
list(REMOVE_ITEM f90flags -O2)
list(REMOVE_ITEM f90flags -O2 -fstack-arrays)
list(APPEND f90flags -g -Wall -Wno-unused-dummy-argument -pedantic
-fbounds-check -ffpe-trap=invalid,overflow,underflow)
list(APPEND ldflags -g)

66
src/cross_section.F90
View file

@ -42,7 +42,6 @@ contains
! Set all material macroscopic cross sections to zero
material_xs % total = ZERO
material_xs % elastic = ZERO
material_xs % absorption = ZERO
material_xs % fission = ZERO
material_xs % nu_fission = ZERO
@ -115,10 +114,6 @@ contains
material_xs % total = material_xs % total + &
atom_density * micro_xs(i_nuclide) % total
! Add contributions to material macroscopic scattering cross section
material_xs % elastic = material_xs % elastic + &
atom_density * micro_xs(i_nuclide) % elastic
! Add contributions to material macroscopic absorption cross section
material_xs % absorption = material_xs % absorption + &
atom_density * micro_xs(i_nuclide) % absorption
@ -162,6 +157,7 @@ contains
real(8) :: sig_t, sig_a, sig_f ! Intermediate multipole variables
! Initialize cached cross sections to zero
micro_xs(i_nuclide) % elastic = CACHE_INVALID
micro_xs(i_nuclide) % thermal = ZERO
micro_xs(i_nuclide) % thermal_elastic = ZERO
@ -182,13 +178,11 @@ contains
micro_xs(i_nuclide) % total = sig_t
micro_xs(i_nuclide) % absorption = sig_a
micro_xs(i_nuclide) % elastic = sig_t - sig_a
micro_xs(i_nuclide) % fission = sig_f
if (nuc % fissionable) then
micro_xs(i_nuclide) % fission = sig_f
micro_xs(i_nuclide) % nu_fission = sig_f * nuc % nu(E, EMISSION_TOTAL)
else
micro_xs(i_nuclide) % fission = ZERO
micro_xs(i_nuclide) % nu_fission = ZERO
end if
@ -233,7 +227,7 @@ contains
if (f > prn()) i_temp = i_temp + 1
end select
associate (grid => nuc % grid(i_temp), xs => nuc % sum_xs(i_temp))
associate (grid => nuc % grid(i_temp), xs => nuc % xs(i_temp))
! Determine the energy grid index using a logarithmic mapping to
! reduce the energy range over which a binary search needs to be
! performed
@ -266,25 +260,21 @@ contains
micro_xs(i_nuclide) % interp_factor = f
! Calculate microscopic nuclide total cross section
micro_xs(i_nuclide) % total = (ONE - f) * xs % total(i_grid) &
+ f * xs % total(i_grid + 1)
! Calculate microscopic nuclide elastic cross section
micro_xs(i_nuclide) % elastic = (ONE - f) * xs % elastic(i_grid) &
+ f * xs % elastic(i_grid + 1)
micro_xs(i_nuclide) % total = (ONE - f) * xs % value(XS_TOTAL,i_grid) &
+ f * xs % value(XS_TOTAL,i_grid + 1)
! Calculate microscopic nuclide absorption cross section
micro_xs(i_nuclide) % absorption = (ONE - f) * xs % absorption( &
i_grid) + f * xs % absorption(i_grid + 1)
micro_xs(i_nuclide) % absorption = (ONE - f) * xs % value(XS_ABSORPTION, &
i_grid) + f * xs % value(XS_ABSORPTION,i_grid + 1)
if (nuc % fissionable) then
! Calculate microscopic nuclide total cross section
micro_xs(i_nuclide) % fission = (ONE - f) * xs % fission(i_grid) &
+ f * xs % fission(i_grid + 1)
micro_xs(i_nuclide) % fission = (ONE - f) * xs % value(XS_FISSION,i_grid) &
+ f * xs % value(XS_FISSION,i_grid + 1)
! Calculate microscopic nuclide nu-fission cross section
micro_xs(i_nuclide) % nu_fission = (ONE - f) * xs % nu_fission( &
i_grid) + f * xs % nu_fission(i_grid + 1)
micro_xs(i_nuclide) % nu_fission = (ONE - f) * xs % value(XS_NU_FISSION, &
i_grid) + f * xs % value(XS_NU_FISSION,i_grid + 1)
else
micro_xs(i_nuclide) % fission = ZERO
micro_xs(i_nuclide) % nu_fission = ZERO
@ -451,6 +441,9 @@ contains
micro_xs(i_nuclide) % thermal = sab_frac * (elastic + inelastic)
micro_xs(i_nuclide) % thermal_elastic = sab_frac * elastic
! Calculate free atom elastic cross section
call calculate_elastic_xs(i_nuclide)
! Correct total and elastic cross sections
micro_xs(i_nuclide) % total = micro_xs(i_nuclide) % total &
+ micro_xs(i_nuclide) % thermal &
@ -579,6 +572,7 @@ contains
! Multiply by smooth cross-section if needed
if (urr % multiply_smooth) then
call calculate_elastic_xs(i_nuclide)
elastic = elastic * micro_xs(i_nuclide) % elastic
capture = capture * (micro_xs(i_nuclide) % absorption - &
micro_xs(i_nuclide) % fission)
@ -607,6 +601,36 @@ contains
end subroutine calculate_urr_xs
!===============================================================================
! CALCULATE_ELASTIC_XS precalculates the free atom elastic scattering cross
! section. Normally it is not needed until a collision actually occurs in a
! material. However, in the thermal and unresolved resonance regions, we have to
! calculate it early to adjust the total cross section correctly.
!===============================================================================
subroutine calculate_elastic_xs(i_nuclide)
integer, intent(in) :: i_nuclide
integer :: i_temp
integer :: i_grid
real(8) :: f
! Get temperature index, grid index, and interpolation factor
i_temp = micro_xs(i_nuclide) % index_temp
i_grid = micro_xs(i_nuclide) % index_grid
f = micro_xs(i_nuclide) % interp_factor
if (i_temp > 0) then
associate (xs => nuclides(i_nuclide) % reactions(1) % xs(i_temp) % value)
micro_xs(i_nuclide) % elastic = (ONE - f)*xs(i_grid) + f*xs(i_grid + 1)
end associate
else
! For multipole, elastic is total - absorption
micro_xs(i_nuclide) % elastic = micro_xs(i_nuclide) % total - &
micro_xs(i_nuclide) % absorption
end if
end subroutine calculate_elastic_xs
!===============================================================================
! MULTIPOLE_EVAL evaluates the windowed multipole equations for cross
! sections in the resolved resonance regions

71
src/nuclide_header.F90
View file

@ -36,13 +36,18 @@ module nuclide_header
real(8), allocatable :: energy(:) ! energy values corresponding to xs
end type EnergyGrid
! Positions for first dimension of Nuclide % xs
integer, parameter :: &
XS_TOTAL = 1, &
XS_ABSORPTION = 2, &
XS_FISSION = 3, &
XS_NU_FISSION = 4
! The array within SumXS is of shape (4, n_energy) where the first dimension
! corresponds to the following values: 1) total, 2) absorption (MT > 100), 3)
! fission, 4) neutron production
type SumXS
real(8), allocatable :: total(:) ! total cross section
real(8), allocatable :: elastic(:) ! elastic scattering
real(8), allocatable :: fission(:) ! fission
real(8), allocatable :: nu_fission(:) ! neutron production
real(8), allocatable :: absorption(:) ! absorption (MT > 100)
real(8), allocatable :: heating(:) ! heating
real(8), allocatable :: value(:,:)
end type SumXS
type :: Nuclide
@ -61,7 +66,7 @@ module nuclide_header
type(EnergyGrid), allocatable :: grid(:)
! Microscopic cross sections
type(SumXS), allocatable :: sum_xs(:)
type(SumXS), allocatable :: xs(:)
! Resonance scattering info
logical :: resonant = .false. ! resonant scatterer?
@ -110,14 +115,19 @@ module nuclide_header
! nuclide at the current energy
!===============================================================================
! Arbitrary value to indicate invalid cache state for elastic scattering
! (NuclideMicroXS % elastic)
real(8), parameter :: CACHE_INVALID = dble(Z"FFE0000000000000")
type NuclideMicroXS
! Microscopic cross sections in barns
real(8) :: total
real(8) :: elastic ! If sab_frac is not 1 or 0, then this value is
! averaged over bound and non-bound nuclei
real(8) :: absorption ! absorption (disappearance)
real(8) :: fission ! fission
real(8) :: nu_fission ! neutron production from fission
real(8) :: elastic ! If sab_frac is not 1 or 0, then this value is
! averaged over bound and non-bound nuclei
real(8) :: thermal ! Bound thermal elastic & inelastic scattering
real(8) :: thermal_elastic ! Bound thermal elastic scattering
@ -148,7 +158,6 @@ module nuclide_header
type MaterialMacroXS
real(8) :: total ! macroscopic total xs
real(8) :: elastic ! macroscopic elastic scattering xs
real(8) :: absorption ! macroscopic absorption xs
real(8) :: fission ! macroscopic fission xs
real(8) :: nu_fission ! macroscopic production xs
@ -566,21 +575,13 @@ contains
type(VectorInt) :: MTs
n_temperature = size(this % kTs)
allocate(this % sum_xs(n_temperature))
allocate(this % xs(n_temperature))
this % reaction_index(:) = 0
do i = 1, n_temperature
! Allocate and initialize derived cross sections
n_grid = size(this % grid(i) % energy)
allocate(this % sum_xs(i) % total(n_grid))
allocate(this % sum_xs(i) % elastic(n_grid))
allocate(this % sum_xs(i) % fission(n_grid))
allocate(this % sum_xs(i) % nu_fission(n_grid))
allocate(this % sum_xs(i) % absorption(n_grid))
this % sum_xs(i) % total(:) = ZERO
this % sum_xs(i) % elastic(:) = ZERO
this % sum_xs(i) % fission(:) = ZERO
this % sum_xs(i) % nu_fission(:) = ZERO
this % sum_xs(i) % absorption(:) = ZERO
allocate(this % xs(i) % value(4,n_grid))
this % xs(i) % value(:,:) = ZERO
end do
i_fission = 0
@ -607,17 +608,14 @@ contains
j = rx % xs(t) % threshold
n = size(rx % xs(t) % value)
! Copy elastic
if (rx % MT == ELASTIC) this % sum_xs(t) % elastic(:) = rx % xs(t) % value
! Add contribution to total cross section
this % sum_xs(t) % total(j:j+n-1) = this % sum_xs(t) % total(j:j+n-1) + &
rx % xs(t) % value
this % xs(t) % value(XS_TOTAL,j:j+n-1) = this % xs(t) % &
value(XS_TOTAL,j:j+n-1) + rx % xs(t) % value
! Add contribution to absorption cross section
if (is_disappearance(rx % MT)) then
this % sum_xs(t) % absorption(j:j+n-1) = this % sum_xs(t) % &
absorption(j:j+n-1) + rx % xs(t) % value
this % xs(t) % value(XS_ABSORPTION,j:j+n-1) = this % xs(t) % &
value(XS_ABSORPTION,j:j+n-1) + rx % xs(t) % value
end if
! Information about fission reactions
@ -633,12 +631,12 @@ contains
! Add contribution to fission cross section
if (is_fission(rx % MT)) then
this % fissionable = .true.
this % sum_xs(t) % fission(j:j+n-1) = this % sum_xs(t) % &
fission(j:j+n-1) + rx % xs(t) % value
this % xs(t) % value(XS_FISSION,j:j+n-1) = this % xs(t) % &
value(XS_FISSION,j:j+n-1) + rx % xs(t) % value
! Also need to add fission cross sections to absorption
this % sum_xs(t) % absorption(j:j+n-1) = this % sum_xs(t) % &
absorption(j:j+n-1) + rx % xs(t) % value
this % xs(t) % value(XS_ABSORPTION,j:j+n-1) = this % xs(t) % &
value(XS_ABSORPTION,j:j+n-1) + rx % xs(t) % value
! If total fission reaction is present, there's no need to store the
! reaction cross-section since it was copied to this % fission
@ -689,12 +687,11 @@ contains
! Calculate nu-fission cross section
do t = 1, n_temperature
if (this % fissionable) then
do i = 1, size(this % sum_xs(t) % fission)
this % sum_xs(t) % nu_fission(i) = this % nu(this % grid(t) % energy(i), &
EMISSION_TOTAL) * this % sum_xs(t) % fission(i)
do i = 1, n_grid
this % xs(t) % value(XS_NU_FISSION,i) = &
this % nu(this % grid(t) % energy(i), EMISSION_TOTAL) * &
this % xs(t) % value(XS_FISSION,i)
end do
else
this % sum_xs(t) % nu_fission(:) = ZERO
end if
end do
end subroutine nuclide_create_derived

7
src/physics.F90
View file

@ -2,7 +2,7 @@ module physics
use algorithm, only: binary_search
use constants
use cross_section, only: elastic_xs_0K
use cross_section, only: elastic_xs_0K, calculate_elastic_xs
use endf, only: reaction_name
use error, only: fatal_error, warning, write_message
use material_header, only: Material, materials
@ -338,6 +338,11 @@ contains
micro_xs(i_nuclide) % absorption)
sampled = .false.
! Calculate elastic cross section if it wasn't precalculated
if (micro_xs(i_nuclide) % elastic == CACHE_INVALID) then
call calculate_elastic_xs(i_nuclide)
end if
prob = micro_xs(i_nuclide) % elastic - micro_xs(i_nuclide) % thermal
if (prob > cutoff) then
! =======================================================================

21
src/tallies/tally.F90
View file

@ -4,7 +4,7 @@ module tally
use algorithm, only: binary_search
use constants
use cross_section, only: multipole_deriv_eval
use cross_section, only: multipole_deriv_eval, calculate_elastic_xs
use dict_header, only: EMPTY
use error, only: fatal_error
use geometry_header
@ -1017,9 +1017,26 @@ contains
else
if (i_nuclide > 0) then
if (micro_xs(i_nuclide) % elastic == CACHE_INVALID) then
call calculate_elastic_xs(i_nuclide)
end if
score = micro_xs(i_nuclide) % elastic * atom_density * flux
else
score = material_xs % elastic * flux
score = ZERO
if (p % material /= MATERIAL_VOID) then
do l = 1, materials(p % material) % n_nuclides
! Get atom density
atom_density_ = materials(p % material) % atom_density(l)
! Get index in nuclides array
i_nuc = materials(p % material) % nuclide(l)
if (micro_xs(i_nuc) % elastic == CACHE_INVALID) then
call calculate_elastic_xs(i_nuc)
end if
score = score + micro_xs(i_nuc) % elastic * atom_density_ * flux
end do
end if
end if
end if