For the case where there are more processes than burnable materials,
e.g. single pin cell or assembly with one burnable material, then
the local_indexes attribute on some FissionYieldHelpers will
be an empty array. This causes errors trying to use this as an
indexer into the tally results.
Instead, simply return from the unpacking if this is the case.
The compute_yields function is called inside the Operator's
unpacking of it's local materials, so that __should__ not be
in a position to call compute_yields for a process with no
burnable materials.
Create an empty matrix to fit energy-dependent fission
rates for [materials, energies, nuclides].
Use the nonzero method on total fission rate to find
indicies along the first and last axis for non-zero total
fission rate. Populate corresponding fractional fission
rates by dividing group fission rates by total fission rate,
using the indices for materials and nuclides with non-zero
total fission rate.
Use numpy.where to find where total fission rate is zero,
and directly set these locations to zero in the final
result matrix.
Also clean up some lint in openmc.deplete.nuclide related to
old variable names and unused imports.
Operator now has a new helper for tallying fission rates
and computing fission yields using energy dependency.
The FissionYieldHelper creates fission rate tallies across
burnable materials with an energy filter. The energy filter
corresponds to each energy group for which there are fission
yields.
During the operator's unpacking of tallies, this helper
computes a relative fission contribution by normalizing
fission rates, such that the sum of all group fission rates
in each material and each nuclide is unity.
A single fission yield for each parent nuclide is computed
by weighting the energy-dependent fission yields by
these relative contributions.
A nested dictionary {parent: {product: yield}} is added to
a library, and then attached to the depletion chain
at the end of the unpacking. These libraries are used
during the depletion process instead of the only-thermal
yields from before.
A test was added that creates a proxy helper that doesn't
use the C-API, but unpacks the "fission tally" and
computes the fission yield library in an identical manner
as the actual class. Since the proxy helper is a subclass
of the one used in real transport, the functionality
is consistent.
Add two new classes for working with fission yield data
on the Chain. The primary class is FissionYieldDistribution.
This class stores distributions for one nuclide with
potentially many energies, with the assumption that the
products don't change __too__ much across energy.
Looking at the CASL chain and one generated by ENDF data,
this appears to be the case. Most distributions are
"full" in the sense that energies produce the same products.
There are some cases where one or two products may not
exist for a given energy.
The FissionYieldDistribution retains a dictionary-like behavior,
e.g. d[0.0253]["Xe135"] is a valid command and returns the
yield for Xe-135 at a "thermal" spectrum, due to yields provided
at 0.0253 eV. This is done with a helper class, _FissionYield,
implemented first to support this behavior, but also to
allow vector-operations in combining fission yields.
These _FissionYield objects store the same product vector and
a view into the underlying yield_matrix for a single energy.
To support simple weighting of yields, the __mull__ and
__iadd__ methods are implemented. A copy method is provided
to remove the chance of modifying the original yield data.
test_deplete_chain and test_deplete_nuclide have been modified
in order to utilize these classes, without ruining the validity
of the tests. Tests for the view / copy methods and vector
operations on _FissionYield objects are included.
Takes a single set of fission yields and passes them as
an additional argument to matrix_func:
>>> A = matrix_func(chain, rates, fission_yields)
Applied to cf4, epc_rk4, celi, and leqi functions.
Assumes that fission yields will not change during a
depletion event. This change is probably overshadowed
by how much the reaction rates may change, but still worth
pointing out.
In support of using energy dependent fission yields, this
commit allows a dictionary of fission yields to be passed
into Chain.form_matrix. The dictionary is expected to
be of the form
{source_name: {target_name: yield_fraction}}
where source_name and target_name are string GND names
for fissionable nuclide and fission yield products, respectively.
Currently, the fission yield dictionary does not have to be passed,
defaulting to using the thermal yield values. This is done to
make testing easier, and for back compatability. The goal of this
feature, however, is to pass material [and thus spectrum] specific
yields into this method.
The test test_deplete_chain::test_form_matrix has been modified
to pass the exact fission yield dictionary into Chain.form_matrix.
The resulting matrix is compared against the matrix when
no yields are passed.
The LEQI-type integrators are the only ones that need the current
depletion index. All other integrators now display the same
_i=None for this index in the call signature, and denote the
variable as optional and unused in the docstring
Much of the previous API is intact, with the major change
being CRAM functions are imported from openmc.deplete.cram
in the test_deplete_cram.
The integrator abstract classes are placed in
openmc/deplete/abc.py with all concrete classes going in to
openmc/deplete/integrators.py
Documentation updated accordingly
Called with a single integer corresponding to the current
depletion step. This method is intended to document the current
beginning-of-step calculations prior to any depletion at each step
and at the final simulation.
The Operator uses this method to write a statepoint using the
C API.
Also rename the Integrator _get_bos_data_from_openmc to
_get_bos_data_from_operator for generality. It is here that the
write_bos_data method is called.
Removed the private method for saving to avoid unnecessary
redirection. Also renamed some private attributes:
1) Integrator._ires -> Integrator._i_res
2) Integrator._N_STAGES -> Integrator._num_stages
Closes#1275 by passing reaction rates and numbers from
a previous result according to what materials are tracked
on an operator. A new method, Results.distribute, creates a new
Results object by
- directly copying over "global" data, like time, keff, and
maps describing where each material lives in the depletion_results
hdf5 file
- mapping volumes for local materials
- Extracting slices of numbers and reaction rates corresponding
to local material
This allows the Operator to create a unique ResultsList instance
containing reaction rates and compositions pertaining to
the materials tracked on this process (local).
This has been tested by comparing depletion_result files from restarts
1) using a serial restart run and 2) using two MPI process.
A quarter PWR assembly with 71 burnable materials was used.
Comparing with h5diff -r -p 0.01 <old> <new> <group> where
group=["eigenvalues", "number", "reaction rates"] and
-p 0.01 reports errors over 0.01% revealed no differences.
The new class method ResultsList.from_hdf5 should be
used to load in results data, rather than passing
the file name into the __init__ method. __init__
passes directly to list.__init__ now.
The motivation for this is to resolve the depletion
restart with MPI bug #1275. To do this, each Operator
will create it's own ResultsList and distribute reaction
rates and densities according to what materials are burned
on this process. Rather than use a normal list, this Operator
expects the various accessor methods like get_atoms to be
present on self.prev_res
- SI_Integrator must be passed positive integer.
- Number of powers computed must be equal to the
number of time steps
Added unit test to ensure that the right errors are raised.
Closes#1283 by writing a statepoint file at the beginning
of every depletion step [BOS] and at the final transport simulation,
[EOL]. These files are named "openmc_simulation_n<N>.h5", where
<N> represents the current depletion step, including restart
steps. The new files are written using the C API.
It is worth noting that OpenMC will still write statepoint files
after each intermediate transport simulation used by various
depletion schemes. This commit makes no effort to move or delete
these files, nor stop their creation. It is not clear how,
using the python nor C API, one can stop OpenMC from
writing a statepoint file.
The source bank is not written to these files as this causes
a segmentation fault.
Comparisons of k_combined pull from all statepoints,
reference depletion_results file, and test depletion_results
file have been added to depletion regression test.
Each scheme is capable of different amounts of concentration
and reaction data, depending on the number of intermediate
transport solutions of interest. The PredictorIntegrator
stores data from a single solution, while the CF4Integrator and
EPC_RK4_Integrator store four stages. This discrepancy can cause
issues when saving data in the depletion output file, as
either some data is calculated and thrown away [more stages going
to fewer stages] or excess space is not utilized [fewer stages
going in to a restart that previously used more stages].
Each subclass of Integrator is expected to declare a
_N_STAGES attribute that provides information on the number of
stages expected. Iterative schemes like SI_CELI_Integrator
only report two stages, for the initial concentration and
first prediction, rather than at each iteration point.
Logic that checks for consistent sizes is included in the
initialization of the Integrators, and has been removed
from Results.save
Add three functions that can be used to modify EnergyFunctionFilters
through the C-API:
- openmc_energyfunc_filter_set_data: set energy and y data
- openmc_energyfunc_filter_get_energy: obtain energies used in
interpolation
- openmc_energyfunc_filter_get_y: obtain ordinate values
These functions are modeled after openmc_energy_filter_[get|set]_bins.
The set_data function relies upon the new
EnergyFunctionFilter::set_data function, which is analogous
to EnergyFilter::set_bins function. Checks are performed
to make sure the energy and ordinate vectors are of equal
size before resetting and populating energy and y
private members.
An EnergyFunctionFilter did exist in openmc.capi, and has
now been flushed out to provide a better __init__ method,
as well as properties for retrieving energies and ordinates
for interpolation.
Changes in d32139b993 highlighted
the repeated data that was saved during non-predictor restart runs.
All Integrator subclasses now have the same
_get_bos_data_from_openmc and _get_bos_data_from_restart methods,
where the restart index ``self._ires`` is one less than the number of
previous results stored on the Operator. This allows all schemes to
write their data in the correct index during restart calculations.
A side benefit of this is that now, the Predictor implementation
can rely entirely on the updated Integrator, without the
need for re-writing the integrate method nor __init__ method.
Many repetitions of the following change:
- assert y1[3] == approx(s2[0])
- assert y2[3] == approx(s2[1])
+ assert y1[2] == approx(s2[0])
+ assert y2[2] == approx(s2[1])
This causes many of the restart tests to fail, because all
schemes but the predictor save repeated data under restart conditions.
Performing one normal depletion event, non-restarted, and then
one restart event, as these tests do, should produce vectors for
time and densities like:
t = [t0, t1, t2]
a = [a0, a1, a2]
| ^ after restart
^ non-restart run
With this, the current tests would fail due to an IndexError, as
the atom vectors should only have three elements. Instead,
non-predictor schemes return vectors
t = [t0, t1, t1, t2]
a = [a0, a1, a1, a2]
and one can access a[3] to obtain the desired EOS concentrations.
The depletion function openmc.deplete.si_leqi has been
removed in favor of the SI_LEQI_Integrator class. The same
depletion scheme can be obtained with the following commands:
>>> leqi = openmc.deplete.SI_LEQI_Integrator(operator, dt, power)
>>> leqi.integrate()
The expression can be onlined for compactness.
The si_celi_inner function has been removed completely now,
as the SI_CELI iteration is performed by directly calling
SI_CELI_Integrator.__call__ through the SI_LEQI_Integrator.
This is similar to how the LEQIIntegrator handles the initial steps.
Tests have been updated to use this class, and the class has been
added to the documentation. No pure-function integration
schemes exist anymore.
Operator.initial_condition return list of vectors representing
the initial composition for a burnable material as an element
in the list. This format is expected from depletion solvers as well.
Previously, when pulling from the previous results on an operator,
the initial concentration was returned simply as a 2D array of
compositions. This commit converts the 2D array to a list of
vectors to be consistent across the integration schemes.
The si_celi depletion function has been removed. The
SI-CELI method can be implemented with
>>> from openmc.deplete import SI_CELI_Integrator
>>> SI_CELI_Integrator(op, dt, power, stages).integrate()
The stages parameter is optional, and defaults to 10 to
be consistent with the previous default for si_celi.
The SI_CELI_Integrator inherits from the new abstract base class
SI_Integrator. There are a few differences in how the SI-based
methods perform the integration.
- The initial, t=0.0, i=0, no restart transport solution is scaled
by the number of stages.
- The SI methods also do not perform a new transport solution at
each successive iteration [t>0, i>0]. Instead, the reaction
rates are pulled from the last stage of the previous step.
- There is no final transport solution once all the depletion
steps have been taken. The final reaction rates are saved as
those from the last stage of the last depletion step.
The si_celi function has been removed from documentation and testing,
and replaced with the SI_CELI_Integrator.
This allows the propagation of the uncertainty through
the SIE iterations, and alows the expressions for updating
k to be preserved. The OperatorResult.k was previous stored
as a tuple of two floats. The interface into the depletion_results
file is not changed, as the Results.save method, which writes the
OperatorResult data, breaks the ufloat into k and uncertainties.
Documentation has been updated on the OperatorResult.
The function openmc.deplete.epc_rk4 has been removed
in favor of openmc.deplete.EPC_RK4_Integrator. The scheme
can be used with::
>>> from openmc.deplete import EPC_RK4_Integrator
>>> EPC_RK4_Integrator(op, dt, power).integrate()
epc_rk4 has been removed from the documentation, and the
EPC_RK4_Integrator class has been added to the depletion
API documentation