Resonace Covariance Module #1024

icmeyer · 2018-07-17T21:11:35Z

This PR incorporates into the API the ability to parse resonance covariance data from ENDF File 32. Its intended functionality can be seen in the new example notebook: openmc/examples/jupyter/nuclear-data-resonance-covariance.ipynb
Namely:

Parse File 32 data
Allow for selection of a subset of the parameters and corresponding covariance data
Sampling of resonance parameters (currently uses numpy.random.multivariate_normal)
Reconstruction of the cross section using sampled resonance parameters

My first PR so comments on style encouraged!

…ples

…as previous range

icmeyer · 2018-07-20T15:18:05Z

Addressed the comments by @paulromano. Added in warning on the sampling method that it can produce nonphysical negative parameters which may cause reconstruction to fail. This is because the covariance matrix assumes a multivariate normal distribution from which negative numbers may be drawn. Further discussion and possible solutions to this issue can be found in a Rochman paper here.

paulromano

Here are a few more comments on the PR. I would also recommend running pycodestyle on your resonance_covariance.py to address style issues.

paulromano · 2018-07-23T12:00:28Z

openmc/data/endf.py

@@ -71,6 +71,24 @@ def float_endf(s):
    return float(_ENDF_FLOAT_RE.sub(r'\1e\2', s))


+def int_endf(s):
+    """Conver string to int. Used for INTG records where blank entries


Typo: convert

paulromano · 2018-07-23T12:01:07Z

openmc/data/endf.py

@@ -71,6 +71,24 @@ def float_endf(s):
    return float(_ENDF_FLOAT_RE.sub(r'\1e\2', s))


+def int_endf(s):


Change to _int_endf since it is unlikely a user would ever need this

paulromano · 2018-07-23T12:01:37Z

openmc/data/endf.py

@@ -249,6 +267,50 @@ def get_tab2_record(file_obj):

    return params, Tabulated2D(breakpoints, interpolation)

+def get_intg_record(file_obj):


Add another blank line before this function (PEP8: 2 lines between functions/classes at top level)

paulromano · 2018-07-23T12:02:18Z

openmc/data/endf.py

+    corr = corr + corr.T - np.diag(corr.diagonal())
+    return corr
+
+


delete one blank line here

paulromano · 2018-07-23T12:15:36Z

openmc/data/resonance_covariance.py

+        self.cov_subset = cov_subset
+
+    def sample_resonance_parameters(self, n_samples, resonances, use_subset=False):
+        """Return a list  size 'n_samples' of openmc.data.ResonanceRange objects.


I would say something more like "Sample resonance parameters based on the covariances." The type of the object returned (and its length) are described in the "Returns" section so probably not necessary to mention them here too.

paulromano · 2018-07-23T12:29:06Z

openmc/data/resonance_covariance.py

+                    res_range.parameters = sample_params
+                    samples.append(res_range)
+
+        self.samples = samples


As we discussed before, the samples object should be returned directly rather than being set as an attribute of the current instance.

paulromano · 2018-07-23T12:30:22Z

openmc/data/resonance_covariance.py

+        cov_subset[tri_indices] = cov_subset_vals
+
+        self.parameters_subset = parameters_subset
+        self.cov_subset = cov_subset


As with the sampling method, this method should return an object rather than setting attributes of the current instance. I think what would make the most sense is to return a new instance of ResonanceCovarianceRange with the specified subset of parameters/covariance matrix. This would allow you to get rid of the use_subset argument on the sampling method since calling sample_resonance_parameters on the new instance would use the subset by default.

paulromano · 2018-07-23T12:34:29Z

examples/jupyter/nuclear-data-resonance-covariance.ipynb

+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "The covariance module also has the ability to sample a new set of parameters using the covariance matrix. Currently the sampling uses np.multivariate_normal(). Because parameters are assumed to have a multivariate normal distribution this method doesn't not currently guarantee that sampled parameters will be positive."


np.multivariate_normal -> numpy.multivariate_normal

paulromano · 2018-07-23T12:41:15Z

openmc/data/resonance_covariance.py

+        items : list
+            Items from the CONT record at the start of the resonance range
+            subsection
+        file2params : openmc.data.ResonanceRange object


It looks like this is still actually taking a DataFrame of resonance parameters rather than the ResonanceRange object that owns it. I think it would be better to use the ResonanceRange object directly and store it as an attribute of the ResonanceCovarianceRange instance. That way, when sample_resonance_parameters is called, there would be no need to pass the parameters as they could be retrieved from the attribute.

Same comment applies to ReichMooreCovariance

paulromano · 2018-07-23T12:42:49Z

examples/jupyter/nuclear-data-resonance-covariance.ipynb

+   "cell_type": "markdown",
+   "metadata": {},
+   "source": [
+    "The sampling routine requires the incorpotation of the `openmc.data.ResonanceRange` for the same resonance range object. This allows each sample itself to be its own `openmc.data.ResonanceRange` with a new set of parameters. Looking at some of the sampled parameters below:"


Typo: incorpotation

…file2 data as an attribute

paulromano · 2018-07-25T12:08:38Z

openmc/data/resonance_covariance.py

+        parameters = self.parameters
+        cov = self.covariance
+        # Copy ResonanceRange object
+        res_range = copy.copy(self.file2res)


Right now, you're creating a single copy of the resonance range. If multiple samples are requested, this method actually returns a list with the same object repeated N times (this is evident in the Jupyter notebook where you can see that two different samples have the same parameters). You should be creating N copies of the original object and then setting the parameters appropriately for each one.

Saw similar code being used for each formalism/mpar so thought I could pull it out front but guess copy.copy() and .append() behaved differently than I thought. Should be better now.

paulromano · 2018-07-25T12:09:52Z

openmc/data/resonance_covariance.py

+                               'captureWidth', 'fissionWidth', 'competitiveWidth']
+                    sample_params = pd.DataFrame.from_records(records,
+                                                              columns=columns)
+                    res_range.parameters = sample_params


Related to the previous comment, this line overwrites the parameters on each iteration of the for loop. (also, I don't see a similar line updating the parameters for the mpar=4/5 and the RM cases)

paulromano

Ok, hopefully this is the last round!

paulromano · 2018-07-30T19:29:09Z

openmc/data/resonance_covariance.py

@@ -0,0 +1,743 @@
+from collections import defaultdict, MutableSequence, Iterable


Remove defaultdict, Iterable (not used)

paulromano · 2018-07-30T19:29:44Z

openmc/data/resonance_covariance.py

+    mpar : int
+        Number of parameters in covariance matrix for each individual resonance
+    formalism : str
+        String descriptor of formalism


Need to add description of file2res here (since it is assumed to be present for any subclass)

paulromano · 2018-07-30T19:30:48Z

openmc/data/resonance_covariance.py

+        self.energy_min = energy_min
+        self.energy_max = energy_max
+
+    def res_subset(self, parameter_str, bounds):


What do you think about just naming this subset? The res seems kind of superfluous to me.

paulromano · 2018-07-30T19:31:40Z

openmc/data/resonance_covariance.py

+                gn = sample[1::mpar]
+                gg = sample[2::mpar]
+                gfa = sample[3::mpar] if mpar > 3 else parameters['fissionWidthA'].values
+                gfb = sample[3::mpar] if mpar > 3 else parameters['fissionWidthB'].values


paulromano · 2018-07-30T19:37:31Z

openmc/data/resonance_covariance.py

+
+            # Add parameters from File 2
+            parameters = _add_file2_contributions(parameters,
+                                                  resonance.parameters)


It looks like the last few lines here (starting with columns = ...) are the same for each branch of the if -- I think those could be moved out of the if altogether (so that they're not repeated). The only thing that's different is in the lcomp == 1 case, mpar is not set, but that looks like a bug since it is used in the call below.

paulromano · 2018-07-30T19:38:23Z

openmc/data/resonance_covariance.py

+
+            # Add parameters from File 2
+            parameters = _add_file2_contributions(parameters,
+                                                  resonance.parameters)


Same comment here -- the lines starting from columns = ... looks redundant on the branches of the if.

paulromano · 2018-07-30T19:38:39Z

tests/unit_tests/test_data_neutron.py

@@ -39,7 +39,7 @@ def sm150():
 def gd154():
    """Gd154 ENDF data (contains Reich Moore resonance range)"""
    filename = os.path.join(_ENDF_DATA, 'neutrons', 'n-064_Gd_154.endf')
-    return openmc.data.IncidentNeutron.from_endf(filename)
+    return openmc.data.IncidentNeutron.from_endf(filename, covariance = True)


covariance=True

paulromano · 2018-07-30T19:41:39Z

openmc/data/resonance_covariance.py

+        res_cov_range.covariance = cov_subset
+        return res_cov_range
+
+    def sample_resonance_parameters(self, n_samples):


What do you think about naming this simply sample()?

Makes sense to me. Wasn't sure how descriptive the names needed to be but as with res_subset I think that having them as methods of the objects provides enough context

…elop

icmeyer · 2018-08-01T20:57:07Z

It looks like the build is failing with Python 3.4 on the following test.
https://github.com/mit-crpg/openmc/blob/ac57dd929841026bb91cd1f89d67b5556b4ad157/tests/unit_tests/test_capi.py#L89-L100
Running on my machine with 3.4 gives me
> assert m.densities == pytest.approx(test_dens)
E ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()
Does this test need to be changed? Also odd that on some of the builds Travis isn't showing the final report of pytest

icmeyer added 30 commits December 15, 2017 13:03

limited functionality for RM covariance only

f741c21

Merge remote-tracking branch 'upstream/develop' into develop

468f48c

Limited functionality for MLBW/SLBW covariances

8bd02e7

Merge remote-tracking branch 'upstream/develop' into develop

03c7f8b

Added some capability for LCOMP=2 format

b6082ac

Added capability for LCOMP=0

27b6abe

Merge remote-tracking branch 'upstream/develop' into develop

8687875

fixed an issue causing LCOMP=2 to fail

c27aeb3

Merge remote-tracking branch 'remotes/upstream/develop' into develop

c298f2e

Sampling covariance working for one nuclide

756690c

Sampling for almost all possible File 32

00b12a1

added function to handle dicts of values

9156410

removed print statements

1d8ac57

Added more resonance sampling capability

77b5149

Sampling and reconstructing working for large files

0a1408d

Fixed l-values, added subset capability

9e321f7

more l-value changes

a85829e

Restructuring of classes, better handling of file 2 contribution

3647306

Sampling and subset functions working

a0182e1

Fixed a lot of doc strings, jupyter example

1845edb

Pull from upstream before testing

0bcb6f8

added reconstruction ability for resonance samples

84ef96c

Small change to endf parser

e88f8cd

Fixed MLBW sample reconstruction

6b836a1

MLBW test added

27c9d3c

Added tests for RM covariance

ebd66a4

some spelling

447d85a

Some cleanup of dependencies

b197b2b

style

219ea89

Added to sphinx autodocs

f69ab27

icmeyer added 9 commits July 18, 2018 16:59

Style

15ffa56

even more style

1beab30

Style and changed sampling to produce ResonanceRange objects

c22b36e

More reconstruction fixes

e485211

Fixes to tests

9ad8dab

Take advantage of np.random.multivariate size option for multiple sam…

d148c24

…ples

Style and change of __init__ methods

fde13a9

Stopped loop from erroneously adding range for unresolved paramaters …

81a15a9

…as previous range

added warning for sampling/reconstruction

010acb8

paulromano reviewed Jul 23, 2018

View reviewed changes

icmeyer added 4 commits July 23, 2018 13:43

Style, changed sampling/subset methods to return new objects

5381ad4

Restructured ResonanceCovarianceRange class to contain corresponding …

e2a27f2

…file2 data as an attribute

Changed tests to reflect changes in module

8c320cf

Removed print statements from test

6b82852

paulromano reviewed Jul 25, 2018

View reviewed changes

icmeyer added 5 commits July 25, 2018 09:41

Fix of sampling routine, change dataframe .as_matrix to .values

78a411e

more changes to sampling

1bf5e7b

Condensed sampling method

3626d8e

Changed __copy__ method of ResonanceRange to mark parameters unprepared

28414ef

remove unused variable

5918f5b

paulromano reviewed Jul 30, 2018

View reviewed changes

icmeyer added 3 commits July 31, 2018 08:29

Some name changes, removed redundandancies

c6c47e1

Fixed naming change in tests

c4b28e4

Merge branch 'develop' of https://github.com/mit-crpg/openmc into dev…

7190792

…elop

Added more tests for covariance module

46e58ef

paulromano approved these changes Aug 3, 2018

View reviewed changes

paulromano merged commit cbb9ace into openmc-dev:develop Aug 3, 2018

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Resonace Covariance Module #1024

Resonace Covariance Module #1024

icmeyer commented Jul 17, 2018 •

edited

Loading

icmeyer commented Jul 20, 2018

paulromano left a comment

paulromano Jul 23, 2018

paulromano Jul 23, 2018

paulromano Jul 23, 2018

paulromano Jul 23, 2018

paulromano Jul 23, 2018

paulromano Jul 23, 2018

paulromano Jul 23, 2018

paulromano Jul 23, 2018

paulromano Jul 23, 2018

paulromano Jul 23, 2018

paulromano Jul 23, 2018

paulromano Jul 25, 2018

icmeyer Jul 25, 2018

paulromano Jul 25, 2018

paulromano left a comment

paulromano Jul 30, 2018

paulromano Jul 30, 2018

paulromano Jul 30, 2018

paulromano Jul 30, 2018

paulromano Jul 30, 2018

paulromano Jul 30, 2018

paulromano Jul 30, 2018

paulromano Jul 30, 2018

icmeyer Jul 31, 2018

icmeyer commented Aug 1, 2018

		@@ -71,6 +71,24 @@ def float_endf(s):
		return float(_ENDF_FLOAT_RE.sub(r'\1e\2', s))


		def int_endf(s):

		@@ -249,6 +267,50 @@ def get_tab2_record(file_obj):

		return params, Tabulated2D(breakpoints, interpolation)

		def get_intg_record(file_obj):

		@@ -0,0 +1,743 @@
		from collections import defaultdict, MutableSequence, Iterable

Resonace Covariance Module #1024

Resonace Covariance Module #1024

Conversation

icmeyer commented Jul 17, 2018 • edited Loading

icmeyer commented Jul 20, 2018

paulromano left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

paulromano left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

icmeyer commented Aug 1, 2018

icmeyer commented Jul 17, 2018 •

edited

Loading