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

Add openmc.capi.global_tallies()

Browse source
This commit is contained in:
Paul Romano 2017-11-03 21:45:49 -05:00
parent ac565724dc
commit f14957940d
2 changed files with 43 additions and 5 deletions
Download patch file Download diff file Expand all files Collapse all files

33
openmc/capi/tally.py
View file

@ -2,6 +2,7 @@ from collections import Mapping
from ctypes import c_int, c_int32, c_double, c_char_p, POINTER
from weakref import WeakValueDictionary
import numpy as np
from numpy.ctypeslib import as_array
from openmc.data.reaction import REACTION_NAME
@ -11,15 +12,18 @@ from .error import _error_handler, AllocationError, InvalidIDError
from .filter import _get_filter
__all__ = ['Tally', 'tallies']
__all__ = ['Tally', 'tallies', 'global_tallies']
# Tally functions
_dll.openmc_get_tally_index.argtypes = [c_int32, POINTER(c_int32)]
_dll.openmc_get_tally_index.restype = c_int
_dll.openmc_get_tally_index.errcheck = _error_handler
_dll.openmc_extend_tallies.argtypes = [c_int32, POINTER(c_int32), POINTER(c_int32)]
_dll.openmc_extend_tallies.restype = c_int
_dll.openmc_extend_tallies.errcheck = _error_handler
_dll.openmc_get_tally_index.argtypes = [c_int32, POINTER(c_int32)]
_dll.openmc_get_tally_index.restype = c_int
_dll.openmc_get_tally_index.errcheck = _error_handler
_dll.openmc_global_tallies.argtypes = [POINTER(POINTER(c_double))]
_dll.openmc_global_tallies.restype = c_int
_dll.openmc_global_tallies.errcheck = _error_handler
_dll.openmc_tally_get_id.argtypes = [c_int32, POINTER(c_int32)]
_dll.openmc_tally_get_id.restype = c_int
_dll.openmc_tally_get_id.errcheck = _error_handler
@ -65,6 +69,27 @@ _SCORES = {
}
def global_tallies():
ptr = POINTER(c_double)()
_dll.openmc_global_tallies(ptr)
array = as_array(ptr, (4, 3))
# Get sum, sum-of-squares, and number of realizations
sum_ = array[:, 1]
sum_sq = array[:, 2]
n = c_int32.in_dll(_dll, 'n_realizations').value
# Determine mean
mean = sum_ / n
# Determine standard deviation
nonzero = np.abs(mean) > 0
stdev = np.zeros_like(mean)
stdev[nonzero] = np.sqrt((sum_sq[nonzero]/n - mean[nonzero]**2)/(n - 1))
return list(zip(mean, stdev))
class Tally(_FortranObjectWithID):
"""Tally stored internally.

15
src/tallies/tally_header.F90
View file

@ -146,7 +146,7 @@ module tally_header
type(VectorInt), public :: active_surface_tallies
! Normalization for statistics
integer, public :: n_realizations = 0 ! # of independent realizations
integer(C_INT32_T), public, bind(C) :: n_realizations = 0 ! # of independent realizations
real(8), public :: total_weight ! total starting particle weight in realization
contains
@ -470,6 +470,19 @@ contains
end function openmc_get_tally_index
function openmc_global_tallies(ptr) result(err) bind(C)
type(C_PTR), intent(out) :: ptr
integer(C_INT) :: err
if (.not. allocated(global_tallies)) then
err = E_ALLOCATE
else
err = 0
ptr = C_LOC(global_tallies)
end if
end function openmc_global_tallies
function openmc_tally_get_id(index, id) result(err) bind(C)
! Return the ID of a tally
integer(C_INT32_T), value :: index