Commit graph

9978 commits

Author SHA1 Message Date
Andrew Johnson
c732336536
Avert bad numpy slices unpacking fission yields with no materials
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.
2019-08-08 08:14:32 -05:00
Andrew Johnson
27f2c6b15f
Guard against divide by zero in fission yield helper
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.
2019-08-08 08:14:27 -05:00
Andrew Johnson
599d473888
Add FissionYieldHelper for energy-dependent fission yields
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.
2019-08-08 08:13:06 -05:00
Andrew Johnson
8cd71be864
Store depletion fission yields with common yield matrix
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.
2019-08-07 19:54:49 -05:00
Andrew Johnson
72bf884d9b
Pass fission yields to depletion matrix_func
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.
2019-08-07 19:54:49 -05:00
Andrew Johnson
f3e7c7c623
Support for passing fission yields to Chain.form_matrix
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.
2019-08-07 19:54:45 -05:00
Paul Romano
f5c421d24e
Merge pull request #1312 from drewejohnson/bug-restart-mpi
Distribute reaction rates and densities when restarting under MPI
2019-08-07 16:46:31 -05:00
Andrew Johnson
becda0c2c2
Merge branch 'develop' into bug-restart-mpi 2019-08-07 13:50:42 -05:00
Paul Romano
21fd9a8a13
Merge pull request #1287 from drewejohnson/integrator-class
Class-based time integration for depletion
2019-08-07 13:18:14 -05:00
Andrew Johnson
144aa3a2e3
Document purpose for ranges in Results.distribute 2019-08-07 11:32:20 -05:00
Andrew Johnson
4043e02fc7
Apply suggestions from code review
Co-Authored-By: Paul Romano <paul.k.romano@gmail.com>
2019-08-07 11:20:34 -05:00
Andrew Johnson
2514918591
Consistently document unused depletion index
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
2019-08-07 11:09:51 -05:00
Andrew Johnson
0bc7800a92
Move openmc.deplete.integrator into openmc.deplete
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
2019-08-07 11:04:37 -05:00
Andrew Johnson
5adc3b5816
Document special and abstract methods for Integrators 2019-08-07 09:13:14 -05:00
Andrew Johnson
18c4d5f432
Document init parameters and attributes for Integrators 2019-08-07 09:01:43 -05:00
Andrew Johnson
1bd663daf6
Rename integrator classes to be PEP8 complaint
Changes:
- EPC_RK4_Integrator -> EPCRK4Integrator
- SI_Integrator -> SIIntegrator
- SI_CELI_Integrator -> SICELIIntegrator
- SI_LEQI_Integrator -> SICELIIntegrator
2019-08-07 08:55:40 -05:00
Andrew Johnson
a4905c03af
Add TransportOperator.write_bos_data abstract method
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.
2019-08-06 17:36:22 -05:00
Andrew Johnson
e286940f8e
Pass prev_res onto TransportOperator.__init__
Document that this attribute will either be a
ResultsList or None, both on TransportOperator and
Operator.
2019-08-06 17:14:44 -05:00
Andrew Johnson
78943d1086
Add power as argument to TransportOperator.__call__ 2019-08-06 16:57:21 -05:00
Andrew Johnson
49ba67d480
Use Results.save directly during Integrator.integrate
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
2019-08-06 16:53:09 -05:00
Andrew Johnson
a7e1bca6da
Apply suggestions from code review
Co-Authored-By: Paul Romano <paul.k.romano@gmail.com>
2019-08-06 16:26:23 -05:00
Andrew Johnson
00095b8b74
Distribute process-specific previous results to operator
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.
2019-08-05 16:21:51 -05:00
Andrew Johnson
0fdef55d42
Create ResultsList from file using from_hdf5
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
2019-08-05 16:09:49 -05:00
Andrew Johnson
f633e84822
Merge pull request #1311 from zxkjack123/fix_docs_userguide_typo
fix a typo
2019-08-05 15:51:04 -05:00
zxkjack123
f0727c654a fix a typo 2019-08-05 11:22:16 -05:00
Andrew Johnson
b9683dbdba
Merge branch 'develop' into integrator-class
Need new travis.yml to pass CI
2019-07-30 15:58:49 -05:00
Andrew Johnson
b1ea5b8942
Add type and value checking to Integrators
- 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.
2019-07-30 15:56:24 -05:00
Andrew Johnson
2ee14a4a0f
Write new step-specific statepoint files at BOS, EOL
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.
2019-07-30 15:23:30 -05:00
Paul Romano
d6cbacf4b3
Merge pull request #1303 from drewejohnson/capi-energy-function-filter
Expose/improve EnergyFunctionFilter through C-API
2019-07-30 15:22:05 -05:00
Andrew Johnson
a02feed171
Restrict depletion restart to schemes with common stages
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
2019-07-30 14:19:56 -05:00
Andrew Johnson
51add8f6e9
Clean up changes introduced by Integrators
Remove some unused imports, redeclared __init__,
long lines, etc.
2019-07-30 13:21:04 -05:00
Andrew Johnson
8d920b792a
set_interp_data -> set_data for capi.EnergyFunctionFilter
Add doctring to set_data method. Update test_capi.py with
this change
2019-07-30 11:49:31 -05:00
Andrew Johnson
9d952c3c58
Use a.all() for equality in capi EnergyFunctionFilter tests 2019-07-30 11:49:30 -05:00
Andrew Johnson
9b7ab0bacb
Expose/improve EnergyFunctionFilter through C-API
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.
2019-07-30 11:49:21 -05:00
Paul Romano
d7bcf60eb4 Store min/max data temperatures. Use in Cell::set_temperature to check bounds 2019-07-30 10:18:17 -05:00
Andrew Johnson
c4567b6d9d
Remove print_out options to Operator.__call__, deplete
Time spent running openmc is already presented in after each
transport run. Time spent performing depletion is written in
the depletion file.
2019-07-30 09:49:54 -05:00
Paul Romano
38a5f4a24f
Merge pull request #1302 from drewejohnson/reaction-rate-bug
Fix bug when reaction rate nuclides change through depletion
2019-07-30 09:09:33 -05:00
Andrew Johnson
afe1b5b68c
Fix restart integration offset bugs
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.
2019-07-30 08:46:23 -05:00
Andrew Johnson
52f42ce76c
Add NUMPY_EXPERIMENTAL_ARRAY_FUNCTION to .travis.yml
Instruct travis to use pre-numpy 1.17 behavior.
Workaround issue for uncertainties package causing CI
builds to fail

Related: https://github.com/lebigot/uncertainties/issues/89
2019-07-30 07:50:50 -05:00
Andrew Johnson
62865b3a55
Documentation fixes for pin func, depletion helpers 2019-07-29 15:57:29 -05:00
Andrew Johnson
d5fa7b5c00
Set correct index for final nuclides in restart tests
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.
2019-07-29 14:46:02 -05:00
Andrew Johnson
1f3e7fdac1
Remove si_leqi in favor of SI_LEQI_Integrator
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.
2019-07-29 10:34:59 -05:00
Andrew Johnson
791f7c78b7
Return BOS concentration as list of vectors from restarts
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.
2019-07-29 08:36:09 -05:00
Andrew Johnson
d28546bcd6
Add SI_Integrator, SI_CELI_Integrator classes; remove si_celi
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.
2019-07-26 16:58:51 -05:00
Patrick Shriwise
232a29a2ab
Merge pull request #1301 from drewejohnson/feat-op-dilute
Allow user control over initial concentration for trace isotopes
2019-07-26 09:39:53 -08:00
Andrew Johnson
98957fb473
Store OperatorResults.k as uncertainties.ufloat
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.
2019-07-26 10:31:24 -05:00
Andrew Johnson
58117628b8
Fully remove openmc.deplete.cf4 function from API 2019-07-26 09:53:14 -05:00
Andrew Johnson
e30e6e3f99
Use CF4Integrator in restart tests 2019-07-26 09:52:55 -05:00
Andrew Johnson
39b4d8a1fe
Remove epc_rk4 for EPC_RK4_Integrator class
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
2019-07-26 09:49:58 -05:00
Patrick Shriwise
994221705b
Merge pull request #1297 from paulromano/update-jupyter-notebooks
Update many of the example notebooks
2019-07-25 10:00:43 -05:00