diff --git a/openmc/capi/core.py b/openmc/capi/core.py index c25146fb5b..ced9390435 100644 --- a/openmc/capi/core.py +++ b/openmc/capi/core.py @@ -166,9 +166,17 @@ def keff(): Mean k-eigenvalue and standard deviation of the mean """ - k = (c_double*2)() - _dll.openmc_get_keff(k) - return tuple(k) + n = openmc.capi.num_realizations() + if n > 3: + # Use the combined estimator if there are enough realizations + k = (c_double*2)() + _dll.openmc_get_keff(k) + return tuple(k) + else: + # Otherwise, return the tracklength estimator + mean = c_double.in_dll(_dll, 'keff').value + std_dev = c_double.in_dll(_dll, 'keff_std').value if n > 1 else np.inf + return (mean, std_dev) def next_batch(): diff --git a/openmc/capi/error.py b/openmc/capi/error.py index c312987bc2..a11d6ea87d 100644 --- a/openmc/capi/error.py +++ b/openmc/capi/error.py @@ -4,39 +4,39 @@ from warnings import warn from . import _dll -class Error(Exception): +class OpenMCError(Exception): """Root exception class for OpenMC.""" -class GeometryError(Error): +class GeometryError(OpenMCError): """Geometry-related error""" -class InvalidIDError(Error): +class InvalidIDError(OpenMCError): """Use of an ID that is invalid.""" -class AllocationError(Error): +class AllocationError(OpenMCError): """Error related to memory allocation.""" -class OutOfBoundsError(Error): +class OutOfBoundsError(OpenMCError): """Index in array out of bounds.""" -class DataError(Error): +class DataError(OpenMCError): """Error relating to nuclear data.""" -class PhysicsError(Error): +class PhysicsError(OpenMCError): """Error relating to performing physics.""" -class InvalidArgumentError(Error): +class InvalidArgumentError(OpenMCError): """Argument passed was invalid.""" -class InvalidTypeError(Error): +class InvalidTypeError(OpenMCError): """Tried to perform an operation on the wrong type.""" @@ -71,4 +71,4 @@ def _error_handler(err, func, args): elif err == errcode('e_warning'): warn(msg) elif err < 0: - raise Exception("Unknown error encountered (code {}).".format(err)) + raise OpenMCError("Unknown error encountered (code {}).".format(err)) diff --git a/openmc/capi/tally.py b/openmc/capi/tally.py index 4a76de41d5..751f0e6031 100644 --- a/openmc/capi/tally.py +++ b/openmc/capi/tally.py @@ -92,12 +92,17 @@ def global_tallies(): n = num_realizations() # Determine mean - mean = sum_ / n + if n > 0: + mean = sum_ / n + else: + mean = sum_.copy() # 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)) + stdev = np.empty_like(mean) + stdev.fill(np.inf) + if n > 1: + stdev[nonzero] = np.sqrt((sum_sq[nonzero]/n - mean[nonzero]**2)/(n - 1)) return list(zip(mean, stdev)) @@ -265,7 +270,7 @@ class Tally(_FortranObjectWithID): def std_dev(self): results = self.results std_dev = np.empty(results.shape[:2]) - std_dev.fill(np.nan) + std_dev.fill(np.inf) n = self.num_realizations if n > 1: diff --git a/src/simulation_header.F90 b/src/simulation_header.F90 index 9719cc72e3..62c66e5f98 100644 --- a/src/simulation_header.F90 +++ b/src/simulation_header.F90 @@ -1,5 +1,7 @@ module simulation_header + use, intrinsic :: ISO_C_BINDING + use bank_header use constants use settings, only: gen_per_batch @@ -33,8 +35,8 @@ module simulation_header ! Temporary k-effective values type(VectorReal) :: k_generation ! single-generation estimates of k - real(8) :: keff = ONE ! average k over active batches - real(8) :: keff_std ! standard deviation of average k + real(C_DOUBLE), bind(C) :: keff = ONE ! average k over active batches + real(C_DOUBLE), bind(C) :: keff_std ! standard deviation of average k real(8) :: k_col_abs = ZERO ! sum over batches of k_collision * k_absorption real(8) :: k_col_tra = ZERO ! sum over batches of k_collision * k_tracklength real(8) :: k_abs_tra = ZERO ! sum over batches of k_absorption * k_tracklength