.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/uncertainty_contributions.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_uncertainty_contributions.py: Uncertainty/Variance Contributions ================================== This example demonstrates how to read (nuclear-data induced) uncertainty/variance contributions from a set of TSUNAMI-B SDF files. .. GENERATED FROM PYTHON SOURCE LINES 6-19 .. code-block:: Python # Getting the Contributions # ------------------------- # The first step is to get the contributions to the uncertainty. This can be done using the # :func:`tsunami_ip_utils.integral_indices.get_uncertainty_contributions` function. from tsunami_ip_utils.integral_indices import get_uncertainty_contributions from paths import EXAMPLES application_filenames = [ f"{EXAMPLES}/data/example_sdfs/HMF/HEU-MET-FAST-003-00{i}.sdf" for i in range(1, 3) ] experiment_filenames = [ f"{EXAMPLES}/data/example_sdfs/HMF/HEU-MET-FAST-003-00{i}.sdf" for i in range(3, 6) ] uncertainty_contributions_nuclide, uncertainty_contributions_nuclide_reaction = \ get_uncertainty_contributions(application_filenames, experiment_filenames) .. GENERATED FROM PYTHON SOURCE LINES 20-24 As explained in the API documentation, the above function returns two dictionaries. The first dictionary, ``uncertainty_contributions_nuclide``, contains the contributions to the uncertainty for each nuclide. The second dictionary, ``uncertainty_contributions_nuclide_reaction``, contains the contributions to the uncertainty for each nuclide-reaction pair. Each dictionary has keys ``'application'`` and ``'experiment'``. .. GENERATED FROM PYTHON SOURCE LINES 24-28 .. code-block:: Python print(uncertainty_contributions_nuclide.keys()) application_nuclide_contributions = uncertainty_contributions_nuclide['application'] .. rst-class:: sphx-glr-script-out .. code-block:: none dict_keys(['application', 'experiment', 'filenames']) .. GENERATED FROM PYTHON SOURCE LINES 29-34 Isotope-Wise Contributions -------------------------- The contents of the dictionary is a list of contributions for each application/experiment. I.e. the (nuclide-wise) contributions to the nuclear data induced uncertainty for application 1 (i.e. ``application_filenames[0]``) can be accessed as follows: .. GENERATED FROM PYTHON SOURCE LINES 34-37 .. code-block:: Python print(application_nuclide_contributions[1]) .. rst-class:: sphx-glr-script-out .. code-block:: none [{'isotope': 'u-235 - u-235', 'contribution': 1.1869831817995358+/-0.0009513935878905702}, {'isotope': 'u-238 - u-238', 'contribution': 0.5665193328849952+/-0.004831934771704849}, {'isotope': 'u-234 - u-234', 'contribution': 0.16184626300599841+/-0.00012274171911905328}, {'isotope': 'u-235 - u-238', 'contribution': 0.050524528884493325+/-3.775832291069502e-06}] .. GENERATED FROM PYTHON SOURCE LINES 38-41 The ouput is a list of dictionaries with the keys ``'isotope'`` and ``'contribution'``. The contributions are the nuclide-wise contributions to the nuclear data induced uncertainty. For TSUNAMI-B formatted SDFs, these values have an associated (monte carlo) uncertainty, and so are represented as :func:`uncertainties.ufloat` objects. .. GENERATED FROM PYTHON SOURCE LINES 41-44 .. code-block:: Python print(application_nuclide_contributions[1][0]['contribution']) .. rst-class:: sphx-glr-script-out .. code-block:: none 1.1870+/-0.0010 .. GENERATED FROM PYTHON SOURCE LINES 45-47 These objects automatically handle uncertainty propagation, and there use is further documented in the :mod:`uncertainties` package. .. GENERATED FROM PYTHON SOURCE LINES 49-52 Isotope-Reaction-Wise Contributions ----------------------------------- Isotope-reaction-wise contributions are accessed similarly. .. GENERATED FROM PYTHON SOURCE LINES 52-56 .. code-block:: Python application_nuclide_contributions = uncertainty_contributions_nuclide_reaction['application'] print(application_nuclide_contributions[1]) .. rst-class:: sphx-glr-script-out .. code-block:: none [{'isotope': 'u-235 - u-235', 'reaction_type': 'n,gamma - n,gamma', 'contribution': 1.0762+/-0.0010083}, {'isotope': 'u-238 - u-238', 'reaction_type': "n,n' - n,n'", 'contribution': 0.89654+/-0.0022208}, {'isotope': 'u-238 - u-238', 'reaction_type': "n,n' - elastic", 'contribution': -0.82886+/-0.0022158}, {'isotope': 'u-238 - u-238', 'reaction_type': 'elastic - elastic', 'contribution': 0.43714+/-0.00090315}, {'isotope': 'u-235 - u-235', 'reaction_type': "n,n' - n,n'", 'contribution': 0.31603+/-0.00044944}, {'isotope': 'u-235 - u-235', 'reaction_type': 'chi - chi', 'contribution': 0.28691+/-0.00081625}, {'isotope': 'u-235 - u-235', 'reaction_type': "n,n' - elastic", 'contribution': -0.25501+/-0.00029919}, {'isotope': 'u-235 - u-235', 'reaction_type': 'fission - fission', 'contribution': 0.24183+/-2.7144e-05}, {'isotope': 'u-235 - u-235', 'reaction_type': 'elastic - n,gamma', 'contribution': 0.23559+/-0.00054821}, {'isotope': 'u-234 - u-234', 'reaction_type': 'fission - fission', 'contribution': 0.15883+/-0.00012507}, {'isotope': 'u-235 - u-235', 'reaction_type': 'elastic - elastic', 'contribution': 0.12347+/-0.00011876}, {'isotope': 'u-238 - u-238', 'reaction_type': 'chi - chi', 'contribution': 0.082855+/-0.00013436}, {'isotope': 'u-235 - u-235', 'reaction_type': 'nubar - nubar', 'contribution': 0.080153+/-3.2557e-06}, {'isotope': 'u-238 - u-238', 'reaction_type': 'nubar - nubar', 'contribution': 0.070821+/-9.4547e-06}, {'isotope': 'u-235 - u-238', 'reaction_type': 'fission - fission', 'contribution': 0.059862+/-3.1781e-06}, {'isotope': 'u-235 - u-235', 'reaction_type': 'elastic - fission', 'contribution': -0.046224+/-7.8172e-06}, {'isotope': 'u-235 - u-238', 'reaction_type': 'fission - n,gamma', 'contribution': -0.032105+/-4.4042e-07}, {'isotope': 'u-238 - u-238', 'reaction_type': 'fission - fission', 'contribution': 0.023509+/-1.3028e-06}, {'isotope': 'u-238 - u-238', 'reaction_type': 'n,gamma - n,gamma', 'contribution': 0.022796+/-1.0321e-06}, {'isotope': 'u-234 - u-234', 'reaction_type': 'nubar - nubar', 'contribution': 0.020862+/-1.5998e-06}, {'isotope': 'u-234 - u-234', 'reaction_type': 'n,gamma - n,gamma', 'contribution': 0.015236+/-1.0777e-07}, {'isotope': 'u-234 - u-234', 'reaction_type': "n,n' - n,n'", 'contribution': 0.013031+/-7.1888e-06}, {'isotope': 'u-238 - u-238', 'reaction_type': 'elastic - n,gamma', 'contribution': 0.010215+/-4.2262e-07}, {'isotope': 'u-234 - u-234', 'reaction_type': 'chi - chi', 'contribution': 0.010136+/-8.1112e-06}, {'isotope': 'u-238 - u-238', 'reaction_type': 'n,2n - n,2n', 'contribution': 0.0090309+/-1.8261e-06}, {'isotope': 'u-235 - u-235', 'reaction_type': 'n,2n - n,2n', 'contribution': 0.0070606+/-6.1117e-07}, {'isotope': 'u-238 - u-238', 'reaction_type': 'elastic - fission', 'contribution': -0.006939+/-2.54e-07}, {'isotope': 'u-234 - u-234', 'reaction_type': 'elastic - elastic', 'contribution': 0.0052285+/-1.3973e-06}, {'isotope': 'u-238 - u-238', 'reaction_type': 'n,2n - elastic', 'contribution': -0.0035981+/-8.4954e-07}, {'isotope': 'u-235 - u-235', 'reaction_type': 'n,2n - elastic', 'contribution': -0.0025521+/-3.5859e-07}, {'isotope': 'u-235 - u-235', 'reaction_type': 'fission - n,gamma', 'contribution': 5.4626e-05+/-3.3895e-10}, {'isotope': 'u-234 - u-234', 'reaction_type': 'n,2n - n,2n', 'contribution': 5.0396e-05+/-1.6122e-09}, {'isotope': 'u-238 - u-238', 'reaction_type': 'fission - n,gamma', 'contribution': -1.2874e-05+/-1.0338e-11}] .. GENERATED FROM PYTHON SOURCE LINES 57-59 The output is a list of dictionaries with the keys ``'isotope'``, ``'reaction'``, and ``'contribution'``. The contributions are (like before) :func:`uncertainties.ufloat` objects. A specific contribution can be accessed via .. GENERATED FROM PYTHON SOURCE LINES 59-63 .. code-block:: Python print(application_nuclide_contributions[1][0]) print(application_nuclide_contributions[1][0]['contribution']) .. rst-class:: sphx-glr-script-out .. code-block:: none {'isotope': 'u-235 - u-235', 'reaction_type': 'n,gamma - n,gamma', 'contribution': 1.0762+/-0.0010083} 1.0762+/-0.0010 .. GENERATED FROM PYTHON SOURCE LINES 64-69 Getting Contributions for a Single File --------------------------------------- As discussed in the API documentation, the application_filenames, and experiment_filenames arugments are both optional - to streamline the process of getting contributions for a single file, the function can be called with only one argument. For example, to get the contributions for the first application file, we can call the function as follows: .. GENERATED FROM PYTHON SOURCE LINES 69-75 .. code-block:: Python application_filenames = [ f"{EXAMPLES}/data/example_sdfs/HMF/HEU-MET-FAST-003-001.sdf" ] uncertainty_contributions_nuclide, uncertainty_contributions_nuclide_reaction = \ get_uncertainty_contributions(application_filenames) print(uncertainty_contributions_nuclide['application']) .. rst-class:: sphx-glr-script-out .. code-block:: none [[{'isotope': 'u-235 - u-235', 'contribution': 1.1547475051845695+/-0.000896155210593727}, {'isotope': 'u-238 - u-238', 'contribution': 0.556418115404617+/-0.004135126369373132}, {'isotope': 'u-234 - u-234', 'contribution': 0.1570725325098565+/-0.00011949411284427181}, {'isotope': 'u-235 - u-238', 'contribution': 0.04953563217523321+/-3.3690884432420322e-06}]] .. GENERATED FROM PYTHON SOURCE LINES 76-78 Note that since we only passed application filenames, the function only returns the contributions for the application. The experiment contributions are an empty list .. GENERATED FROM PYTHON SOURCE LINES 78-81 .. code-block:: Python print(uncertainty_contributions_nuclide['experiment']) .. rst-class:: sphx-glr-script-out .. code-block:: none [] .. GENERATED FROM PYTHON SOURCE LINES 82-86 Variance Contributions ---------------------- The contributions to the nuclear data induced variance can also be obtained by passing the ``variance=True`` argument to the function. This will return the contributions to the variance, rather than the uncertainty. .. GENERATED FROM PYTHON SOURCE LINES 86-92 .. code-block:: Python variance_contributions_nuclide, variance_contributions_nuclide_reaction = \ get_uncertainty_contributions(application_filenames, variance=True) print(variance_contributions_nuclide['application']) .. rst-class:: sphx-glr-script-out .. code-block:: none [[{'isotope': 'u-235 - u-235', 'contribution': 1.3334418007299873+/-0.0020696659873825174}, {'isotope': 'u-238 - u-238', 'contribution': 0.30960111915042576+/-0.004601718442813069}, {'isotope': 'u-234 - u-234', 'contribution': 0.024671780469059927+/-3.75384858489367e-05}, {'isotope': 'u-235 - u-238', 'contribution': 0.002453778855+/-3.337798517805328e-07}]] .. GENERATED FROM PYTHON SOURCE LINES 93-95 The variances are just the squares of the uncertainties, except when the are negative, in which case the variance is the negative of the square of the uncertainty. .. GENERATED FROM PYTHON SOURCE LINES 95-103 .. code-block:: Python variance_contribution = variance_contributions_nuclide['application'][0][0]['contribution'] uncertainty_contribution = uncertainty_contributions_nuclide['application'][0][0]['contribution'] squared_uncertainty = uncertainty_contribution**2 if uncertainty_contribution >= 0 else -(uncertainty_contribution)**2 print(variance_contribution, squared_uncertainty) assert ( variance_contribution.n == squared_uncertainty.n ) and ( variance_contribution.s == squared_uncertainty.s ) .. rst-class:: sphx-glr-script-out .. code-block:: none 1.3334+/-0.0021 1.3334+/-0.0021 .. GENERATED FROM PYTHON SOURCE LINES 104-107 Note that even though the variance contribution and the squared uncertainty contribution have the same nominal value and standard deviation, they are not equal, i.e. ``variance_contribution == squared_uncertainty`` evaluates to ``False``. This is because the uncertainties are represented as :func:`uncertainties.ufloat` which represent distinct random variables, as discussed in `the uncertainties package documentation `_. .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 10.449 seconds) .. _sphx_glr_download_auto_examples_uncertainty_contributions.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: uncertainty_contributions.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: uncertainty_contributions.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: uncertainty_contributions.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_