crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-28 06:05:58 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Merge pull request #876 from paulromano/assign-sab-fix

Browse source 
Fix assignment of S(a,b) tables
This commit is contained in:
Will Boyd 2017-05-10 13:51:18 -04:00 • committed by GitHub
commit fa32b036fc
2 changed files with 45 additions and 34 deletions
Download patch file Download diff file Expand all files Collapse all files

77
src/input_xml.F90
View file

@ -5038,6 +5038,7 @@ contains
integer :: m ! position for sorting
integer :: temp_nuclide ! temporary value for sorting
integer :: temp_table ! temporary value for sorting
logical :: found
type(VectorInt) :: i_sab_tables
type(VectorInt) :: i_sab_nuclides
@ -5051,17 +5052,19 @@ contains
! In order to know which nuclide the S(a,b) table applies to, we need
! to search through the list of nuclides for one which has a matching
! name
found = .false.
associate (sab => sab_tables(mat % i_sab_tables(k)))
FIND_NUCLIDE: do j = 1, size(mat % nuclide)
if (any(sab % nuclides == nuclides(mat % nuclide(j)) % name)) then
call i_sab_tables % push_back(k)
call i_sab_tables % push_back(mat % i_sab_tables(k))
call i_sab_nuclides % push_back(j)
found = .true.
end if
end do FIND_NUCLIDE
end associate
! Check to make sure S(a,b) table matched a nuclide
if (find(i_sab_tables, k) == -1) then
if (.not. found) then
call fatal_error("S(a,b) table " // trim(mat % &
sab_names(k)) // " did not match any nuclide on material " &
// trim(to_str(mat % id)))
@ -5298,43 +5301,51 @@ contains
subroutine assign_0K_elastic_scattering(nuc)
type(Nuclide), intent(inout) :: nuc
integer :: i, j
integer :: i
real(8) :: xs_cdf_sum
do i = 1, size(res_scat_nuclides)
if (nuc % name == res_scat_nuclides(i)) then
! Make sure nuclide has 0K data
if (.not. allocated(nuc % energy_0K)) then
call fatal_error("Cannot treat " // trim(nuc % name) // " as a &
&resonant scatterer because 0 K elastic scattering data is not &
&present.")
nuc % resonant = .false.
if (allocated(res_scat_nuclides)) then
! If resonant nuclides were specified, check the list explicitly
do i = 1, size(res_scat_nuclides)
if (nuc % name == res_scat_nuclides(i)) then
nuc % resonant = .true.
! Make sure nuclide has 0K data
if (.not. allocated(nuc % energy_0K)) then
call fatal_error("Cannot treat " // trim(nuc % name) // " as a &
&resonant scatterer because 0 K elastic scattering data is &
&not present.")
end if
exit
end if
end do
else
! Otherwise, assume that any that have 0 K elastic scattering data are
! resonant
nuc % resonant = allocated(nuc % energy_0K)
end if
! Set nuclide to be resonant
nuc % resonant = .true.
if (nuc % resonant) then
! Build CDF for 0K elastic scattering
xs_cdf_sum = ZERO
allocate(nuc % xs_cdf(0:size(nuc % energy_0K)))
nuc % xs_cdf(0) = ZERO
! Build CDF for 0K elastic scattering
xs_cdf_sum = ZERO
allocate(nuc % xs_cdf(0:size(nuc % energy_0K)))
nuc % xs_cdf(0) = ZERO
associate (E => nuc % energy_0K, xs => nuc % elastic_0K)
do j = 1, size(E) - 1
! Negative cross sections result in a CDF that is not monotonically
! increasing. Set all negative xs values to zero.
if (xs(j) < ZERO) xs(j) = ZERO
! build xs cdf
xs_cdf_sum = xs_cdf_sum + (sqrt(E(j))*xs(j) + sqrt(E(j+1))*xs(j+1))&
/ TWO * (E(j+1) - E(j))
nuc % xs_cdf(j) = xs_cdf_sum
end do
end associate
exit
end if
end do
associate (E => nuc % energy_0K, xs => nuc % elastic_0K)
do i = 1, size(E) - 1
! Negative cross sections result in a CDF that is not monotonically
! increasing. Set all negative xs values to zero.
if (xs(i) < ZERO) xs(i) = ZERO
! build xs cdf
xs_cdf_sum = xs_cdf_sum + (sqrt(E(i))*xs(i) + sqrt(E(i+1))*xs(i+1))&
/ TWO * (E(i+1) - E(i))
nuc % xs_cdf(i) = xs_cdf_sum
end do
end associate
end if
end subroutine assign_0K_elastic_scattering
!===============================================================================

2
tests/test_salphabeta/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
8.465448E-01 5.743842E-03
8.544160E-01 1.274133E-02