Make better use of C API functions in cmfd_input and input_xml modules

This commit is contained in:
Paul Romano 2017-08-16 14:06:11 -05:00
parent 2a14926661
commit 4998e1f832
7 changed files with 131 additions and 171 deletions

View file

@ -107,6 +107,14 @@ class TallyView(object):
_dll.openmc_tally_results(self._index, data, shape)
return as_array(data, tuple(shape[::-1]))
@filters.setter
def filters(self, filters):
# Get filter indices as int32_t[]
n = len(filters)
indices = (c_int32*n)(*(f._index for f in filters))
_dll.openmc_tally_set_filters(self._index, n, indices)
@nuclides.setter
def nuclides(self, nuclides):
nucs = (c_char_p * len(nuclides))()

View file

@ -261,14 +261,16 @@ contains
logical :: energy_filters
integer :: i ! loop counter
integer :: n ! size of arrays in mesh specification
integer :: ng ! number of energy groups (default 1)
integer(C_INT32_T) :: ng ! number of energy groups (default 1)
integer :: n_filter ! number of filters
integer :: i_start, i_end
integer :: i_filt_start, i_filt_end
integer(C_INT32_T), allocatable :: filter_indices(:)
integer(C_INT) :: err
integer :: i_filt ! index in filters array
integer :: iarray3(3) ! temp integer array
real(8) :: rarray3(3) ! temp double array
real(C_DOUBLE), allocatable :: energies(:)
type(TallyObject), pointer :: t
type(RegularMesh), pointer :: m
type(XMLNode) :: node_mesh
@ -385,60 +387,39 @@ contains
! Set up mesh filter
i_filt = i_filt_start
allocate(MeshFilter :: filters(i_filt) % obj)
select type (filt => filters(i_filt) % obj)
type is (MeshFilter)
filt % id = i_filt
filt % n_bins = product(m % dimension)
filt % mesh = i_start
! Add filter to dictionary
call filter_dict % add_key(filt % id, i_filt)
end select
err = openmc_filter_set_type(i_filt, C_CHAR_'mesh' // C_NULL_CHAR)
err = openmc_filter_set_id(i_filt, i_filt)
err = openmc_mesh_filter_set_mesh(i_filt, i_start)
if (energy_filters) then
! Read and set incoming energy mesh filter
i_filt = i_filt + 1
allocate(EnergyFilter :: filters(i_filt) % obj)
select type (filt => filters(i_filt) % obj)
type is (EnergyFilter)
filt % id = i_filt
ng = node_word_count(node_mesh, "energy")
filt % n_bins = ng - 1
allocate(filt % bins(ng))
call get_node_array(node_mesh, "energy", filt % bins)
! Add filter to dictionary
call filter_dict % add_key(filt % id, i_filt)
end select
err = openmc_filter_set_type(i_filt, C_CHAR_'energy' // C_NULL_CHAR)
err = openmc_filter_set_id(i_filt, i_filt)
! Get energies and set bins
ng = node_word_count(node_mesh, "energy")
allocate(energies(ng))
call get_node_array(node_mesh, "energy", energies)
err = openmc_energy_filter_set_bins(i_filt, ng, energies)
! Read and set outgoing energy mesh filter
i_filt = i_filt + 1
allocate(EnergyoutFilter :: filters(i_filt) % obj)
select type (filt => filters(i_filt) % obj)
type is (EnergyoutFilter)
filt % id = i_filt
ng = node_word_count(node_mesh, "energy")
filt % n_bins = ng - 1
allocate(filt % bins(ng))
call get_node_array(node_mesh, "energy", filt % bins)
! Add filter to dictionary
call filter_dict % add_key(filt % id, i_filt)
end select
err = openmc_filter_set_type(i_filt, C_CHAR_'energyout' // C_NULL_CHAR)
err = openmc_filter_set_id(i_filt, i_filt)
err = openmc_energy_filter_set_bins(i_filt, ng, energies)
end if
! Duplicate the mesh filter for the mesh current tally since other
! tallies use this filter and we need to change the dimension
i_filt = i_filt + 1
allocate(MeshFilter :: filters(i_filt) % obj)
select type (filt => filters(i_filt) % obj)
type is (MeshFilter)
filt % id = i_filt
! We need to increase the dimension by one since we also need
! currents coming into and out of the boundary mesh cells.
filt % n_bins = product(m % dimension + 1)
filt % mesh = i_start
! Add filter to dictionary
call filter_dict % add_key(filt % id, i_filt)
end select
err = openmc_filter_set_type(i_filt, C_CHAR_'mesh' // C_NULL_CHAR)
err = openmc_filter_set_id(i_filt, i_filt)
err = openmc_mesh_filter_set_mesh(i_filt, i_start)
! We need to increase the dimension by one since we also need
! currents coming into and out of the boundary mesh cells.
filters(i_filt) % obj % n_bins = product(m % dimension + 1)
! Set up surface filter
i_filt = i_filt + 1
@ -476,15 +457,11 @@ contains
call get_node_value(root, "reset", t % reset)
end if
! Set the mesh filter index in the tally find_filter array
n_filter = 1
t % find_filter(FILTER_MESH) = n_filter
! Set the incoming energy mesh filter index in the tally find_filter
! array
n_filter = 1
if (energy_filters) then
n_filter = n_filter + 1
t % find_filter(FILTER_ENERGYIN) = n_filter
end if
! Set number of nucilde bins
@ -507,11 +484,13 @@ contains
t % type = TALLY_VOLUME
! Allocate and set filters
allocate(t % filter(n_filter))
t % filter(1) = i_filt_start
allocate(filter_indices(n_filter))
filter_indices(1) = i_filt_start
if (energy_filters) then
t % filter(2) = i_filt_start + 1
filter_indices(2) = i_filt_start + 1
end if
err = openmc_tally_set_filters(i_start + i - 1, n_filter, filter_indices)
deallocate(filter_indices)
! Allocate scoring bins
allocate(t % score_bins(3))
@ -543,16 +522,17 @@ contains
! array
if (energy_filters) then
n_filter = n_filter + 1
t % find_filter(FILTER_ENERGYOUT) = n_filter
end if
! Allocate and set indices in filters array
allocate(t % filter(n_filter))
t % filter(1) = i_filt_start
allocate(filter_indices(n_filter))
filter_indices(1) = i_filt_start
if (energy_filters) then
t % filter(2) = i_filt_start + 1
t % filter(3) = i_filt_start + 2
filter_indices(2) = i_filt_start + 1
filter_indices(3) = i_filt_start + 2
end if
err = openmc_tally_set_filters(i_start + i - 1, n_filter, filter_indices)
deallocate(filter_indices)
! Allocate macro reactions
allocate(t % score_bins(2))
@ -577,15 +557,16 @@ contains
! Set the surface filter index in the tally find_filter array
n_filter = n_filter + 1
t % find_filter(FILTER_SURFACE) = n_filter
! Allocate and set filters
allocate(t % filter(n_filter))
t % filter(1) = i_filt_end - 1
t % filter(n_filter) = i_filt_end
allocate(filter_indices(n_filter))
filter_indices(1) = i_filt_end - 1
filter_indices(n_filter) = i_filt_end
if (energy_filters) then
t % filter(2) = i_filt_start + 1
filter_indices(2) = i_filt_start + 1
end if
err = openmc_tally_set_filters(i_start + i - 1, n_filter, filter_indices)
deallocate(filter_indices)
! Allocate macro reactions
allocate(t % score_bins(1))

View file

@ -29,7 +29,7 @@ module initialize
use state_point, only: load_state_point
use string, only: to_str, starts_with, ends_with, str_to_int
use summary, only: write_summary
use tally_header, only: TallyObject, configure_tallies
use tally_header, only: TallyObject
use tally_filter
use tally, only: init_tally_routines
@ -116,9 +116,6 @@ contains
end if
if (run_mode /= MODE_PLOTTING) then
! Allocate and setup tally stride, filter_matches, and tally maps
call configure_tallies()
! Set up tally procedure pointers
call init_tally_routines()

View file

@ -3315,11 +3315,7 @@ contains
end select
! Set filter id
f % obj % id = filter_id
! Add filter to dictionary
call filter_dict % add_key(filter_id, i)
err = openmc_filter_set_id(i_start + i - 1, filter_id)
end do READ_FILTERS
! ==========================================================================
@ -3394,64 +3390,28 @@ contains
allocate(temp_filter(n_filter))
if (n_filter > 0) then
call get_node_array(node_tal, "filters", temp_filter)
do j = 1, n_filter
! Get pointer to filter
if (filter_dict % has_key(temp_filter(j))) then
i_filt = filter_dict % get_key(temp_filter(j))
f => filters(i_filt)
else
call fatal_error("Could not find filter " &
// trim(to_str(temp_filter(j))) // " specified on tally " &
// trim(to_str(t % id)))
end if
! Store the index of the filter
temp_filter(j) = i_filt
end do
! Set the filters
err = openmc_tally_set_filters(i_start + i - 1, n_filter, temp_filter)
else
allocate(t % filter(n_filter))
end if
do j = 1, n_filter
! Get pointer to filter
if (filter_dict % has_key(temp_filter(j))) then
i_filt = filter_dict % get_key(temp_filter(j))
f => filters(i_filt)
else
call fatal_error("Could not find filter " &
// trim(to_str(temp_filter(j))) // " specified on tally " &
// trim(to_str(t % id)))
end if
! Set the filter index in the tally find_filter array
select type (filt => f % obj)
type is (DistribcellFilter)
t % find_filter(FILTER_DISTRIBCELL) = j
type is (CellFilter)
t % find_filter(FILTER_CELL) = j
type is (CellFromFilter)
t % find_filter(FILTER_CELLFROM) = j
type is (CellbornFilter)
t % find_filter(FILTER_CELLBORN) = j
type is (MaterialFilter)
t % find_filter(FILTER_MATERIAL) = j
type is (UniverseFilter)
t % find_filter(FILTER_UNIVERSE) = j
type is (SurfaceFilter)
t % find_filter(FILTER_SURFACE) = j
type is (MeshFilter)
t % find_filter(FILTER_MESH) = j
type is (EnergyFilter)
t % find_filter(FILTER_ENERGYIN) = j
type is (EnergyoutFilter)
t % find_filter(FILTER_ENERGYOUT) = j
! Set to analog estimator
t % estimator = ESTIMATOR_ANALOG
type is (DelayedGroupFilter)
t % find_filter(FILTER_DELAYEDGROUP) = j
type is (MuFilter)
t % find_filter(FILTER_MU) = j
! Set to analog estimator
t % estimator = ESTIMATOR_ANALOG
type is (PolarFilter)
t % find_filter(FILTER_POLAR) = j
type is (AzimuthalFilter)
t % find_filter(FILTER_AZIMUTHAL) = j
type is (EnergyFunctionFilter)
t % find_filter(FILTER_ENERGYFUNCTION) = j
end select
! Store the index of the filter
temp_filter(j) = i_filt
end do
! Store the filter indices
call move_alloc(FROM=temp_filter, TO=t % filter)
! =======================================================================
! READ DATA FOR NUCLIDES

View file

@ -21,6 +21,7 @@ module simulation
use state_point, only: write_state_point, write_source_point
use string, only: to_str
use tally, only: accumulate_tallies, setup_active_tallies
use tally_header, only: configure_tallies
use trigger, only: check_triggers
use tracking, only: transport
@ -371,7 +372,11 @@ contains
subroutine initialize_simulation()
! Allocate tally results arrays if they're not allocated yet
call configure_tallies()
!$omp parallel
! Allocate array for microscopic cross section cache
allocate(micro_xs(n_nuclides))
! Allocate array for matching filter bins

View file

@ -182,6 +182,8 @@ contains
if (index >= 1 .and. index <= n_filters) then
if (allocated(filters(index) % obj)) then
filters(index) % obj % id = id
call filter_dict % add_key(id, index)
err = 0
else
err = E_FILTER_NOT_ALLOCATED

View file

@ -108,7 +108,7 @@ module tally_header
contains
procedure :: read_results_hdf5 => tally_read_results_hdf5
procedure :: write_results_hdf5 => tally_write_results_hdf5
procedure :: setup_arrays => tally_setup_arrays
procedure :: allocate_results => tally_allocate_results
end type TallyObject
integer(C_INT32_T), public, bind(C) :: n_tallies = 0 ! # of tallies
@ -210,45 +210,32 @@ contains
end subroutine tally_read_results_hdf5
!===============================================================================
! SETUP_ARRAYS allocates and populates several member arrays of the TallyObject
! derived type, including stride, filter_matches, and results.
! ALLOCATE_RESULTS allocates and initializes the results component of the
! TallyObject derived type
!===============================================================================
subroutine tally_setup_arrays(this)
subroutine tally_allocate_results(this)
class(TallyObject), intent(inout) :: this
integer :: j ! loop index for filters
integer :: n ! temporary stride
integer :: i_filt ! filter index
! Allocate stride
if (allocated(this % filter)) then
if (allocated(this % stride)) deallocate(this % stride)
allocate(this % stride(size(this % filter)))
! The filters are traversed in opposite order so that the last filter has
! the shortest stride in memory and the first filter has the largest
! stride
n = 1
STRIDE: do j = size(this % filter), 1, -1
i_filt = this % filter(j)
this % stride(j) = n
n = n * filters(i_filt) % obj % n_bins
end do STRIDE
this % n_filter_bins = n
else
this % n_filter_bins = 1
end if
! Set total number of filter and scoring bins
this % total_score_bins = this % n_score_bins * this % n_nuclide_bins
! Allocate results array
if (allocated(this % results)) deallocate(this % results)
allocate(this % results(3, this % total_score_bins, this % n_filter_bins))
if (allocated(this % results)) then
! If results was already allocated but shape is wrong, then reallocate it
! to the correct shape
if (this % total_score_bins /= size(this % results, 2) .or. &
this % n_filter_bins /= size(this % results, 3)) then
deallocate(this % results)
allocate(this % results(3, this % total_score_bins, this % n_filter_bins))
end if
else
allocate(this % results(3, this % total_score_bins, this % n_filter_bins))
end if
! Initialize results array to zero
this % results(:,:,:) = ZERO
end subroutine tally_setup_arrays
end subroutine tally_allocate_results
!===============================================================================
! CONFIGURE_TALLIES initializes several data structures related to tallies. This
@ -260,12 +247,14 @@ contains
integer :: i
! Allocate global tallies
allocate(global_tallies(3, N_GLOBAL_TALLIES))
! Allocate and initialize global tallies
if (.not. allocated(global_tallies)) then
allocate(global_tallies(3, N_GLOBAL_TALLIES))
end if
global_tallies(:,:) = ZERO
do i = 1, n_tallies
call tallies(i) % setup_arrays()
call tallies(i) % allocate_results()
end do
end subroutine configure_tallies
@ -412,24 +401,25 @@ contains
integer(C_INT32_T), intent(in) :: filter_indices(n)
integer(C_INT) :: err
integer :: i
integer :: j
integer :: i ! index in t % filter/stride
integer :: j ! index in t % find_filter
integer :: k ! index in global filters array
integer :: stride ! filter stride
err = 0
if (index >= 1 .and. index <= n_tallies) then
associate (t => tallies(index))
if (allocated(t % filter)) deallocate(t % filter)
allocate(t % filter(n))
t % find_filter(:) = 0
do i = 1, n
if (filter_indices(i) < 1 .or. filter_indices(i) > n_filters) then
k = filter_indices(i)
if (k < 1 .or. k > n_filters) then
err = E_OUT_OF_BOUNDS
exit
end if
! Set the filter index in the tally find_filter array
select type (filt => filters(i) % obj)
select type (filt => filters(k) % obj)
type is (DistribcellFilter)
j = FILTER_DISTRIBCELL
type is (CellFilter)
@ -465,7 +455,28 @@ contains
end select
t % find_filter(j) = i
end do
if (err == 0) t % filter(:) = filter_indices
if (err == 0) then
if (allocated(t % filter)) deallocate(t % filter)
if (allocated(t % stride)) deallocate(t % stride)
allocate(t % filter(n), t % stride(n))
! Filters are traversed in reverse so that the last filter has the
! shortest stride in memory and the first filter has the largest stride
stride = 1
do i = n, 1, -1
! Set filter and stride
k = filter_indices(i)
t % filter(i) = k
t % stride(i) = stride
! Multiply stride by number of bins in this filter
stride = stride * filters(k) % obj % n_bins
end do
! Set total number of filter bins
t % n_filter_bins = stride
end if
end associate
else
err = E_OUT_OF_BOUNDS
@ -509,8 +520,6 @@ contains
end select
end do
call t % setup_arrays()
err = 0
end associate
else
@ -674,8 +683,6 @@ contains
end select
end do
call t % setup_arrays()
err = 0
end associate
else