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

Use interoperable TallyDerivative structs

Browse source
This commit is contained in:
Sterling Harper 2019-01-26 13:53:01 -05:00
parent 06caf67ccf
commit 9f6a7466e2
11 changed files with 139 additions and 158 deletions
Download patch file Download diff file Expand all files Collapse all files

2
include/openmc/tallies/derivative.h
View file

@ -12,7 +12,7 @@
namespace openmc {
struct TallyDerivative {
extern "C" struct TallyDerivative {
int id; //!< User-defined identifier
int variable; //!< Independent variable (like temperature)
int diff_material; //!< Material this derivative is applied to

4
include/openmc/tallies/tally.h
View file

@ -39,8 +39,10 @@ public:
int surface_filter_ {C_NONE};
int mesh_filter_ {C_NONE};
int deriv_ {F90_NONE}; //!< Index of a TallyDerivative object for diff tallies.
private:
std::vector<int32_t> filters_; //< Filter indices in global filters array
std::vector<int32_t> filters_; //!< Filter indices in global filters array
//! Index strides assigned to each filter to support 1D indexing.
std::vector<int32_t> strides_;

2
src/api.F90
View file

@ -66,7 +66,6 @@ contains
use settings
use simulation_header
use surface_header
use tally_derivative_header
use tally_filter_header
use tally_header
use trigger_header
@ -102,7 +101,6 @@ contains
call free_memory_mesh()
call free_memory_tally()
call free_memory_tally_filter()
call free_memory_tally_derivative()
call free_memory_bank()
call free_memory_cmfd()
#ifdef DAGMC

49
src/input_xml.F90
View file

@ -913,6 +913,7 @@ contains
integer :: l ! loop over bins
integer :: filter_id ! user-specified identifier for filter
integer :: tally_id ! user-specified identifier for filter
integer :: deriv_id
integer :: i_filt ! index in filters array
integer :: i_elem ! index of entry in dictionary
integer :: n ! size of arrays in mesh specification
@ -948,6 +949,7 @@ contains
type(XMLNode), allocatable :: node_trigger_list(:)
type(XMLNode), allocatable :: node_deriv_list(:)
type(DictEntryCI) :: elem
type(TallyDerivative), pointer :: deriv
interface
subroutine read_tally_derivatives(node_ptr) bind(C)
@ -960,12 +962,6 @@ contains
filename = trim(path_input) // "tallies.xml"
inquire(FILE=filename, EXIST=file_exists)
if (.not. file_exists) then
! We need to allocate tally_derivs to avoid segfaults. Also needs to be
! done in parallel because tally derivs are threadprivate.
!$omp parallel
allocate(tally_derivs(0))
!$omp end parallel
! Since a tallies.xml file is optional, no error is issued here
return
end if
@ -1008,28 +1004,10 @@ contains
call read_tally_derivatives(root % ptr)
! Get pointer list to XML <derivative> nodes and allocate global array.
! The array is threadprivate so it must be allocated in parallel.
call get_node_list(root, "derivative", node_deriv_list)
!$omp parallel
allocate(tally_derivs(size(node_deriv_list)))
!$omp end parallel
! Make sure this is not an MG run.
if (.not. run_CE .and. size(node_deriv_list) > 0) then
call fatal_error("Differential tallies not supported in multi-group mode")
end if
! Read derivative attributes.
do i = 1, size(node_deriv_list)
!$omp parallel
!$omp critical (ReadTallyDeriv)
call tally_derivs(i) % from_xml(node_deriv_list(i))
!$omp end critical (ReadTallyDeriv)
!$omp end parallel
! Update tally derivative dictionary
call tally_deriv_dict % set(tally_derivs(i) % id, i)
deriv => tally_deriv_c(i - 1)
call tally_deriv_dict % set(deriv % id, i)
end do
! ==========================================================================
@ -1618,13 +1596,13 @@ contains
! Check for a tally derivative.
if (check_for_node(node_tal, "derivative")) then
! Temporarily store the derivative id.
call get_node_value(node_tal, "derivative", t % deriv)
call get_node_value(node_tal, "derivative", deriv_id)
! Find the derivative with the given id, and store it's index.
do j = 1, size(tally_derivs)
if (tally_derivs(j) % id == t % deriv) then
t % deriv = j
do j = 1, n_tally_derivs()
deriv => tally_deriv_c(j - 1)
if (deriv % id == deriv_id) then
call t % set_deriv(j)
! Only analog or collision estimators are supported for differential
! tallies.
if (t % estimator == ESTIMATOR_TRACKLENGTH) then
@ -1633,15 +1611,16 @@ contains
! We found the derivative we were looking for; exit the do loop.
exit
end if
if (j == size(tally_derivs)) then
if (j == n_tally_derivs()) then
call fatal_error("Could not find derivative " &
// trim(to_str(t % deriv)) // " specified on tally " &
// trim(to_str(deriv_id)) // " specified on tally " &
// trim(to_str(t % id)))
end if
end do
if (tally_derivs(t % deriv) % variable == DIFF_NUCLIDE_DENSITY &
.or. tally_derivs(t % deriv) % variable == DIFF_TEMPERATURE) then
deriv => tally_deriv_c(t % deriv() - 1)
if (deriv % variable == DIFF_NUCLIDE_DENSITY &
.or. deriv % variable == DIFF_TEMPERATURE) then
if (any(t % nuclide_bins == -1)) then
if (has_energyout) then
call fatal_error("Error on tally " // trim(to_str(t % id)) &

8
src/output.F90
View file

@ -471,6 +471,7 @@ contains
! to be applied at write-time
integer, allocatable :: filter_bins(:)
character(MAX_WORD_LEN) :: temp_name
type(TallyDerivative), pointer :: deriv
! Skip if there are no tallies
if (n_tallies == 0) return
@ -530,8 +531,9 @@ contains
endif
! Write derivative information.
if (t % deriv /= NONE) then
associate(deriv => tally_derivs(t % deriv))
if (t % deriv() /= NONE) then
!associate(deriv => tally_derivs(t % deriv()))
deriv => tally_deriv_c(t % deriv() - 1)
select case (deriv % variable)
case (DIFF_DENSITY)
write(unit=unit_tally, fmt="(' Density derivative Material ',A)") &
@ -548,7 +550,7 @@ contains
call fatal_error("Differential tally dependent variable for tally "&
// trim(to_str(t % id)) // " not defined in output.F90.")
end select
end associate
!end associate
end if
! WARNING: Admittedly, the logic for moving for printing results is

17
src/state_point.F90
View file

@ -28,7 +28,7 @@ module state_point
use string, only: to_str, count_digits, zero_padded, to_f_string
use tally_header
use tally_filter_header
use tally_derivative_header, only: tally_derivs
use tally_derivative_header
use timer_header
implicit none
@ -69,6 +69,7 @@ contains
character(len=:, kind=C_CHAR), allocatable :: filename_
character(MAX_WORD_LEN, kind=C_CHAR) :: temp_name
logical :: parallel
type(TallyDerivative), pointer :: deriv
interface
subroutine meshes_to_hdf5(group) bind(C)
@ -174,10 +175,11 @@ contains
call meshes_to_hdf5(tallies_group)
! Write information for derivatives.
if (size(tally_derivs) > 0) then
if (n_tally_derivs() > 0) then
derivs_group = create_group(tallies_group, "derivatives")
do i = 1, size(tally_derivs)
associate(deriv => tally_derivs(i))
do i = 1, n_tally_derivs()
!associate(deriv => tally_derivs(i))
deriv => tally_deriv_c(i - 1)
deriv_group = create_group(derivs_group, "derivative " &
// trim(to_str(deriv % id)))
select case (deriv % variable)
@ -200,7 +202,7 @@ contains
&state_point.F90.")
end select
call close_group(deriv_group)
end associate
!end associate
end do
call close_group(derivs_group)
end if
@ -303,9 +305,10 @@ contains
deallocate(str_array)
! Write derivative information.
if (tally % deriv /= NONE) then
if (tally % deriv() /= NONE) then
deriv => tally_deriv_c(tally % deriv() - 1)
call write_dataset(tally_group, "derivative", &
tally_derivs(tally % deriv) % id)
deriv % id)
end if
! Write scores.

35
src/tallies/derivative.cpp
View file

@ -33,12 +33,12 @@ TallyDerivative::TallyDerivative(pugi::xml_node node)
if (id <= 0)
fatal_error("<derivative> IDs must be an integer greater than zero");
std::string variable = get_node_value(node, "variable");
std::string variable_str = get_node_value(node, "variable");
if (variable == "density") {
if (variable_str == "density") {
variable = DIFF_DENSITY;
} else if (variable == "nuclide_density") {
} else if (variable_str == "nuclide_density") {
variable = DIFF_NUCLIDE_DENSITY;
std::string nuclide_name = get_node_value(node, "nuclide");
@ -46,7 +46,8 @@ TallyDerivative::TallyDerivative(pugi::xml_node node)
for (auto i = 0; i < data::nuclides.size(); ++i) {
if (data::nuclides[i]->name_ == nuclide_name) {
found = true;
diff_nuclide = i;
//TODO: off-by-one
diff_nuclide = i + 1;
}
}
if (!found) {
@ -56,12 +57,13 @@ TallyDerivative::TallyDerivative(pugi::xml_node node)
fatal_error(out);
}
} else if (variable == "temperature") {
} else if (variable_str == "temperature") {
variable = DIFF_TEMPERATURE;
} else {
std::stringstream out;
out << "Unrecognized variable \"" << variable << "\" on derivative " << id;
out << "Unrecognized variable \"" << variable_str
<< "\" on derivative " << id;
fatal_error(out);
}
@ -100,4 +102,25 @@ read_tally_derivatives(pugi::xml_node* node)
fatal_error("Differential tallies not supported in multi-group mode");
}
/*
//! Set the flux derivatives on differential tallies to zero.
extern "C" void
zero_flux_derivs_c()
{
for (auto& deriv : model::tally_derivs) deriv.flux_deriv = 0.;
}
*/
//==============================================================================
// Fortran interop
//==============================================================================
extern "C" int n_tally_derivs() {return model::tally_derivs.size();}
extern "C" TallyDerivative*
tally_deriv_c(int i)
{
return &(model::tally_derivs[i]);
}
}// namespace openmc

37
src/tallies/tally.F90
View file

@ -17,7 +17,7 @@ module tally
use settings
use simulation_header
use string, only: to_str
use tally_derivative_header, only: tally_derivs
use tally_derivative_header
use tally_filter
use tally_header
@ -1184,7 +1184,7 @@ contains
!#########################################################################
! Add derivative information on score for differential tallies.
if (t % deriv /= NONE) then
if (t % deriv() /= NONE) then
call apply_derivative_to_score(p, t, i_nuclide, atom_density, &
score_bin, score)
end if
@ -2375,7 +2375,7 @@ contains
! Add derivative information for differential tallies. Note that the
! i_nuclide and atom_density arguments do not matter since this is an
! analog estimator.
if (t % deriv /= NONE) then
if (t % deriv() /= NONE) then
call apply_derivative_to_score(p, t, 0, ZERO, SCORE_NU_FISSION, score)
end if
@ -2992,6 +2992,7 @@ contains
integer, intent(in) :: score_bin
real(8), intent(inout) :: score
type(TallyDerivative), pointer :: deriv
integer :: l
logical :: scoring_diff_nuclide
real(8) :: flux_deriv
@ -3004,10 +3005,12 @@ contains
! where (1/f * d_f/d_p) is the (logarithmic) flux derivative and p is the
! perturbated variable.
associate(deriv => tally_derivs(t % deriv))
!associate(deriv => tally_derivs(t % deriv()))
deriv => tally_deriv_c(t % deriv() - 1)
flux_deriv = deriv % flux_deriv
select case (tally_derivs(t % deriv) % variable)
!select case (tally_derivs(t % deriv()) % variable)
select case (deriv % variable)
!=========================================================================
! Density derivative:
@ -3546,7 +3549,7 @@ contains
&analog and collision estimators.")
end select
end select
end associate
!end associate
end subroutine apply_derivative_to_score
!===============================================================================
@ -3558,14 +3561,16 @@ contains
type(Particle), intent(in) :: p
real(8), intent(in) :: distance ! Neutron flight distance
type(TallyDerivative), pointer :: deriv
integer :: i, l
real(8) :: dsig_s, dsig_a, dsig_f
! A void material cannot be perturbed so it will not affect flux derivatives
if (p % material == MATERIAL_VOID) return
do i = 1, size(tally_derivs)
associate(deriv => tally_derivs(i))
do i = 1, n_tally_derivs()
!associate(deriv => tally_derivs(i))
deriv => tally_deriv_c(i - 1)
select case (deriv % variable)
case (DIFF_DENSITY)
@ -3609,7 +3614,7 @@ contains
end if
end associate
end select
end associate
!end associate
end do
end subroutine score_track_derivative
@ -3632,14 +3637,16 @@ contains
subroutine score_collision_derivative(p)
type(Particle), intent(in) :: p
type(TallyDerivative), pointer :: deriv
integer :: i, j, l
real(8) :: dsig_s, dsig_a, dsig_f
! A void material cannot be perturbed so it will not affect flux derivatives
if (p % material == MATERIAL_VOID) return
do i = 1, size(tally_derivs)
associate(deriv => tally_derivs(i))
do i = 1, n_tally_derivs()
!associate(deriv => tally_derivs(i))
deriv => tally_deriv_c(i - 1)
select case (deriv % variable)
case (DIFF_DENSITY)
@ -3702,7 +3709,7 @@ contains
end if
end associate
end select
end associate
!end associate
end do
end subroutine score_collision_derivative
@ -3712,8 +3719,10 @@ contains
subroutine zero_flux_derivs()
integer :: i
do i = 1, size(tally_derivs)
tally_derivs(i) % flux_deriv = ZERO
type(TallyDerivative), pointer :: deriv
do i = 1, n_tally_derivs()
deriv => tally_deriv_c(i - 1)
deriv % flux_deriv = ZERO
end do
end subroutine zero_flux_derivs

6
src/tallies/tally.cpp
View file

@ -3,6 +3,7 @@
#include "openmc/capi.h"
#include "openmc/constants.h"
#include "openmc/error.h"
#include "openmc/tallies/derivative.h"
#include "openmc/tallies/filter.h"
#include "openmc/tallies/filter_energy.h"
#include "openmc/tallies/filter_delayedgroup.h"
@ -175,6 +176,7 @@ free_memory_tally_c()
#pragma omp parallel
{
simulation::filter_matches.clear();
model::tally_derivs.clear();
}
model::tally_filters.clear();
@ -260,6 +262,10 @@ extern "C" {
int tally_get_mesh_filter_c(Tally* tally)
{return tally->mesh_filter_;}
int tally_get_deriv_c(Tally* tally) {return tally->deriv_;}
int tally_set_deriv_c(Tally* tally, int deriv) {tally->deriv_ = deriv;}
}
} // namespace openmc

106
src/tallies/tally_derivative_header.F90
View file

@ -1,104 +1,38 @@
module tally_derivative_header
use constants
use dict_header, only: DictIntInt, EMPTY
use error, only: fatal_error
use nuclide_header, only: nuclide_dict
use string, only: to_str, to_lower
use xml_interface
use, intrinsic :: ISO_C_BINDING
use dict_header, only: DictIntInt
implicit none
private
public :: free_memory_tally_derivative
!===============================================================================
! TALLYDERIVATIVE describes a first-order derivative that can be applied to
! tallies.
!===============================================================================
type, public :: TallyDerivative
integer :: id
integer :: variable
integer :: diff_material
integer :: diff_nuclide
real(8) :: flux_deriv
contains
procedure :: from_xml
type, bind(C), public :: TallyDerivative
integer(C_INT) :: id
integer(C_INT) :: variable
integer(C_INT) :: diff_material
integer(C_INT) :: diff_nuclide
real(C_DOUBLE) :: flux_deriv
end type TallyDerivative
type(TallyDerivative), public, allocatable :: tally_derivs(:)
!$omp threadprivate(tally_derivs)
! Dictionary that maps user IDs to indices in 'tally_derivs'
type(DictIntInt), public :: tally_deriv_dict
contains
interface
function n_tally_derivs() result(n) bind(C)
import C_INT
integer(C_INT) :: n
end function
subroutine from_xml(this, node)
class(TallyDerivative), intent(inout) :: this
type(XMLNode), intent(in) :: node
character(MAX_WORD_LEN) :: temp_str
character(MAX_WORD_LEN) :: word
integer :: val
! Copy the derivative id.
if (check_for_node(node, "id")) then
call get_node_value(node, "id", this % id)
else
call fatal_error("Must specify an ID for <derivative> elements in the&
& tally XML file")
end if
! Make sure the id is > 0.
if (this % id <= 0) then
call fatal_error("<derivative> IDs must be an integer greater than &
&zero")
end if
! Make sure this id has not already been used.
if (tally_deriv_dict % has(this % id)) then
call fatal_error("Two or more <derivative>'s use the same unique &
&ID: " // trim(to_str(this % id)))
end if
! Read the independent variable name.
call get_node_value(node, "variable", temp_str)
temp_str = to_lower(temp_str)
select case(temp_str)
case("density")
this % variable = DIFF_DENSITY
case("nuclide_density")
this % variable = DIFF_NUCLIDE_DENSITY
call get_node_value(node, "nuclide", word)
word = trim(to_lower(word))
val = nuclide_dict % get(word)
if (val == EMPTY) then
call fatal_error("Could not find the nuclide " &
// trim(word) // " specified in derivative " &
// trim(to_str(this % id)) // " in any material.")
end if
this % diff_nuclide = val
case("temperature")
this % variable = DIFF_TEMPERATURE
end select
call get_node_value(node, "material", this % diff_material)
end subroutine from_xml
!===============================================================================
! FREE_MEMORY_TALLY_DERIVATIVE deallocates global arrays defined in this module
!===============================================================================
subroutine free_memory_tally_derivative()
!$omp parallel
if (allocated(tally_derivs)) deallocate(tally_derivs)
!$omp end parallel
end subroutine free_memory_tally_derivative
function tally_deriv_c(i) result(deriv) bind(C)
import C_INT, TallyDerivative
integer(C_INT), value :: i
type(TallyDerivative), pointer :: deriv
end function
end interface
end module tally_derivative_header

31
src/tallies/tally_header.F90
View file

@ -93,9 +93,6 @@ module tally_header
integer :: n_triggers = 0 ! # of triggers
type(TriggerObject), allocatable :: triggers(:) ! Array of triggers
! Index for the TallyDerivative for differential tallies.
integer :: deriv = NONE
contains
procedure :: accumulate => tally_accumulate
procedure :: allocate_results => tally_allocate_results
@ -110,6 +107,8 @@ module tally_header
procedure :: delayedgroup_filter => tally_get_delayedgroup_filter
procedure :: surface_filter => tally_get_surface_filter
procedure :: mesh_filter => tally_get_mesh_filter
procedure :: deriv => tally_get_deriv
procedure :: set_deriv => tally_set_deriv
end type TallyObject
type, public :: TallyContainer
@ -404,6 +403,32 @@ contains
filt = tally_get_mesh_filter_c(this % ptr)
end function
function tally_get_deriv(this) result(deriv)
class(TallyObject) :: this
integer(C_INT) :: deriv
interface
function tally_get_deriv_c(tally) result(deriv) bind(C)
import C_PTR, C_INT
type(C_PTR), value :: tally
integer(C_INT) :: deriv
end function
end interface
deriv = tally_get_deriv_c(this % ptr)
end function
subroutine tally_set_deriv(this, deriv)
class(TallyObject) :: this
integer(C_INT) :: deriv
interface
subroutine tally_set_deriv_c(tally, deriv) bind(C)
import C_PTR, C_INT
type(C_PTR), value :: tally
integer(C_INT), value :: deriv
end subroutine
end interface
call tally_set_deriv_c(this % ptr, deriv)
end subroutine
!===============================================================================
! CONFIGURE_TALLIES initializes several data structures related to tallies. This
! is called after the basic tally data has already been read from the