mirror of
https://github.com/openmc-dev/openmc.git
synced 2026-07-27 13:45:36 -04:00
Merge branch 'dlopen_source' of https://github.com/makeclean/openmc into dlopen_source
This commit is contained in:
Andrew Davis
commit
e5e8e1823d
15 changed files with 384 additions and 185 deletions
|
|
@ -1,5 +1,8 @@
|
|||
get_filename_component(OpenMC_CMAKE_DIR "${CMAKE_CURRENT_LIST_FILE}" DIRECTORY)
|
||||
|
||||
find_package(fmt REQUIRED HINTS ${OpenMC_CMAKE_DIR}/../fmt)
|
||||
find_package(gsl-lite REQUIRED HINTS ${OpenMC_CMAKE_DIR}/../gsl-lite)
|
||||
find_package(pugixml REQUIRED HINTS ${OpenMC_CMAKE_DIR}/../pugixml)
|
||||
find_package(xtl REQUIRED HINTS ${OpenMC_CMAKE_DIR}/../xtl)
|
||||
find_package(xtensor REQUIRED HINTS ${OpenMC_CMAKE_DIR}/../xtensor)
|
||||
|
||||
|
|
|
|||
|
|
@ -33,12 +33,11 @@ material compositions over time. Each method appears as a different class.
|
|||
For example, :class:`openmc.deplete.CECMIntegrator` runs a depletion calculation
|
||||
using the CE/CM algorithm (deplete over a timestep using the middle-of-step
|
||||
reaction rates). An instance of :class:`openmc.deplete.Operator` is passed to
|
||||
one of these functions along with the power level and timesteps::
|
||||
one of these functions along with the timesteps and power level::
|
||||
|
||||
power = 1200.0e6
|
||||
days = 24*60*60
|
||||
timesteps = [10.0*days, 10.0*days, 10.0*days]
|
||||
openmc.deplete.CECMIntegrator(op, power, timesteps).integrate()
|
||||
power = 1200.0e6 # watts
|
||||
timesteps = [10.0, 10.0, 10.0] # days
|
||||
openmc.deplete.CECMIntegrator(op, timesteps, power, timestep_units='d').integrate()
|
||||
|
||||
The coupled transport-depletion problem is executed, and once it is done a
|
||||
``depletion_results.h5`` file is written. The results can be analyzed using the
|
||||
|
|
@ -67,7 +66,7 @@ the energy deposited during a transport calculation will be lower than expected.
|
|||
This causes the reaction rates to be over-adjusted to hit the user-specific power,
|
||||
or power density, leading to an over-depletion of burnable materials.
|
||||
|
||||
There are some remedies. First, the fission Q values can be directly set in a
|
||||
There are some remedies. First, the fission Q values can be directly set in a
|
||||
variety of ways. This requires knowing what the total fission energy release should
|
||||
be, including indirect components. Some examples are provided below::
|
||||
|
||||
|
|
@ -99,11 +98,11 @@ Local Spectra and Repeated Materials
|
|||
------------------------------------
|
||||
|
||||
It is not uncommon to explicitly create a single burnable material across many locations.
|
||||
From a pure transport perspective, there is nothing wrong with creating a single
|
||||
From a pure transport perspective, there is nothing wrong with creating a single
|
||||
3.5 wt.% enriched fuel ``fuel_3``, and placing that fuel in every fuel pin in an assembly
|
||||
or even full core problem. This certainly expedites the model making process, but can pose
|
||||
issues with depletion.
|
||||
Under this setup, :mod:`openmc.deplete` will deplete a single ``fuel_3`` material using
|
||||
issues with depletion.
|
||||
Under this setup, :mod:`openmc.deplete` will deplete a single ``fuel_3`` material using
|
||||
a single set of reaction rates, and produce a single new composition for the next time
|
||||
step. This can be problematic if the same ``fuel_3`` is used in very different regions
|
||||
of the problem.
|
||||
|
|
|
|||
|
|
@ -1,25 +1,7 @@
|
|||
import openmc
|
||||
import openmc.deplete
|
||||
import numpy as np
|
||||
import matplotlib.pyplot as plt
|
||||
|
||||
###############################################################################
|
||||
# Simulation Input File Parameters
|
||||
###############################################################################
|
||||
|
||||
# OpenMC simulation parameters
|
||||
batches = 100
|
||||
inactive = 10
|
||||
particles = 1000
|
||||
|
||||
# Depletion simulation parameters
|
||||
time_step = 1*24*60*60 # s
|
||||
final_time = 5*24*60*60 # s
|
||||
time_steps = np.full(final_time // time_step, time_step)
|
||||
|
||||
chain_file = './chain_simple.xml'
|
||||
power = 174 # W/cm, for 2D simulations only (use W for 3D)
|
||||
|
||||
###############################################################################
|
||||
# Load previous simulation results
|
||||
###############################################################################
|
||||
|
|
@ -37,31 +19,34 @@ previous_results = openmc.deplete.ResultsList("depletion_results.h5")
|
|||
###############################################################################
|
||||
|
||||
# Instantiate a Settings object, set all runtime parameters
|
||||
settings_file = openmc.Settings()
|
||||
settings_file.batches = batches
|
||||
settings_file.inactive = inactive
|
||||
settings_file.particles = particles
|
||||
settings = openmc.Settings()
|
||||
settings.batches = 100
|
||||
settings.inactive = 10
|
||||
settings.particles = 10000
|
||||
|
||||
# Create an initial uniform spatial source distribution over fissionable zones
|
||||
bounds = [-0.62992, -0.62992, -1, 0.62992, 0.62992, 1]
|
||||
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:], only_fissionable=True)
|
||||
settings_file.source = openmc.source.Source(space=uniform_dist)
|
||||
settings.source = openmc.source.Source(space=uniform_dist)
|
||||
|
||||
entropy_mesh = openmc.RegularMesh()
|
||||
entropy_mesh.lower_left = [-0.39218, -0.39218, -1.e50]
|
||||
entropy_mesh.upper_right = [0.39218, 0.39218, 1.e50]
|
||||
entropy_mesh.dimension = [10, 10, 1]
|
||||
settings_file.entropy_mesh = entropy_mesh
|
||||
settings.entropy_mesh = entropy_mesh
|
||||
|
||||
###############################################################################
|
||||
# Initialize and run depletion calculation
|
||||
###############################################################################
|
||||
|
||||
op = openmc.deplete.Operator(geometry, settings_file, chain_file,
|
||||
previous_results)
|
||||
# Create depletion "operator"
|
||||
chain_file = './chain_simple.xml'
|
||||
op = openmc.deplete.Operator(geometry, settings, chain_file, previous_results)
|
||||
|
||||
# Perform simulation using the predictor algorithm
|
||||
integrator = openmc.deplete.PredictorIntegrator(op, time_steps, power)
|
||||
time_steps = [1.0, 1.0, 1.0, 1.0, 1.0] # days
|
||||
power = 174 # W/cm, for 2D simulations only (use W for 3D)
|
||||
integrator = openmc.deplete.PredictorIntegrator(op, time_steps, power, timestep_units='d')
|
||||
integrator.integrate()
|
||||
|
||||
###############################################################################
|
||||
|
|
@ -77,27 +62,28 @@ time, keff = results.get_eigenvalue()
|
|||
# Obtain U235 concentration as a function of time
|
||||
time, n_U235 = results.get_atoms('1', 'U235')
|
||||
|
||||
# Obtain Xe135 absorption as a function of time
|
||||
time, Xe_gam = results.get_reaction_rate('1', 'Xe135', '(n,gamma)')
|
||||
# Obtain Xe135 capture reaction rate as a function of time
|
||||
time, Xe_capture = results.get_reaction_rate('1', 'Xe135', '(n,gamma)')
|
||||
|
||||
###############################################################################
|
||||
# Generate plots
|
||||
###############################################################################
|
||||
|
||||
days = 24*60*60
|
||||
plt.figure()
|
||||
plt.plot(time/(24*60*60), keff, label="K-effective")
|
||||
plt.plot(time/days, keff, label="K-effective")
|
||||
plt.xlabel("Time (days)")
|
||||
plt.ylabel("Keff")
|
||||
plt.show()
|
||||
|
||||
plt.figure()
|
||||
plt.plot(time/(24*60*60), n_U235, label="U 235")
|
||||
plt.plot(time/days, n_U235, label="U 235")
|
||||
plt.xlabel("Time (days)")
|
||||
plt.ylabel("n U5 (-)")
|
||||
plt.show()
|
||||
|
||||
plt.figure()
|
||||
plt.plot(time/(24*60*60), Xe_gam, label="Xe135 absorption")
|
||||
plt.plot(time/days, Xe_capture, label="Xe135 capture")
|
||||
plt.xlabel("Time (days)")
|
||||
plt.ylabel("RR (-)")
|
||||
plt.show()
|
||||
|
|
|
|||
|
|
@ -1,47 +1,31 @@
|
|||
from math import pi
|
||||
|
||||
import openmc
|
||||
import openmc.deplete
|
||||
import numpy as np
|
||||
import matplotlib.pyplot as plt
|
||||
|
||||
###############################################################################
|
||||
# Simulation Input File Parameters
|
||||
###############################################################################
|
||||
|
||||
# OpenMC simulation parameters
|
||||
batches = 100
|
||||
inactive = 10
|
||||
particles = 1000
|
||||
|
||||
# Depletion simulation parameters
|
||||
time_step = 1*24*60*60 # s
|
||||
final_time = 5*24*60*60 # s
|
||||
time_steps = np.full(final_time // time_step, time_step)
|
||||
chain_file = './chain_simple.xml'
|
||||
power = 174 # W/cm, for 2D simulations only (use W for 3D)
|
||||
|
||||
###############################################################################
|
||||
# Define materials
|
||||
###############################################################################
|
||||
|
||||
# Instantiate some Materials and register the appropriate Nuclides
|
||||
uo2 = openmc.Material(material_id=1, name='UO2 fuel at 2.4% wt enrichment')
|
||||
uo2 = openmc.Material(name='UO2 fuel at 2.4% wt enrichment')
|
||||
uo2.set_density('g/cm3', 10.29769)
|
||||
uo2.add_element('U', 1., enrichment=2.4)
|
||||
uo2.add_element('O', 2.)
|
||||
uo2.depletable = True
|
||||
|
||||
helium = openmc.Material(material_id=2, name='Helium for gap')
|
||||
helium = openmc.Material(name='Helium for gap')
|
||||
helium.set_density('g/cm3', 0.001598)
|
||||
helium.add_element('He', 2.4044e-4)
|
||||
|
||||
zircaloy = openmc.Material(material_id=3, name='Zircaloy 4')
|
||||
zircaloy = openmc.Material(name='Zircaloy 4')
|
||||
zircaloy.set_density('g/cm3', 6.55)
|
||||
zircaloy.add_element('Sn', 0.014 , 'wo')
|
||||
zircaloy.add_element('Sn', 0.014, 'wo')
|
||||
zircaloy.add_element('Fe', 0.00165, 'wo')
|
||||
zircaloy.add_element('Cr', 0.001 , 'wo')
|
||||
zircaloy.add_element('Cr', 0.001, 'wo')
|
||||
zircaloy.add_element('Zr', 0.98335, 'wo')
|
||||
|
||||
borated_water = openmc.Material(material_id=4, name='Borated water')
|
||||
borated_water = openmc.Material(name='Borated water')
|
||||
borated_water.set_density('g/cm3', 0.740582)
|
||||
borated_water.add_element('B', 4.0e-5)
|
||||
borated_water.add_element('H', 5.0e-2)
|
||||
|
|
@ -52,87 +36,62 @@ borated_water.add_s_alpha_beta('c_H_in_H2O')
|
|||
# Create geometry
|
||||
###############################################################################
|
||||
|
||||
# Instantiate ZCylinder surfaces
|
||||
fuel_or = openmc.ZCylinder(surface_id=1, x0=0, y0=0, r=0.39218, name='Fuel OR')
|
||||
clad_ir = openmc.ZCylinder(surface_id=2, x0=0, y0=0, r=0.40005, name='Clad IR')
|
||||
clad_or = openmc.ZCylinder(surface_id=3, x0=0, y0=0, r=0.45720, name='Clad OR')
|
||||
left = openmc.XPlane(surface_id=4, x0=-0.62992, name='left')
|
||||
right = openmc.XPlane(surface_id=5, x0=0.62992, name='right')
|
||||
bottom = openmc.YPlane(surface_id=6, y0=-0.62992, name='bottom')
|
||||
top = openmc.YPlane(surface_id=7, y0=0.62992, name='top')
|
||||
# Define surfaces
|
||||
pitch = 1.25984
|
||||
fuel_or = openmc.ZCylinder(r=0.39218, name='Fuel OR')
|
||||
clad_ir = openmc.ZCylinder(r=0.40005, name='Clad IR')
|
||||
clad_or = openmc.ZCylinder(r=0.45720, name='Clad OR')
|
||||
box = openmc.model.rectangular_prism(pitch, pitch, boundary_type='reflective')
|
||||
|
||||
left.boundary_type = 'reflective'
|
||||
right.boundary_type = 'reflective'
|
||||
top.boundary_type = 'reflective'
|
||||
bottom.boundary_type = 'reflective'
|
||||
# Define cells
|
||||
fuel = openmc.Cell(fill=uo2, region=-fuel_or)
|
||||
gap = openmc.Cell(fill=helium, region=+fuel_or & -clad_ir)
|
||||
clad = openmc.Cell(fill=zircaloy, region=+clad_ir & -clad_or)
|
||||
water = openmc.Cell(fill=borated_water, region=+clad_or & box)
|
||||
|
||||
# Instantiate Cells
|
||||
fuel = openmc.Cell(cell_id=1, name='cell 1')
|
||||
gap = openmc.Cell(cell_id=2, name='cell 2')
|
||||
clad = openmc.Cell(cell_id=3, name='cell 3')
|
||||
water = openmc.Cell(cell_id=4, name='cell 4')
|
||||
|
||||
# Use surface half-spaces to define regions
|
||||
fuel.region = -fuel_or
|
||||
gap.region = +fuel_or & -clad_ir
|
||||
clad.region = +clad_ir & -clad_or
|
||||
water.region = +clad_or & +left & -right & +bottom & -top
|
||||
|
||||
# Register Materials with Cells
|
||||
fuel.fill = uo2
|
||||
gap.fill = helium
|
||||
clad.fill = zircaloy
|
||||
water.fill = borated_water
|
||||
|
||||
# Instantiate Universe
|
||||
root = openmc.Universe(universe_id=0, name='root universe')
|
||||
|
||||
# Register Cells with Universe
|
||||
root.add_cells([fuel, gap, clad, water])
|
||||
|
||||
# Instantiate a Geometry, register the root Universe
|
||||
geometry = openmc.Geometry(root)
|
||||
# Define overall geometry
|
||||
geometry = openmc.Geometry([fuel, gap, clad, water])
|
||||
|
||||
###############################################################################
|
||||
# Set volumes of depletable materials
|
||||
###############################################################################
|
||||
|
||||
# Compute cell areas
|
||||
area = {}
|
||||
area[fuel] = np.pi * fuel_or.coefficients['r'] ** 2
|
||||
|
||||
# Set materials volume for depletion. Set to an area for 2D simulations
|
||||
uo2.volume = area[fuel]
|
||||
# Set material volume for depletion. For 2D simulations, this should be an area.
|
||||
uo2.volume = pi * fuel_or.r**2
|
||||
|
||||
###############################################################################
|
||||
# Transport calculation settings
|
||||
###############################################################################
|
||||
|
||||
# Instantiate a Settings object, set all runtime parameters, and export to XML
|
||||
settings_file = openmc.Settings()
|
||||
settings_file.batches = batches
|
||||
settings_file.inactive = inactive
|
||||
settings_file.particles = particles
|
||||
settings = openmc.Settings()
|
||||
settings.batches = 100
|
||||
settings.inactive = 10
|
||||
settings.particles = 1000
|
||||
|
||||
# Create an initial uniform spatial source distribution over fissionable zones
|
||||
bounds = [-0.62992, -0.62992, -1, 0.62992, 0.62992, 1]
|
||||
uniform_dist = openmc.stats.Box(bounds[:3], bounds[3:], only_fissionable=True)
|
||||
settings_file.source = openmc.source.Source(space=uniform_dist)
|
||||
settings.source = openmc.source.Source(space=uniform_dist)
|
||||
|
||||
entropy_mesh = openmc.RegularMesh()
|
||||
entropy_mesh.lower_left = [-0.39218, -0.39218, -1.e50]
|
||||
entropy_mesh.upper_right = [0.39218, 0.39218, 1.e50]
|
||||
entropy_mesh.dimension = [10, 10, 1]
|
||||
settings_file.entropy_mesh = entropy_mesh
|
||||
settings.entropy_mesh = entropy_mesh
|
||||
|
||||
###############################################################################
|
||||
# Initialize and run depletion calculation
|
||||
###############################################################################
|
||||
|
||||
op = openmc.deplete.Operator(geometry, settings_file, chain_file)
|
||||
# Create depletion "operator"
|
||||
chain_file = './chain_simple.xml'
|
||||
op = openmc.deplete.Operator(geometry, settings, chain_file)
|
||||
|
||||
# Perform simulation using the predictor algorithm
|
||||
integrator = openmc.deplete.PredictorIntegrator(op, time_steps, power)
|
||||
time_steps = [1.0, 1.0, 1.0, 1.0, 1.0] # days
|
||||
power = 174 # W/cm, for 2D simulations only (use W for 3D)
|
||||
integrator = openmc.deplete.PredictorIntegrator(op, time_steps, power, timestep_units='d')
|
||||
integrator.integrate()
|
||||
|
||||
###############################################################################
|
||||
|
|
@ -148,27 +107,28 @@ time, keff = results.get_eigenvalue()
|
|||
# Obtain U235 concentration as a function of time
|
||||
time, n_U235 = results.get_atoms('1', 'U235')
|
||||
|
||||
# Obtain Xe135 absorption as a function of time
|
||||
time, Xe_gam = results.get_reaction_rate('1', 'Xe135', '(n,gamma)')
|
||||
# Obtain Xe135 capture reaction rate as a function of time
|
||||
time, Xe_capture = results.get_reaction_rate('1', 'Xe135', '(n,gamma)')
|
||||
|
||||
###############################################################################
|
||||
# Generate plots
|
||||
###############################################################################
|
||||
|
||||
days = 24*60*60
|
||||
plt.figure()
|
||||
plt.plot(time/(24*60*60), keff, label="K-effective")
|
||||
plt.plot(time/days, keff, label="K-effective")
|
||||
plt.xlabel("Time (days)")
|
||||
plt.ylabel("Keff")
|
||||
plt.show()
|
||||
|
||||
plt.figure()
|
||||
plt.plot(time/(24*60*60), n_U235, label="U 235")
|
||||
plt.plot(time/days, n_U235, label="U235")
|
||||
plt.xlabel("Time (days)")
|
||||
plt.ylabel("n U5 (-)")
|
||||
plt.show()
|
||||
|
||||
plt.figure()
|
||||
plt.plot(time/(24*60*60), Xe_gam, label="Xe135 absorption")
|
||||
plt.plot(time/days, Xe_capture, label="Xe135 capture")
|
||||
plt.xlabel("Time (days)")
|
||||
plt.ylabel("RR (-)")
|
||||
plt.show()
|
||||
|
|
|
|||
|
|
@ -34,12 +34,6 @@ namespace openmc {
|
|||
// use to store the bins for delayed group tallies.
|
||||
constexpr int MAX_DELAYED_GROUPS {8};
|
||||
|
||||
// Maximum number of lost particles
|
||||
constexpr int MAX_LOST_PARTICLES {10};
|
||||
|
||||
// Maximum number of lost particles, relative to the total number of particles
|
||||
constexpr double REL_MAX_LOST_PARTICLES {1.0e-6};
|
||||
|
||||
constexpr double CACHE_INVALID {-1.0};
|
||||
|
||||
//==============================================================================
|
||||
|
|
|
|||
|
|
@ -64,10 +64,12 @@ extern std::string path_source_library; //!< path to the source shared object
|
|||
extern std::string path_sourcepoint; //!< path to a source file
|
||||
extern "C" std::string path_statepoint; //!< path to a statepoint file
|
||||
|
||||
extern "C" int32_t n_batches; //!< number of (inactive+active) batches
|
||||
extern "C" int32_t n_inactive; //!< number of inactive batches
|
||||
extern "C" int32_t gen_per_batch; //!< number of generations per batch
|
||||
extern "C" int64_t n_particles; //!< number of particles per generation
|
||||
extern "C" int32_t n_batches; //!< number of (inactive+active) batches
|
||||
extern "C" int32_t n_inactive; //!< number of inactive batches
|
||||
extern "C" int32_t max_lost_particles; //!< maximum number of lost particles
|
||||
extern double rel_max_lost_particles; //!< maximum number of lost particles, relative to the total number of particles
|
||||
extern "C" int32_t gen_per_batch; //!< number of generations per batch
|
||||
extern "C" int64_t n_particles; //!< number of particles per generation
|
||||
|
||||
|
||||
extern int64_t max_particles_in_flight; //!< Max num. event-based particles in flight
|
||||
|
|
|
|||
|
|
@ -31,6 +31,10 @@ __all__ = [
|
|||
"Integrator", "SIIntegrator", "DepSystemSolver"]
|
||||
|
||||
|
||||
_SECONDS_PER_MINUTE = 60
|
||||
_SECONDS_PER_HOUR = 60*60
|
||||
_SECONDS_PER_DAY = 24*60*60
|
||||
|
||||
OperatorResult = namedtuple('OperatorResult', ['k', 'rates'])
|
||||
OperatorResult.__doc__ = """\
|
||||
Result of applying transport operator
|
||||
|
|
@ -597,9 +601,11 @@ class Integrator(ABC):
|
|||
----------
|
||||
operator : openmc.deplete.TransportOperator
|
||||
Operator to perform transport simulations
|
||||
timesteps : iterable of float
|
||||
Array of timesteps in units of [s]. Note that values are not
|
||||
cumulative.
|
||||
timesteps : iterable of float or iterable of tuple
|
||||
Array of timesteps. Note that values are not cumulative. The units are
|
||||
specified by the `timestep_units` argument when `timesteps` is an
|
||||
iterable of float. Alternatively, units can be specified for each step
|
||||
by passing an iterable of (value, unit) tuples.
|
||||
power : float or iterable of float, optional
|
||||
Power of the reactor in [W]. A single value indicates that
|
||||
the power is constant over all timesteps. An iterable
|
||||
|
|
@ -612,6 +618,11 @@ class Integrator(ABC):
|
|||
Power density of the reactor in [W/gHM]. It is multiplied by
|
||||
initial heavy metal inventory to get total power if ``power``
|
||||
is not speficied.
|
||||
timestep_units : {'s', 'min', 'h', 'd', 'MWd/kg'}
|
||||
Units for values specified in the `timesteps` argument. 's' means
|
||||
seconds, 'min' means minutes, 'h' means hours, and 'MWd/kg' indicates
|
||||
that the values are given in burnup (MW-d of energy deposited per
|
||||
kilogram of initial heavy metal).
|
||||
|
||||
Attributes
|
||||
----------
|
||||
|
|
@ -625,7 +636,8 @@ class Integrator(ABC):
|
|||
Power of the reactor in [W] for each interval in :attr:`timesteps`
|
||||
"""
|
||||
|
||||
def __init__(self, operator, timesteps, power=None, power_density=None):
|
||||
def __init__(self, operator, timesteps, power=None, power_density=None,
|
||||
timestep_units='s'):
|
||||
# Check number of stages previously used
|
||||
if operator.prev_res is not None:
|
||||
res = operator.prev_res[-1]
|
||||
|
|
@ -638,27 +650,59 @@ class Integrator(ABC):
|
|||
self._num_stages))
|
||||
self.operator = operator
|
||||
self.chain = operator.chain
|
||||
if not isinstance(timesteps, Iterable):
|
||||
self.timesteps = [timesteps]
|
||||
else:
|
||||
self.timesteps = timesteps
|
||||
|
||||
# Determine power and normalize units to W
|
||||
if power is None:
|
||||
if power_density is None:
|
||||
raise ValueError("Either power or power density must be set")
|
||||
if not isinstance(power_density, Iterable):
|
||||
power = power_density * operator.heavy_metal
|
||||
else:
|
||||
power = [p * operator.heavy_metal for p in power_density]
|
||||
|
||||
power = [p*operator.heavy_metal for p in power_density]
|
||||
if not isinstance(power, Iterable):
|
||||
# Ensure that power is single value if that is the case
|
||||
power = [power] * len(self.timesteps)
|
||||
elif len(power) != len(self.timesteps):
|
||||
raise ValueError(
|
||||
"Number of time steps != number of powers. {} vs {}".format(
|
||||
len(self.timesteps), len(power)))
|
||||
power = [power] * len(timesteps)
|
||||
|
||||
self.power = power
|
||||
if len(power) != len(timesteps):
|
||||
raise ValueError(
|
||||
"Number of time steps ({}) != number of powers ({})".format(
|
||||
len(timesteps), len(power)))
|
||||
|
||||
# Get list of times / units
|
||||
if isinstance(timesteps[0], Iterable):
|
||||
times, units = zip(*timesteps)
|
||||
else:
|
||||
times = timesteps
|
||||
units = [timestep_units] * len(timesteps)
|
||||
|
||||
# Determine number of seconds for each timestep
|
||||
seconds = []
|
||||
for time, unit, watts in zip(times, units, power):
|
||||
# Make sure values passed make sense
|
||||
check_type('timestep', time, Real)
|
||||
check_greater_than('timestep', time, 0.0, False)
|
||||
check_type('timestep units', unit, str)
|
||||
check_type('power', watts, Real)
|
||||
check_greater_than('power', watts, 0.0, True)
|
||||
|
||||
if unit in ('s', 'sec'):
|
||||
seconds.append(time)
|
||||
elif unit in ('min', 'minute'):
|
||||
seconds.append(time*_SECONDS_PER_MINUTE)
|
||||
elif unit in ('h', 'hr', 'hour'):
|
||||
seconds.append(time*_SECONDS_PER_HOUR)
|
||||
elif unit in ('d', 'day'):
|
||||
seconds.append(time*_SECONDS_PER_DAY)
|
||||
elif unit.lower() == 'mwd/kg':
|
||||
watt_days_per_kg = 1e6*time
|
||||
kilograms = 1e-3*operator.heavy_metal
|
||||
days = watt_days_per_kg * kilograms / watts
|
||||
seconds.append(days*_SECONDS_PER_DAY)
|
||||
else:
|
||||
raise ValueError("Invalid timestep unit '{}'".format(unit))
|
||||
|
||||
self.timesteps = asarray(seconds)
|
||||
self.power = asarray(power)
|
||||
|
||||
@abstractmethod
|
||||
def __call__(self, conc, rates, dt, power, i):
|
||||
|
|
@ -772,9 +816,11 @@ class SIIntegrator(Integrator):
|
|||
----------
|
||||
operator : openmc.deplete.TransportOperator
|
||||
The operator object to simulate on.
|
||||
timesteps : iterable of float
|
||||
Array of timesteps in units of [s]. Note that values are not
|
||||
cumulative.
|
||||
timesteps : iterable of float or iterable of tuple
|
||||
Array of timesteps. Note that values are not cumulative. The units are
|
||||
specified by the `timestep_units` argument when `timesteps` is an
|
||||
iterable of float. Alternatively, units can be specified for each step
|
||||
by passing an iterable of (value, unit) tuples.
|
||||
power : float or iterable of float, optional
|
||||
Power of the reactor in [W]. A single value indicates that
|
||||
the power is constant over all timesteps. An iterable
|
||||
|
|
@ -787,6 +833,11 @@ class SIIntegrator(Integrator):
|
|||
Power density of the reactor in [W/gHM]. It is multiplied by
|
||||
initial heavy metal inventory to get total power if ``power``
|
||||
is not speficied.
|
||||
timestep_units : {'s', 'min', 'h', 'd', 'MWd/kg'}
|
||||
Units for values specified in the `timesteps` argument. 's' means
|
||||
seconds, 'min' means minutes, 'h' means hours, and 'MWd/kg' indicates
|
||||
that the values are given in burnup (MW-d of energy deposited per
|
||||
kilogram of initial heavy metal).
|
||||
n_steps : int, optional
|
||||
Number of stochastic iterations per depletion interval.
|
||||
Must be greater than zero. Default : 10
|
||||
|
|
@ -805,10 +856,10 @@ class SIIntegrator(Integrator):
|
|||
Number of stochastic iterations per depletion interval
|
||||
"""
|
||||
def __init__(self, operator, timesteps, power=None, power_density=None,
|
||||
n_steps=10):
|
||||
timestep_units='s', n_steps=10):
|
||||
check_type("n_steps", n_steps, Integral)
|
||||
check_greater_than("n_steps", n_steps, 0)
|
||||
super().__init__(operator, timesteps, power, power_density)
|
||||
super().__init__(operator, timesteps, power, power_density, timestep_units)
|
||||
self.n_steps = n_steps
|
||||
|
||||
def _get_bos_data_from_operator(self, step_index, step_power, bos_conc):
|
||||
|
|
|
|||
|
|
@ -31,9 +31,11 @@ class PredictorIntegrator(Integrator):
|
|||
----------
|
||||
operator : openmc.deplete.TransportOperator
|
||||
Operator to perform transport simulations
|
||||
timesteps : iterable of float
|
||||
Array of timesteps in units of [s]. Note that values are not
|
||||
cumulative.
|
||||
timesteps : iterable of float or iterable of tuple
|
||||
Array of timesteps. Note that values are not cumulative. The units are
|
||||
specified by the `timestep_units` argument when `timesteps` is an
|
||||
iterable of float. Alternatively, units can be specified for each step
|
||||
by passing an iterable of (value, unit) tuples.
|
||||
power : float or iterable of float, optional
|
||||
Power of the reactor in [W]. A single value indicates that
|
||||
the power is constant over all timesteps. An iterable
|
||||
|
|
@ -46,6 +48,11 @@ class PredictorIntegrator(Integrator):
|
|||
Power density of the reactor in [W/gHM]. It is multiplied by
|
||||
initial heavy metal inventory to get total power if ``power``
|
||||
is not speficied.
|
||||
timestep_units : {'s', 'min', 'h', 'd', 'MWd/kg'}
|
||||
Units for values specified in the `timesteps` argument. 's' means
|
||||
seconds, 'min' means minutes, 'h' means hours, and 'MWd/kg' indicates
|
||||
that the values are given in burnup (MW-d of energy deposited per
|
||||
kilogram of initial heavy metal).
|
||||
|
||||
Attributes
|
||||
----------
|
||||
|
|
@ -113,9 +120,11 @@ class CECMIntegrator(Integrator):
|
|||
----------
|
||||
operator : openmc.deplete.TransportOperator
|
||||
Operator to perform transport simulations
|
||||
timesteps : iterable of float
|
||||
Array of timesteps in units of [s]. Note that values are not
|
||||
cumulative.
|
||||
timesteps : iterable of float or iterable of tuple
|
||||
Array of timesteps. Note that values are not cumulative. The units are
|
||||
specified by the `timestep_units` argument when `timesteps` is an
|
||||
iterable of float. Alternatively, units can be specified for each step
|
||||
by passing an iterable of (value, unit) tuples.
|
||||
power : float or iterable of float, optional
|
||||
Power of the reactor in [W]. A single value indicates that
|
||||
the power is constant over all timesteps. An iterable
|
||||
|
|
@ -128,6 +137,11 @@ class CECMIntegrator(Integrator):
|
|||
Power density of the reactor in [W/gHM]. It is multiplied by
|
||||
initial heavy metal inventory to get total power if ``power``
|
||||
is not speficied.
|
||||
timestep_units : {'s', 'min', 'h', 'd', 'MWd/kg'}
|
||||
Units for values specified in the `timesteps` argument. 's' means
|
||||
seconds, 'min' means minutes, 'h' means hours, and 'MWd/kg' indicates
|
||||
that the values are given in burnup (MW-d of energy deposited per
|
||||
kilogram of initial heavy metal).
|
||||
|
||||
Attributes
|
||||
----------
|
||||
|
|
@ -203,9 +217,11 @@ class CF4Integrator(Integrator):
|
|||
----------
|
||||
operator : openmc.deplete.TransportOperator
|
||||
Operator to perform transport simulations
|
||||
timesteps : iterable of float
|
||||
Array of timesteps in units of [s]. Note that values are not
|
||||
cumulative.
|
||||
timesteps : iterable of float or iterable of tuple
|
||||
Array of timesteps. Note that values are not cumulative. The units are
|
||||
specified by the `timestep_units` argument when `timesteps` is an
|
||||
iterable of float. Alternatively, units can be specified for each step
|
||||
by passing an iterable of (value, unit) tuples.
|
||||
power : float or iterable of float, optional
|
||||
Power of the reactor in [W]. A single value indicates that
|
||||
the power is constant over all timesteps. An iterable
|
||||
|
|
@ -218,6 +234,11 @@ class CF4Integrator(Integrator):
|
|||
Power density of the reactor in [W/gHM]. It is multiplied by
|
||||
initial heavy metal inventory to get total power if ``power``
|
||||
is not speficied.
|
||||
timestep_units : {'s', 'min', 'h', 'd', 'MWd/kg'}
|
||||
Units for values specified in the `timesteps` argument. 's' means
|
||||
seconds, 'min' means minutes, 'h' means hours, and 'MWd/kg' indicates
|
||||
that the values are given in burnup (MW-d of energy deposited per
|
||||
kilogram of initial heavy metal).
|
||||
|
||||
Attributes
|
||||
----------
|
||||
|
|
@ -310,9 +331,11 @@ class CELIIntegrator(Integrator):
|
|||
----------
|
||||
operator : openmc.deplete.TransportOperator
|
||||
Operator to perform transport simulations
|
||||
timesteps : iterable of float
|
||||
Array of timesteps in units of [s]. Note that values are not
|
||||
cumulative.
|
||||
timesteps : iterable of float or iterable of tuple
|
||||
Array of timesteps. Note that values are not cumulative. The units are
|
||||
specified by the `timestep_units` argument when `timesteps` is an
|
||||
iterable of float. Alternatively, units can be specified for each step
|
||||
by passing an iterable of (value, unit) tuples.
|
||||
power : float or iterable of float, optional
|
||||
Power of the reactor in [W]. A single value indicates that
|
||||
the power is constant over all timesteps. An iterable
|
||||
|
|
@ -325,6 +348,11 @@ class CELIIntegrator(Integrator):
|
|||
Power density of the reactor in [W/gHM]. It is multiplied by
|
||||
initial heavy metal inventory to get total power if ``power``
|
||||
is not speficied.
|
||||
timestep_units : {'s', 'min', 'h', 'd', 'MWd/kg'}
|
||||
Units for values specified in the `timesteps` argument. 's' means
|
||||
seconds, 'min' means minutes, 'h' means hours, and 'MWd/kg' indicates
|
||||
that the values are given in burnup (MW-d of energy deposited per
|
||||
kilogram of initial heavy metal).
|
||||
|
||||
Attributes
|
||||
----------
|
||||
|
|
@ -404,9 +432,11 @@ class EPCRK4Integrator(Integrator):
|
|||
----------
|
||||
operator : openmc.deplete.TransportOperator
|
||||
Operator to perform transport simulations
|
||||
timesteps : iterable of float
|
||||
Array of timesteps in units of [s]. Note that values are not
|
||||
cumulative.
|
||||
timesteps : iterable of float or iterable of tuple
|
||||
Array of timesteps. Note that values are not cumulative. The units are
|
||||
specified by the `timestep_units` argument when `timesteps` is an
|
||||
iterable of float. Alternatively, units can be specified for each step
|
||||
by passing an iterable of (value, unit) tuples.
|
||||
power : float or iterable of float, optional
|
||||
Power of the reactor in [W]. A single value indicates that
|
||||
the power is constant over all timesteps. An iterable
|
||||
|
|
@ -419,6 +449,11 @@ class EPCRK4Integrator(Integrator):
|
|||
Power density of the reactor in [W/gHM]. It is multiplied by
|
||||
initial heavy metal inventory to get total power if ``power``
|
||||
is not speficied.
|
||||
timestep_units : {'s', 'min', 'h', 'd', 'MWd/kg'}
|
||||
Units for values specified in the `timesteps` argument. 's' means
|
||||
seconds, 'min' means minutes, 'h' means hours, and 'MWd/kg' indicates
|
||||
that the values are given in burnup (MW-d of energy deposited per
|
||||
kilogram of initial heavy metal).
|
||||
|
||||
Attributes
|
||||
----------
|
||||
|
|
@ -518,9 +553,11 @@ class LEQIIntegrator(Integrator):
|
|||
----------
|
||||
operator : openmc.deplete.TransportOperator
|
||||
Operator to perform transport simulations
|
||||
timesteps : iterable of float
|
||||
Array of timesteps in units of [s]. Note that values are not
|
||||
cumulative.
|
||||
timesteps : iterable of float or iterable of tuple
|
||||
Array of timesteps. Note that values are not cumulative. The units are
|
||||
specified by the `timestep_units` argument when `timesteps` is an
|
||||
iterable of float. Alternatively, units can be specified for each step
|
||||
by passing an iterable of (value, unit) tuples.
|
||||
power : float or iterable of float, optional
|
||||
Power of the reactor in [W]. A single value indicates that
|
||||
the power is constant over all timesteps. An iterable
|
||||
|
|
@ -533,6 +570,11 @@ class LEQIIntegrator(Integrator):
|
|||
Power density of the reactor in [W/gHM]. It is multiplied by
|
||||
initial heavy metal inventory to get total power if ``power``
|
||||
is not speficied.
|
||||
timestep_units : {'s', 'min', 'h', 'd', 'MWd/kg'}
|
||||
Units for values specified in the `timesteps` argument. 's' means
|
||||
seconds, 'min' means minutes, 'h' means hours, and 'MWd/kg' indicates
|
||||
that the values are given in burnup (MW-d of energy deposited per
|
||||
kilogram of initial heavy metal).
|
||||
|
||||
Attributes
|
||||
----------
|
||||
|
|
@ -629,9 +671,11 @@ class SICELIIntegrator(SIIntegrator):
|
|||
----------
|
||||
operator : openmc.deplete.TransportOperator
|
||||
The operator object to simulate on.
|
||||
timesteps : iterable of float
|
||||
Array of timesteps in units of [s]. Note that values are not
|
||||
cumulative.
|
||||
timesteps : iterable of float or iterable of tuple
|
||||
Array of timesteps. Note that values are not cumulative. The units are
|
||||
specified by the `timestep_units` argument when `timesteps` is an
|
||||
iterable of float. Alternatively, units can be specified for each step
|
||||
by passing an iterable of (value, unit) tuples.
|
||||
power : float or iterable of float, optional
|
||||
Power of the reactor in [W]. A single value indicates that
|
||||
the power is constant over all timesteps. An iterable
|
||||
|
|
@ -644,6 +688,11 @@ class SICELIIntegrator(SIIntegrator):
|
|||
Power density of the reactor in [W/gHM]. It is multiplied by
|
||||
initial heavy metal inventory to get total power if ``power``
|
||||
is not speficied.
|
||||
timestep_units : {'s', 'min', 'h', 'd', 'MWd/kg'}
|
||||
Units for values specified in the `timesteps` argument. 's' means
|
||||
seconds, 'min' means minutes, 'h' means hours, and 'MWd/kg' indicates
|
||||
that the values are given in burnup (MW-d of energy deposited per
|
||||
kilogram of initial heavy metal).
|
||||
n_steps : int, optional
|
||||
Number of stochastic iterations per depletion interval.
|
||||
Must be greater than zero. Default : 10
|
||||
|
|
@ -730,9 +779,11 @@ class SILEQIIntegrator(SIIntegrator):
|
|||
----------
|
||||
operator : openmc.deplete.TransportOperator
|
||||
The operator object to simulate on.
|
||||
timesteps : iterable of float
|
||||
Array of timesteps in units of [s]. Note that values are not
|
||||
cumulative.
|
||||
timesteps : iterable of float or iterable of tuple
|
||||
Array of timesteps. Note that values are not cumulative. The units are
|
||||
specified by the `timestep_units` argument when `timesteps` is an
|
||||
iterable of float. Alternatively, units can be specified for each step
|
||||
by passing an iterable of (value, unit) tuples.
|
||||
power : float or iterable of float, optional
|
||||
Power of the reactor in [W]. A single value indicates that
|
||||
the power is constant over all timesteps. An iterable
|
||||
|
|
@ -745,6 +796,11 @@ class SILEQIIntegrator(SIIntegrator):
|
|||
Power density of the reactor in [W/gHM]. It is multiplied by
|
||||
initial heavy metal inventory to get total power if ``power``
|
||||
is not speficied.
|
||||
timestep_units : {'s', 'min', 'h', 'd', 'MWd/kg'}
|
||||
Units for values specified in the `timesteps` argument. 's' means
|
||||
seconds, 'min' means minutes, 'h' means hours, and 'MWd/kg' indicates
|
||||
that the values are given in burnup (MW-d of energy deposited per
|
||||
kilogram of initial heavy metal).
|
||||
n_steps : int, optional
|
||||
Number of stochastic iterations per depletion interval.
|
||||
Must be greater than zero. Default : 10
|
||||
|
|
|
|||
|
|
@ -22,6 +22,8 @@ class _Settings(object):
|
|||
entropy_on = _DLLGlobal(c_bool, 'entropy_on')
|
||||
generations_per_batch = _DLLGlobal(c_int32, 'gen_per_batch')
|
||||
inactive = _DLLGlobal(c_int32, 'n_inactive')
|
||||
max_lost_particles = _DLLGlobal(c_int32, 'max_lost_particles')
|
||||
rel_max_lost_particles = _DLLGlobal(c_double, 'rel_max_lost_particles')
|
||||
particles = _DLLGlobal(c_int64, 'n_particles')
|
||||
restart_run = _DLLGlobal(c_bool, 'restart_run')
|
||||
run_CE = _DLLGlobal(c_bool, 'run_CE')
|
||||
|
|
|
|||
|
|
@ -58,6 +58,10 @@ class Settings(object):
|
|||
history-based parallelism.
|
||||
generations_per_batch : int
|
||||
Number of generations per batch
|
||||
max_lost_particles : int
|
||||
Maximum number of lost particles
|
||||
rel_max_lost_particles : int
|
||||
Maximum number of lost particles, relative to the total number of particles
|
||||
inactive : int
|
||||
Number of inactive batches
|
||||
keff_trigger : dict
|
||||
|
|
@ -176,6 +180,8 @@ class Settings(object):
|
|||
self._batches = None
|
||||
self._generations_per_batch = None
|
||||
self._inactive = None
|
||||
self._max_lost_particles = None
|
||||
self._rel_max_lost_particles = None
|
||||
self._particles = None
|
||||
self._keff_trigger = None
|
||||
|
||||
|
|
@ -254,6 +260,14 @@ class Settings(object):
|
|||
def inactive(self):
|
||||
return self._inactive
|
||||
|
||||
@property
|
||||
def max_lost_particles(self):
|
||||
return self._max_lost_particles
|
||||
|
||||
@property
|
||||
def rel_max_lost_particles(self):
|
||||
return self._rel_max_lost_particles
|
||||
|
||||
@property
|
||||
def particles(self):
|
||||
return self._particles
|
||||
|
|
@ -417,6 +431,19 @@ class Settings(object):
|
|||
cv.check_greater_than('inactive batches', inactive, 0, True)
|
||||
self._inactive = inactive
|
||||
|
||||
@max_lost_particles.setter
|
||||
def max_lost_particles(self, max_lost_particles):
|
||||
cv.check_type('max_lost_particles', max_lost_particles, Integral)
|
||||
cv.check_greater_than('max_lost_particles', max_lost_particles, 0)
|
||||
self._max_lost_particles = max_lost_particles
|
||||
|
||||
@rel_max_lost_particles.setter
|
||||
def rel_max_lost_particles(self, rel_max_lost_particles):
|
||||
cv.check_type('rel_max_lost_particles', rel_max_lost_particles, Real)
|
||||
cv.check_greater_than('rel_max_lost_particles', rel_max_lost_particles, 0)
|
||||
cv.check_less_than('rel_max_lost_particles', rel_max_lost_particles, 1)
|
||||
self._rel_max_lost_particles = rel_max_lost_particles
|
||||
|
||||
@particles.setter
|
||||
def particles(self, particles):
|
||||
cv.check_type('particles', particles, Integral)
|
||||
|
|
@ -763,6 +790,16 @@ class Settings(object):
|
|||
element = ET.SubElement(root, "inactive")
|
||||
element.text = str(self._inactive)
|
||||
|
||||
def _create_max_lost_particles_subelement(self, root):
|
||||
if self._max_lost_particles is not None:
|
||||
element = ET.SubElement(root, "max_lost_particles")
|
||||
element.text = str(self._max_lost_particles)
|
||||
|
||||
def _create_rel_max_lost_particles_subelement(self, root):
|
||||
if self._rel_max_lost_particles is not None:
|
||||
element = ET.SubElement(root, "rel_max_lost_particles")
|
||||
element.text = str(self._rel_max_lost_particles)
|
||||
|
||||
def _create_particles_subelement(self, root):
|
||||
if self._particles is not None:
|
||||
element = ET.SubElement(root, "particles")
|
||||
|
|
@ -1009,6 +1046,8 @@ class Settings(object):
|
|||
self._particles_from_xml_element(elem)
|
||||
self._batches_from_xml_element(elem)
|
||||
self._inactive_from_xml_element(elem)
|
||||
self._max_lost_particles_from_xml_element(elem)
|
||||
self._rel_max_lost_particles_from_xml_element(elem)
|
||||
self._generations_per_batch_from_xml_element(elem)
|
||||
|
||||
def _run_mode_from_xml_element(self, root):
|
||||
|
|
@ -1031,6 +1070,16 @@ class Settings(object):
|
|||
if text is not None:
|
||||
self.inactive = int(text)
|
||||
|
||||
def _max_lost_particles_from_xml_element(self, root):
|
||||
text = get_text(root, 'max_lost_particles')
|
||||
if text is not None:
|
||||
self.max_lost_particles = int(text)
|
||||
|
||||
def _rel_max_lost_particles_from_xml_element(self, root):
|
||||
text = get_text(root, 'rel_max_lost_particles')
|
||||
if text is not None:
|
||||
self.rel_max_lost_particles = float(text)
|
||||
|
||||
def _generations_per_batch_from_xml_element(self, root):
|
||||
text = get_text(root, 'generations_per_batch')
|
||||
if text is not None:
|
||||
|
|
@ -1270,6 +1319,8 @@ class Settings(object):
|
|||
self._create_particles_subelement(root_element)
|
||||
self._create_batches_subelement(root_element)
|
||||
self._create_inactive_subelement(root_element)
|
||||
self._create_max_lost_particles_subelement(root_element)
|
||||
self._create_rel_max_lost_particles_subelement(root_element)
|
||||
self._create_generations_per_batch_subelement(root_element)
|
||||
self._create_keff_trigger_subelement(root_element)
|
||||
self._create_source_subelement(root_element)
|
||||
|
|
@ -1340,6 +1391,8 @@ class Settings(object):
|
|||
settings._particles_from_xml_element(root)
|
||||
settings._batches_from_xml_element(root)
|
||||
settings._inactive_from_xml_element(root)
|
||||
settings._max_lost_particles_from_xml_element(root)
|
||||
settings._rel_max_lost_particles_from_xml_element(root)
|
||||
settings._generations_per_batch_from_xml_element(root)
|
||||
settings._keff_trigger_from_xml_element(root)
|
||||
settings._source_from_xml_element(root)
|
||||
|
|
|
|||
|
|
@ -310,7 +310,7 @@ class StatePoint(object):
|
|||
if self.run_mode == 'eigenvalue':
|
||||
return self._f['n_inactive'][()]
|
||||
else:
|
||||
return None
|
||||
return None
|
||||
|
||||
@property
|
||||
def n_particles(self):
|
||||
|
|
|
|||
|
|
@ -631,8 +631,8 @@ Particle::mark_as_lost(const char* message)
|
|||
|
||||
// Abort the simulation if the maximum number of lost particles has been
|
||||
// reached
|
||||
if (simulation::n_lost_particles >= MAX_LOST_PARTICLES &&
|
||||
simulation::n_lost_particles >= REL_MAX_LOST_PARTICLES*n) {
|
||||
if (simulation::n_lost_particles >= settings::max_lost_particles &&
|
||||
simulation::n_lost_particles >= settings::rel_max_lost_particles*n) {
|
||||
fatal_error("Maximum number of lost particles has been reached.");
|
||||
}
|
||||
}
|
||||
|
|
|
|||
|
|
@ -80,6 +80,8 @@ std::string path_statepoint;
|
|||
|
||||
int32_t n_batches;
|
||||
int32_t n_inactive {0};
|
||||
int32_t max_lost_particles {10};
|
||||
double rel_max_lost_particles {1.0e-6};
|
||||
int32_t gen_per_batch {1};
|
||||
int64_t n_particles {-1};
|
||||
|
||||
|
|
@ -144,6 +146,16 @@ void get_run_parameters(pugi::xml_node node_base)
|
|||
}
|
||||
if (!trigger_on) n_max_batches = n_batches;
|
||||
|
||||
// Get max number of lost particles
|
||||
if (check_for_node(node_base, "max_lost_particles")) {
|
||||
max_lost_particles = std::stoi(get_node_value(node_base, "max_lost_particles"));
|
||||
}
|
||||
|
||||
// Get relative number of lost particles
|
||||
if (check_for_node(node_base, "rel_max_lost_particles")) {
|
||||
rel_max_lost_particles = std::stod(get_node_value(node_base, "rel_max_lost_particles"));
|
||||
}
|
||||
|
||||
// Get number of inactive batches
|
||||
if (run_mode == RunMode::EIGENVALUE) {
|
||||
if (check_for_node(node_base, "inactive")) {
|
||||
|
|
@ -344,14 +356,18 @@ void read_settings_xml()
|
|||
// Read run parameters
|
||||
get_run_parameters(node_mode);
|
||||
|
||||
// Check number of active batches, inactive batches, and particles
|
||||
// Check number of active batches, inactive batches, max lost particles and particles
|
||||
if (n_batches <= n_inactive) {
|
||||
fatal_error("Number of active batches must be greater than zero.");
|
||||
} else if (n_inactive < 0) {
|
||||
fatal_error("Number of inactive batches must be non-negative.");
|
||||
} else if (n_particles <= 0) {
|
||||
fatal_error("Number of particles must be greater than zero.");
|
||||
}
|
||||
} else if (max_lost_particles <= 0) {
|
||||
fatal_error("Number of max lost particles must be greater than zero.");
|
||||
} else if (rel_max_lost_particles <= 0.0 || rel_max_lost_particles >= 1.0) {
|
||||
fatal_error("Relative max lost particles must be between zero and one.");
|
||||
}
|
||||
}
|
||||
|
||||
// Copy random number seed if specified
|
||||
|
|
|
|||
|
|
@ -7,6 +7,7 @@ will be left unimplemented and testing will be done via regression.
|
|||
"""
|
||||
|
||||
import copy
|
||||
from random import uniform
|
||||
from unittest.mock import MagicMock
|
||||
|
||||
import numpy as np
|
||||
|
|
@ -15,11 +16,24 @@ import pytest
|
|||
|
||||
from openmc.deplete import (
|
||||
ReactionRates, Results, ResultsList, comm, OperatorResult,
|
||||
PredictorIntegrator, SICELIIntegrator)
|
||||
PredictorIntegrator, CECMIntegrator, CF4Integrator, CELIIntegrator,
|
||||
EPCRK4Integrator, LEQIIntegrator, SICELIIntegrator, SILEQIIntegrator)
|
||||
|
||||
from tests import dummy_operator
|
||||
|
||||
|
||||
INTEGRATORS = [
|
||||
PredictorIntegrator,
|
||||
CECMIntegrator,
|
||||
CF4Integrator,
|
||||
CELIIntegrator,
|
||||
EPCRK4Integrator,
|
||||
LEQIIntegrator,
|
||||
SICELIIntegrator,
|
||||
SILEQIIntegrator
|
||||
]
|
||||
|
||||
|
||||
def test_results_save(run_in_tmpdir):
|
||||
"""Test data save module"""
|
||||
|
||||
|
|
@ -105,14 +119,14 @@ def test_results_save(run_in_tmpdir):
|
|||
np.testing.assert_array_equal(res[1].time, t2)
|
||||
|
||||
|
||||
@pytest.mark.parametrize("timesteps", (1, [1]))
|
||||
def test_bad_integrator_inputs(timesteps):
|
||||
def test_bad_integrator_inputs():
|
||||
"""Test failure modes for Integrator inputs"""
|
||||
|
||||
op = MagicMock()
|
||||
op.prev_res = None
|
||||
op.chain = None
|
||||
op.heavy_metal = 1.0
|
||||
timesteps = [1]
|
||||
|
||||
# No power nor power density given
|
||||
with pytest.raises(ValueError, match="Either power or power density"):
|
||||
|
|
@ -159,3 +173,62 @@ def test_integrator(run_in_tmpdir, scheme):
|
|||
dep_time = res.get_depletion_time()
|
||||
assert dep_time.shape == (2, )
|
||||
assert all(dep_time > 0)
|
||||
|
||||
|
||||
@pytest.mark.parametrize("integrator", INTEGRATORS)
|
||||
def test_timesteps(integrator):
|
||||
# Crate fake operator
|
||||
op = MagicMock()
|
||||
op.prev_res = None
|
||||
op.chain = None
|
||||
|
||||
# Set heavy metal mass and power randomly
|
||||
op.heavy_metal = uniform(0, 10000)
|
||||
power = uniform(0, 1e6)
|
||||
|
||||
# Reference timesteps in seconds
|
||||
day = 86400.0
|
||||
ref_timesteps = [1*day, 2*day, 5*day, 10*day]
|
||||
|
||||
# Case 1, timesteps in seconds
|
||||
timesteps = ref_timesteps
|
||||
x = integrator(op, timesteps, power, timestep_units='s')
|
||||
assert np.allclose(x.timesteps, ref_timesteps)
|
||||
|
||||
# Case 2, timesteps in minutes
|
||||
minute = 60
|
||||
timesteps = [t / minute for t in ref_timesteps]
|
||||
x = integrator(op, timesteps, power, timestep_units='min')
|
||||
assert np.allclose(x.timesteps, ref_timesteps)
|
||||
|
||||
# Case 3, timesteps in hours
|
||||
hour = 60*60
|
||||
timesteps = [t / hour for t in ref_timesteps]
|
||||
x = integrator(op, timesteps, power, timestep_units='h')
|
||||
assert np.allclose(x.timesteps, ref_timesteps)
|
||||
|
||||
# Case 4, timesteps in days
|
||||
timesteps = [t / day for t in ref_timesteps]
|
||||
x = integrator(op, timesteps, power, timestep_units='d')
|
||||
assert np.allclose(x.timesteps, ref_timesteps)
|
||||
|
||||
# Case 5, timesteps in MWd/kg
|
||||
kilograms = op.heavy_metal / 1000.0
|
||||
days = [t/day for t in ref_timesteps]
|
||||
megawatts = power / 1000000.0
|
||||
burnup = [t * megawatts / kilograms for t in days]
|
||||
x = integrator(op, burnup, power, timestep_units='MWd/kg')
|
||||
assert np.allclose(x.timesteps, ref_timesteps)
|
||||
|
||||
# Case 6, mixed units
|
||||
burnup_per_day = (1e-6*power) / kilograms
|
||||
timesteps = [(burnup_per_day, 'MWd/kg'), (2*day, 's'), (5, 'd'),
|
||||
(10*burnup_per_day, 'MWd/kg')]
|
||||
x = integrator(op, timesteps, power)
|
||||
assert np.allclose(x.timesteps, ref_timesteps)
|
||||
|
||||
# Bad units should raise an exception
|
||||
with pytest.raises(ValueError, match="unit"):
|
||||
integrator(op, ref_timesteps, power, timestep_units='🐨')
|
||||
with pytest.raises(ValueError, match="unit"):
|
||||
integrator(op, [(800.0, 'gorillas')], power)
|
||||
|
|
|
|||
|
|
@ -9,6 +9,8 @@ def test_export_to_xml(run_in_tmpdir):
|
|||
s.generations_per_batch = 10
|
||||
s.inactive = 100
|
||||
s.particles = 1000000
|
||||
s.max_lost_particles = 5
|
||||
s.rel_max_lost_particles = 1e-4
|
||||
s.keff_trigger = {'type': 'std_dev', 'threshold': 0.001}
|
||||
s.energy_mode = 'continuous-energy'
|
||||
s.max_order = 5
|
||||
|
|
@ -62,6 +64,8 @@ def test_export_to_xml(run_in_tmpdir):
|
|||
assert s.generations_per_batch == 10
|
||||
assert s.inactive == 100
|
||||
assert s.particles == 1000000
|
||||
assert s.max_lost_particles == 5
|
||||
assert s.rel_max_lost_particles == 1e-4
|
||||
assert s.keff_trigger == {'type': 'std_dev', 'threshold': 0.001}
|
||||
assert s.energy_mode == 'continuous-energy'
|
||||
assert s.max_order == 5
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue