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

completed conversion of Urr to C++

Browse source
This commit is contained in:
Adam G Nelson 2018-12-06 19:28:26 -05:00
parent 5eab6dbed0
commit d4c8c3d34e
7 changed files with 25 additions and 140 deletions
Download patch file Download diff file Expand all files Collapse all files

1
CMakeLists.txt
View file

@ -361,7 +361,6 @@ add_library(libopenmc SHARED
src/timer_header.F90
src/tracking.F90
src/track_output.F90
src/urr_header.F90
src/vector_header.F90
src/volume_calc.F90
src/volume_header.F90

3
include/openmc/nuclide.h
View file

@ -175,10 +175,9 @@ extern "C" void set_micro_xs();
extern "C" bool nuclide_wmp_present(int i_nuclide);
extern "C" double nuclide_wmp_emin(int i_nuclide);
extern "C" double nuclide_wmp_emax(int i_nuclide);
extern "C" void nuclide_calculate_urr_xs(const int i_nuclide,
extern "C" void nuclide_calculate_urr_xs(const bool use_mp, const int i_nuclide,
const int i_temp, const double E);
} // namespace openmc
#endif // OPENMC_NUCLIDE_H

1
include/openmc/urr.h
View file

@ -20,6 +20,7 @@ public:
int inelastic_flag_; // inelastic competition flag
int absorption_flag_; // other absorption flag
bool multiply_smooth_; // multiply by smooth cross section?
int n_energy_; // number of energy points
xt::xtensor<double, 1> energy_; // incident energies
xt::xtensor<double, 3> prob_; // Actual probability tables

12
src/nuclide.cpp
View file

@ -576,9 +576,17 @@ void set_micro_xs()
}
extern "C" void
nuclide_calculate_urr_xs(const int i_nuclide, const int i_temp, const double E)
nuclide_calculate_urr_xs(const bool use_mp, const int i_nuclide,
const int i_temp, const double E)
{
data::nuclides[i_nuclide - 1]->calculate_urr_xs(i_temp - 1, E);
Nuclide* nuc = data::nuclides[i_nuclide - 1].get();
if (settings::urr_ptables_on && (nuc->urr_present_ && !use_mp)) {
if ((E > nuc->urr_data_[i_temp - 1].energy_(0)) &&
(E < nuc->urr_data_[i_temp - 1].energy_(
nuc->urr_data_[i_temp - 1].n_energy_ - 1))) {
nuc->calculate_urr_xs(i_temp - 1, E);
}
}
}
} // namespace openmc

73
src/nuclide_header.F90
View file

@ -23,7 +23,6 @@ module nuclide_header
use stl_vector, only: VectorInt, VectorReal
use string
use simulation_header, only: need_depletion_rx
use urr_header, only: UrrData
implicit none
@ -76,11 +75,6 @@ module nuclide_header
integer, allocatable :: index_fission(:) ! indices in reactions
class(Function1D), allocatable :: total_nu
! Unresolved resonance data
logical :: urr_present = .false.
integer :: urr_inelastic
type(UrrData), allocatable :: urr_data(:)
! Multipole data
logical :: mp_present = .false.
type(MultipoleArray), pointer :: multipole => null()
@ -199,12 +193,13 @@ module nuclide_header
type(C_PTR) :: path
end function
subroutine nuclide_calculate_urr_xs_c(i_nuclide, i_temp, E) &
subroutine nuclide_calculate_urr_xs_c(use_mp, i_nuclide, i_temp, E) &
bind(C, name='nuclide_calculate_urr_xs')
import C_INT, C_DOUBLE
integer(C_INT), value, intent(in) :: i_nuclide
integer(C_INT), value, intent(in) :: i_temp
real(C_DOUBLE), value, intent(in) :: E
import C_BOOL, C_INT, C_DOUBLE
logical(C_BOOL), value, intent(in) :: use_mp
integer(C_INT), value, intent(in) :: i_nuclide
integer(C_INT), value, intent(in) :: i_temp
real(C_DOUBLE), value, intent(in) :: E
end subroutine nuclide_calculate_urr_xs_c
end interface
@ -261,7 +256,7 @@ contains
integer :: i
integer :: i_closest
integer :: n_temperature
integer(HID_T) :: urr_group, nu_group
integer(HID_T) :: nu_group
integer(HID_T) :: energy_group, energy_dset
integer(HID_T) :: kT_group
integer(HID_T) :: rxs_group
@ -446,47 +441,6 @@ contains
end do
call close_group(rxs_group)
! Read unresolved resonance probability tables if present
if (object_exists(group_id, 'urr')) then
this % urr_present = .true.
allocate(this % urr_data(n_temperature))
do i = 1, n_temperature
! Get temperature as a string
temp_str = trim(to_str(temps_to_read % data(i))) // "K"
! Read probability tables for i-th temperature
urr_group = open_group(group_id, 'urr/' // trim(temp_str))
call this % urr_data(i) % from_hdf5(urr_group)
call close_group(urr_group)
! Check for negative values
if (any(this % urr_data(i) % prob < ZERO) .and. master) then
call warning("Negative value(s) found on probability table &
&for nuclide " // this % name // " at " // trim(temp_str))
end if
end do
! if the inelastic competition flag indicates that the inelastic cross
! section should be determined from a normal reaction cross section, we
! need to get the index of the reaction
if (n_temperature > 0) then
if (this % urr_data(1) % inelastic_flag > 0) then
do i = 1, size(this % reactions)
if (this % reactions(i) % MT == this % urr_data(1) % inelastic_flag) then
this % urr_inelastic = i
end if
end do
! Abort if no corresponding inelastic reaction was found
if (this % urr_inelastic == NONE) then
call fatal_error("Could not find inelastic reaction specified on &
&unresolved resonance probability table.")
end if
end if
end if
end if
! Check for nu-total
if (object_exists(group_id, 'total_nu')) then
nu_group = open_group(group_id, 'total_nu')
@ -784,7 +738,7 @@ contains
real(8), intent(in) :: sab_frac ! fraction of atoms affected by S(a,b)
type(NuclideMicroXS), intent(inout) :: micro_xs ! Cross section cache
logical :: use_mp ! true if XS can be calculated with windowed multipole
logical(C_BOOL) :: use_mp ! true if XS can be calculated with windowed multipole
integer :: i_temp ! index for temperature
integer :: i_grid ! index on nuclide energy grid
integer :: i_low ! lower logarithmic mapping index
@ -839,7 +793,7 @@ contains
! 1. physics.F90 - scatter - For inelastic scatter.
! 2. physics.F90 - sample_fission - For partial fissions.
! 3. tally.F90 - score_general - For tallying on MTxxx reactions.
! 4. nuclide_header.F90 - calculate_urr_xs - For unresolved purposes.
! 4. nuclide.h - calculate_urr_xs - For unresolved purposes.
! It is worth noting that none of these occur in the resolved
! resonance range, so the value here does not matter. index_temp is
! set to -1 to force a segfault in case a developer messes up and tries
@ -982,15 +936,10 @@ contains
call calculate_sab_xs(this, i_sab, E, sqrtkT, sab_frac, micro_xs)
end if
! If the particle is in the unresolved resonance range and there are
! probability tables, we need to determine cross sections from the table
if (urr_ptables_on .and. this % urr_present .and. .not. use_mp) then
if (E > this % urr_data(i_temp) % energy(1) .and. E < this % &
urr_data(i_temp) % energy(this % urr_data(i_temp) % n_energy)) then
call nuclide_calculate_urr_xs_c(this % i_nuclide, i_temp, E)
end if
end if
call nuclide_calculate_urr_xs_c(use_mp, this % i_nuclide, i_temp, E)
micro_xs % last_E = E
micro_xs % last_sqrtkT = sqrtkT

3
src/urr.cpp
View file

@ -27,6 +27,7 @@ UrrData::from_hdf5(hid_t group_id)
interp_ = Interpolation::log_log;
}
// read the metadata
read_attribute(group_id, "inelastic", inelastic_flag_);
read_attribute(group_id, "absorption", absorption_flag_);
int temp_multiply_smooth;
@ -38,7 +39,7 @@ UrrData::from_hdf5(hid_t group_id)
hsize_t dims[1];
get_shape(dset, dims);
close_dataset(dset);
// n_energy_ = static_cast<int>(dims[0]);
n_energy_ = static_cast<int>(dims[0]);
energy_ = xt::xtensor<double, 1>({dims[0]}, 0.);
read_dataset_as_shape(group_id, "energy", energy_);

72
src/urr_header.F90
View file

@ -1,72 +0,0 @@
module urr_header
use hdf5_interface, only: read_attribute, open_dataset, read_dataset, &
close_dataset, get_shape, HID_T, HSIZE_T
implicit none
!===============================================================================
! URRDATA contains probability tables for the unresolved resonance range.
!===============================================================================
type UrrData
integer :: n_energy ! # of incident neutron energies
integer :: n_prob ! # of probabilities
integer :: interp ! inteprolation (2=lin-lin, 5=log-log)
integer :: inelastic_flag ! inelastic competition flag
integer :: absorption_flag ! other absorption flag
logical :: multiply_smooth ! multiply by smooth cross section?
real(8), allocatable :: energy(:) ! incident energies
real(8), allocatable :: prob(:,:,:) ! actual probabibility tables
contains
procedure :: from_hdf5 => urr_from_hdf5
end type UrrData
contains
subroutine urr_from_hdf5(this, group_id)
class(UrrData), intent(inout) :: this
integer(HID_T), intent(in) :: group_id
integer :: i, j, k
integer(HID_T) :: energy
integer(HID_T) :: table
integer(HSIZE_T) :: dims(1)
integer(HSIZE_T) :: dims3(3)
real(8), allocatable :: temp(:,:,:)
! Read interpolation and other flags
call read_attribute(this % interp, group_id, 'interpolation')
call read_attribute(this % inelastic_flag, group_id, 'inelastic')
call read_attribute(this % absorption_flag, group_id, 'absorption')
call read_attribute(i, group_id, 'multiply_smooth')
this % multiply_smooth = (i == 1)
! Read energies at which tables exist
energy = open_dataset(group_id, 'energy')
call get_shape(energy, dims)
this % n_energy = int(dims(1), 4)
allocate(this % energy(this % n_energy))
call read_dataset(this % energy, energy)
call close_dataset(energy)
! Read URR tables
table = open_dataset(group_id, 'table')
call get_shape(table, dims3)
this % n_prob = int(dims3(1), 4)
allocate(temp(this % n_prob, 6, this % n_energy))
call read_dataset(temp, table)
call close_dataset(table)
! Swap first and last indices
allocate(this % prob(this % n_energy, 6, this % n_prob))
do i = 1, this % n_energy
do j = 1, 6
do k = 1, this % n_prob
this % prob(i, j, k) = temp(k, j, i)
end do
end do
end do
end subroutine urr_from_hdf5
end module urr_header