diff --git a/notebooks/notebook_Edistr_custom_perturbation.ipynb b/notebooks/notebook_Edistr_custom_perturbation.ipynb new file mode 100644 index 00000000..83adc1bf --- /dev/null +++ b/notebooks/notebook_Edistr_custom_perturbation.ipynb @@ -0,0 +1,133 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import sandy" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import pandas as pd\n", + "import numpy as np\n", + "import matplotlib.pyplot as plt" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "%load_ext autoreload\n", + "%autoreload 2" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "col = columns = pd.MultiIndex.from_product([[9237], [18], [0], [1.00000e-05], [2.00000e-05]], names=['MAT', 'MT', 'K', 'ELO', 'EHI'])\n", + "pert = pd.DataFrame([1.3, 1.5, 1.2, 1.1, 1, 0.7, 1.6, 0.5, 1.8, 1.5, 0.6, 1.5, 0.8], \n", + " index=sandy.energy_grids.CASMO12,\n", + " columns=col)\n", + "pert = sandy.Pert(pert)\n", + "pert" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Custom perturbation of U-238:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "extra_points = np.logspace(-5, 7, 49)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "edistr = sandy.Edistr.from_endf6(sandy.get_endf6_file(\"jeff_33\", \"xs\", 922380)).reshape(extra_points)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "edistr_pert = edistr.custom_perturbation(pert)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "mask = (9237 == edistr.data.MAT) & \\\n", + " (18 == edistr.data.MT) & \\\n", + " (0 == edistr.data.K) & \\\n", + " (1.00000e-05 == edistr.data.EIN)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(figsize=(7, 3.5), dpi=100)\n", + "edistr.data[mask][['EOUT', \"VALUE\"]].set_index('EOUT').plot(logx=True, logy=True, alpha=0.8, linewidth=0.9, ax=ax)\n", + "edistr_pert.data[mask][['EOUT', \"VALUE\"]].set_index('EOUT').plot(logx=True, logy=True, alpha=0.8, linewidth=0.9, ax=ax)\n", + "ax.set_ylabel(\"Legendre polynomial coefficients for P=1\")\n", + "ax.set_xlabel(\"energy / eV\")\n", + "fig.tight_layout()" + ] + } + ], + "metadata": { + "interpreter": { + "hash": "d167b81dab4cfdd7e83bd7f0413ce20e956b35802f4873808516fc813094dcfb" + }, + "kernelspec": { + "display_name": "Python 3.7.10 ('sandy-devel')", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.7.10" + }, + "orig_nbformat": 4 + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/notebooks/notebook_Lpc_custom_perturbation.ipynb b/notebooks/notebook_Lpc_custom_perturbation.ipynb new file mode 100644 index 00000000..448b6357 --- /dev/null +++ b/notebooks/notebook_Lpc_custom_perturbation.ipynb @@ -0,0 +1,111 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import sandy" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import pandas as pd\n", + "import numpy as np\n", + "import matplotlib.pyplot as plt" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "col = pd.MultiIndex.from_arrays([[2], [9237], [1]], names=('MT', 'MAT', 'P'))\n", + "pert = pd.DataFrame([1.3, 1.5, 1.2, 1.1, 1, 0.7, 1.6, 0.5, 1.8, 1.5, 0.6, 1.5, 0.8], \n", + " index=sandy.energy_grids.CASMO12,\n", + " columns=col)\n", + "pert = sandy.Pert(pert)\n", + "pert" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Custom perturbation of U-238:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "extra_points = np.logspace(-5, 7, 49)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "lpc = sandy.Lpc.from_endf6(sandy.get_endf6_file(\"jeff_33\", \"xs\", 922380)).reshape(extra_points)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "lpc_pert = lpc.custom_perturbation(pert)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(figsize=(7, 3.5), dpi=100)\n", + "lpc.data.loc[(9237, 2), (1)].plot(logx=True, logy=True, alpha=0.8, linewidth=0.9, ax=ax)\n", + "lpc_pert.data.loc[(9237, 2), (1)].plot(logx=True, logy=True, alpha=0.8, linewidth=0.9, ax=ax)\n", + "ax.set_ylabel(\"Legendre polynomial coefficients for P=1\")\n", + "ax.set_xlabel(\"energy / eV\")\n", + "fig.tight_layout()" + ] + } + ], + "metadata": { + "interpreter": { + "hash": "d167b81dab4cfdd7e83bd7f0413ce20e956b35802f4873808516fc813094dcfb" + }, + "kernelspec": { + "display_name": "Python 3.7.10 ('sandy-devel')", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.7.10" + }, + "orig_nbformat": 4 + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/notebooks/notebook_Xs_custom_perturbation.ipynb b/notebooks/notebook_Xs_custom_perturbation.ipynb new file mode 100644 index 00000000..c538e342 --- /dev/null +++ b/notebooks/notebook_Xs_custom_perturbation.ipynb @@ -0,0 +1,100 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import sandy" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import pandas as pd\n", + "import matplotlib.pyplot as plt" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "col = pd.MultiIndex.from_arrays([[2], [125]], names=('MT', 'MAT'))\n", + "pert = pd.DataFrame([1.3, 1.5, 1.2, 1.1, 1, 0.7, 1.6, 0.5, 1.8, 1.5, 0.6, 1.5, 0.8], \n", + " index=sandy.energy_grids.CASMO12,\n", + " columns=col)\n", + "pert = sandy.Pert(pert)\n", + "pert" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Custom perturbation of H:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "xs = sandy.Xs.from_endf6(sandy.get_endf6_file(\"jeff_33\", \"xs\", 10010))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "xs_pert = xs.custom_perturbation(pert)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(figsize=(7, 3.5), dpi=100)\n", + "xs.data.plot(logx=True, logy=True, alpha=0.8, linewidth=0.9, ax=ax)\n", + "xs_pert.data.plot(logx=True, logy=True, alpha=0.8, linewidth=0.9, ax=ax)\n", + "ax.set_ylabel(\"cross section / b\")\n", + "ax.set_xlabel(\"energy / eV\")\n", + "ax.set_xlim([1e-5, 2e7])\n", + "ax.set_ylim([1e-1, 1e4])\n", + "fig.tight_layout()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python3 (sandy-devel)", + "language": "python", + "name": "sandy-devel" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.7.10" + }, + "orig_nbformat": 4 + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/notebooks/notebook_from_sandy.Samples_to_sandy.Pert.ipynb b/notebooks/notebook_from_sandy.Samples_to_sandy.Pert.ipynb new file mode 100644 index 00000000..da8c3b21 --- /dev/null +++ b/notebooks/notebook_from_sandy.Samples_to_sandy.Pert.ipynb @@ -0,0 +1,179 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import sandy" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import pandas as pd" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "XS samples:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "samples = pd.DataFrame([[0.5, 1.5], [0.75, 1.25]] , index=[0, 1])\n", + "samples.columns = pd.MultiIndex.from_arrays([[125, 125], [2, 2], [1e-05, 19970500.0]], names=['MAT', 'MT', 'E'])\n", + "samples = sandy.Samples(samples)\n", + "samples" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "samples.to_pert()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Single sample:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "samples.to_pert(smp=0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "LPC samples:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "samples = pd.DataFrame([[0.5, 1.5], [0.75, 1.25]] , index=[0, 1])\n", + "samples.columns = pd.MultiIndex.from_arrays([[125, 125], [2, 2], [1.00000e-05, 1.00000e-05], [1.00000e-05, \t1.00000e-05], [1e-05, 19970500.0]], names=['MAT', 'MT', 'ELO', 'EHI', 'E'])\n", + "samples = sandy.Samples(samples)\n", + "samples" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "samples.to_pert()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Single sample:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "samples.to_pert(smp=0)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Energy distribution samples:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "samples = pd.DataFrame([[0.5, 1.5], [0.75, 1.25]] , index=[0, 1])\n", + "samples.columns = pd.MultiIndex.from_arrays([[125, 125], [2, 2], [1, 1], [1e-05, 19970500.0]], names=['MAT', 'MT', 'L', 'E'])\n", + "samples = sandy.Samples(samples)\n", + "samples" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "samples.to_pert()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Single sample:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "samples.to_pert(smp=0)\n" + ] + } + ], + "metadata": { + "interpreter": { + "hash": "d167b81dab4cfdd7e83bd7f0413ce20e956b35802f4873808516fc813094dcfb" + }, + "kernelspec": { + "display_name": "Python 3.7.10 ('sandy-devel')", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.7.10" + }, + "orig_nbformat": 4 + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/notebooks/notebook_reorder_of_sandy.Pert.ipynb b/notebooks/notebook_reorder_of_sandy.Pert.ipynb new file mode 100644 index 00000000..896a1d1e --- /dev/null +++ b/notebooks/notebook_reorder_of_sandy.Pert.ipynb @@ -0,0 +1,127 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import sandy" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import pandas as pd" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Rearrange `sandy.Pert` columns levels using input order" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "col = pd.MultiIndex.from_arrays([[5, 2, 2, 3], [2631, 2631, 125, 125]], names=('MT', 'MAT'))\n", + "pert = pd.DataFrame([[1, 1.05, 0.95, 1.30], [1.05, 1, 1.15, 0.7]], \n", + " index=pd.IntervalIndex.from_breaks(pd.Index([1.94000e+08, 1.96000e+08+1]).insert(0, 0)),\n", + " columns=col)\n", + "pert = sandy.Pert(pert)\n", + "pert" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Reorder for custom perturbation of H:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "xs = sandy.Xs.from_endf6(sandy.get_endf6_file(\"jeff_33\", \"xs\", 10010))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "pert.reorder(xs.data.columns)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Reorder for custom perturbation of Fe:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "xs = sandy.Xs.from_endf6(sandy.get_endf6_file(\"jeff_33\", \"xs\", 260560))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "pert.reorder(xs.data.columns)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "\n" + ] + } + ], + "metadata": { + "interpreter": { + "hash": "d167b81dab4cfdd7e83bd7f0413ce20e956b35802f4873808516fc813094dcfb" + }, + "kernelspec": { + "display_name": "Python 3.7.10 ('sandy-devel')", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.7.10" + }, + "orig_nbformat": 4 + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/notebooks/notebook_reshape_of_sandy.Pert.ipynb b/notebooks/notebook_reshape_of_sandy.Pert.ipynb new file mode 100644 index 00000000..d3dabd14 --- /dev/null +++ b/notebooks/notebook_reshape_of_sandy.Pert.ipynb @@ -0,0 +1,105 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import sandy" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "import pandas as pd\n", + "import matplotlib.pyplot as plt\n", + "import seaborn as sns \n", + "import numpy as np\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Reshape of `sandy.Pert` object" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "pert = pd.Series([1, 1.05], index=pd.IntervalIndex.from_breaks(pd.Index([10, 100]).insert(0, 0)))\n", + "pert = sandy.Pert(pert)\n", + "plot = pert.right.reset_index().rename(columns={\"index\": \"E\", 0:'Perturbation'})" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(figsize=(4, 4), dpi=100)\n", + "plot.plot(kind=\"scatter\", x=\"E\", y='Perturbation', ax=ax, grid=True, legend=None)\n", + "plt.xticks(np.arange(0, 130, 10))\n", + "plt.show()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Reshape:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "pert_reshape = pert.reshape(np.arange(0, 120)).right.reset_index().rename(columns={\"index\": \"E\", 0:'Perturbation'})" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "fig, ax = plt.subplots(figsize=(4, 4), dpi=100)\n", + "pert_reshape.plot(kind=\"scatter\", x=\"E\", y='Perturbation', ax=ax, grid=True, legend=None)\n", + "plt.xticks(np.arange(0, 130, 10))\n", + "plt.show()" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python3 (sandy-devel)", + "language": "python", + "name": "sandy-devel" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.7.10" + }, + "orig_nbformat": 4 + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/sandy/core/cov.py b/sandy/core/cov.py index 91dfcea5..27eb0469 100644 --- a/sandy/core/cov.py +++ b/sandy/core/cov.py @@ -683,16 +683,18 @@ def sampling(self, nsmp, seed=None, rows=None, pdf='normal', -------- Draw 3 sets of samples using custom seed: >>> sandy.CategoryCov([[1, 0.4],[0.4, 1]]).sampling(3, seed=11) - 0 1 - 0 -7.49455e-01 -2.13159e+00 - 1 1.28607e+00 1.10684e+00 - 2 1.48457e+00 9.00879e-01 + 0 1 + SMP + 0 -7.49455e-01 -2.13159e+00 + 1 1.28607e+00 1.10684e+00 + 2 1.48457e+00 9.00879e-01 >>> sandy.CategoryCov([[1, 0.4],[0.4, 1]]).sampling(3, seed=11, rows=1) - 0 1 - 0 -7.49455e-01 -2.13159e+00 - 1 1.28607e+00 1.10684e+00 - 2 1.48457e+00 9.00879e-01 + 0 1 + SMP + 0 -7.49455e-01 -2.13159e+00 + 1 1.28607e+00 1.10684e+00 + 2 1.48457e+00 9.00879e-01 >>> sample = sandy.CategoryCov([[1, 0.4],[0.4, 1]]).sampling(1000000, seed=11) >>> sample.data.cov() @@ -702,10 +704,11 @@ def sampling(self, nsmp, seed=None, rows=None, pdf='normal', Small negative eigenvalue: >>> sandy.CategoryCov([[1, -.2],[-.2, 3]]).sampling(3, seed=11, tolerance=0) - 0 1 - 0 2.74945e+00 5.21505e+00 - 1 7.13927e-01 1.07147e+00 - 2 5.15435e-01 1.64683e+00 + 0 1 + SMP + 0 2.74945e+00 5.21505e+00 + 1 7.13927e-01 1.07147e+00 + 2 5.15435e-01 1.64683e+00 >>> sandy.CategoryCov([[1, -.2],[-.2, 3]]).sampling(1000000, seed=11, tolerance=0).data.cov() 0 1 @@ -714,16 +717,18 @@ def sampling(self, nsmp, seed=None, rows=None, pdf='normal', Sampling with different `pdf`: >>> sandy.CategoryCov([[1, -.2],[-.2, 3]]).sampling(3, seed=11, pdf='uniform', tolerance=0) - 0 1 - 0 -1.07578e-01 2.34960e+00 - 1 -6.64587e-01 5.21222e-01 - 2 8.72585e-01 9.12563e-01 + 0 1 + SMP + 0 -1.07578e-01 2.34960e+00 + 1 -6.64587e-01 5.21222e-01 + 2 8.72585e-01 9.12563e-01 >>> sandy.CategoryCov([[1, .2],[.2, 3]]).sampling(3, seed=11, pdf='lognormal', tolerance=0) - 0 1 - 0 3.03419e+00 1.57919e+01 - 1 5.57248e-01 4.74160e-01 - 2 4.72366e-01 6.50840e-01 + 0 1 + SMP + 0 3.03419e+00 1.57919e+01 + 1 5.57248e-01 4.74160e-01 + 2 4.72366e-01 6.50840e-01 >>> sandy.CategoryCov([[1, -.2],[-.2, 3]]).sampling(1000000, seed=11, pdf='uniform', tolerance=0).data.cov() 0 1 diff --git a/sandy/core/lpc.py b/sandy/core/lpc.py index 7470010f..a94c6e0d 100644 --- a/sandy/core/lpc.py +++ b/sandy/core/lpc.py @@ -140,9 +140,13 @@ def filter_by(self, key, value): >>> comp = LPC.filter_by('E', 1e-05).data.index.get_level_values(2) == 1e-05 >>> assert comp.all() == True """ - info = {'MAT': 0, 'MT': 1, 'E': 2} - condition = self.data.index.get_level_values(info[key]) == value - out = self.data.copy()[condition] + value_ = [value] if isinstance(value, int) or isinstance(value, float) else value + if key == 'P': + condition = self.data.columns.get_level_values(key).isin(value_) + out = self.data.copy().loc[:, condition] + else: + condition = self.data.index.get_level_values(key).isin(value_) + out = self.data.copy()[condition] if out.empty: raise sandy.Error("applied filter returned empty dataframe") return self.__class__(out) @@ -177,21 +181,15 @@ def _filters(self, conditions): out = out.filter_by(keys, values) return out - def custom_perturbation(self, mat, mt, p, pert): + def custom_perturbation(self, pert): """ Apply a custom perturbation (energy dependent) to a given Legendre polynomial coefficient. Parameters ---------- - mat : `int` - MAT material number - mt : `int` - MT reaction number - p : `int` - order of the Legendre polynomial coefficient - pert : `sandy.Pert` - tabulated perturbations + pert : `sandy.Pert` or `pd.DataFrame` + tabulated perturbation Returns ------- @@ -205,24 +203,42 @@ def custom_perturbation(self, mat, mt, p, pert): >>> mat= 9228 >>> mt= 2 >>> p = 1 - >>> pert= sandy.Pert([1, 1.05], index=[1.00000e+03, 5.00000e+03]) - >>> pert = LPC.custom_perturbation(mat, mt, p, pert)._filters({'MAT': mat, 'MT':2}).data.loc[:,1].iloc[0:5] + >>> columns = pd.MultiIndex.from_product([[mat], [mt], [p]], names=['MAT', 'MT', 'P']) + >>> pert= sandy.Pert(pd.DataFrame([1, 1.05], index=[1.00000e+03, 5.00000e+03], columns=columns)) + >>> pert = LPC.custom_perturbation(pert)._filters({'MAT': mat, 'MT':2}).data.loc[:,1].iloc[0:5] >>> data = LPC._filters({'MAT': mat, 'MT':2}).data.loc[:,1].iloc[0:5] >>> (pert/data).fillna(0) MAT MT E 9228 2 1.00000e-05 0.00000e+00 + 7.73304e+01 0.00000e+00 1.00000e+03 1.00000e+00 2.00000e+03 1.05000e+00 5.00000e+03 1.05000e+00 - 1.00000e+04 1.00000e+00 + 1.00000e+04 0.00000e+00 Name: 1, dtype: float64 """ - eg = self.data.loc[(mat, mt)].index.values - enew = np.union1d(eg, pert.right.index.values) - u_lpc = self.reshape(enew, selected_mat=mat, selected_mt=mt) - u_pert = pert.reshape(enew) - u_lpc.data.loc[(mat, mt)][p] *= u_pert.right.values - return Lpc(u_lpc.data) + pert_ = sandy.Pert(pert) if not isinstance(pert, sandy.Pert) else pert + if 'L' in pert_.data.columns.names: + pert_.data.columns = pert_.data.columns.set_names('P', level='L') + + # New energy grid: + energy_grid = self.data.index.get_level_values('E').unique() + enew = energy_grid.union(pert.right.index) + enew = enew[(enew <= energy_grid.max()) & (enew >= energy_grid.min())] + + # Reshape to new energy grid and columns estructure: + mat = pert.data.columns.get_level_values('MAT').unique().values + mt = pert.data.columns.get_level_values('MT').unique().values + u_lpc = self.reshape(enew, selected_mat=mat, selected_mt=mt)\ + .data.unstack(level=['MAT', 'MT']) + u_pert = pert_.reorder(u_lpc.columns).reshape(enew).right + + # Aply perturbation: + u_lpc.loc[:, u_pert.columns] = u_lpc.loc[:, u_pert.columns]\ + .multiply(u_pert, axis='columns') + u_lpc = u_lpc.stack(level=['MAT', 'MT'])\ + .reorder_levels(self.data.index.names, axis=0) + return self.__class__(u_lpc) def reshape(self, eg, selected_mat=None, selected_mt=None): """ @@ -258,34 +274,31 @@ def reshape(self, eg, selected_mat=None, selected_mt=None): 0 1.00000e-05 0.00000e+00 0.00000e+00 1 1.00000e+00 1.13380e-06 2.53239e-09 2 2.00000e+00 2.26761e-06 5.06481e-09 - 3 1.00000e+03 1.13381e-03 2.53242e-06 - 4 2.00000e+03 2.93552e-03 1.59183e-05 + 3 7.73304e+01 8.76780e-05 1.95833e-07 + 4 1.00000e+03 1.13381e-03 2.53242e-06 """ - listdf = [] - for (mat, mt), df in self.data.groupby(["MAT", "MT"]): - if selected_mat: - if mat != selected_mat: - listdf.append(df) - continue - if selected_mt: - if mt != selected_mt: - listdf.append(df) - continue - df = df.T[mat][mt].T - enew = df.index.union(eg).astype("float").rename("E") + lpc = self + if selected_mat: + selected_mat_ = [selected_mat] if isinstance(selected_mat, int) else selected_mat + lpc = lpc.filter_by('MAT', selected_mat_) + if selected_mt: + selected_mt_ = [selected_mt] if isinstance(selected_mt, int) else selected_mt + lpc = lpc.filter_by('MT', selected_mt_) + energy_grid = lpc.data.index.get_level_values('E').unique() + enew = energy_grid.union(pd.Index(eg)).astype("float").rename("E") + + def reshape_lpc(df, eg): valsnew = sandy.shared.reshape_differential( - df.index.values, + df.index.get_level_values('E').values, df.values, enew, ) - dfnew = pd.DataFrame(valsnew, index=enew, columns=df.columns) \ - .reset_index(drop=False) - dfnew.insert(0, "MAT", mat) - dfnew.insert(0, "MT", mt) - dfnew.set_index(self._indexnames, inplace=True) - listdf.append(dfnew) - data = pd.concat(listdf, axis=0) - return Lpc(data) + dfnew = pd.DataFrame(valsnew, index=enew, columns=df.columns) + return dfnew.loc[(dfnew.index >= energy_grid.min()) & + (dfnew.index <= energy_grid.max()), :] + + data = lpc.data.groupby(["MAT", "MT"]).apply(reshape_lpc, enew) + return self.__class__(data) def to_endf6(self, endf6): """ @@ -435,7 +448,7 @@ def _add_points(self, extra_points): rdf["MT"] = mt rdf = rdf.set_index(["MAT", "MT", "E"]) List.append(rdf) - return Lpc(pd.concat(List, axis=0)) + return self.__class__(pd.concat(List, axis=0)) def _perturb(self, pert, method=2, **kwargs): """ diff --git a/sandy/core/samples.py b/sandy/core/samples.py index a7021b50..4b8da7e9 100644 --- a/sandy/core/samples.py +++ b/sandy/core/samples.py @@ -89,6 +89,7 @@ def data(self): @data.setter def data(self, data): self._data = data + self._data.index.name = 'SMP' @property def condition_number(self): @@ -262,3 +263,98 @@ def from_csv(cls, file, **kwargs): """ df = pd.read_csv(file, **kwargs) return cls(df) + + def to_pert(self, smp=None): + """ + Samples in the format of `Sandy.Pert`. If the sample number is not + entered, a `pd.DataFrame` is created with a multiindex, the first level + being the sample number and the second the energy grid. + + Parameters + ---------- + smp : `int`, optional + Number of the sample. The default is None. + + Returns + ------- + `sandy.Pert` or `pd.DataFrame` + Samples in `sandy.Pert` format. + + Examples + -------- + >>> samples = pd.DataFrame([[0.5, 1.5], [0.75, 1.25]] , index=[0, 1]) + >>> samples.columns = pd.MultiIndex.from_arrays([[125, 125], [2, 2], [1e-05, 19970500.0]], names=['MAT', 'MT', 'E']) + >>> sandy.Samples(samples).to_pert() + MAT 125 + MT 2 + SMP E + 0 1.00000e-05 5.00000e-01 + 1.99705e+07 1.50000e+00 + 1 1.00000e-05 7.50000e-01 + 1.99705e+07 1.25000e+00 + + >>> sandy.Samples(samples).to_pert(smp=0) + MAT 125 + MT 2 + ENERGY + (0.0, 1e-05] 5.00000e-01 + (1e-05, 19970500.0] 1.50000e+00 + + >>> samples = pd.DataFrame([[0.5, 1.5], [0.75, 1.25]] , index=[0, 1]) + >>> samples.columns = pd.MultiIndex.from_arrays([[125, 125], [2, 2], [1.00000e-05, 1.00000e-05], [1.00000e-05, 1.00000e-05], [1e-05, 19970500.0]], names=['MAT', 'MT', 'ELO', 'EHI', 'E']) + >>> sandy.Samples(samples).to_pert() + MAT 125 + MT 2 + ELO 1.00000e-05 + EHI 1.00000e-05 + SMP E + 0 1.00000e-05 5.00000e-01 + 1.99705e+07 1.50000e+00 + 1 1.00000e-05 7.50000e-01 + 1.99705e+07 1.25000e+00 + + >>> sandy.Samples(samples).to_pert(smp=0) + MAT 125 + MT 2 + ELO 1.00000e-05 + EHI 1.00000e-05 + ENERGY + (0.0, 1e-05] 5.00000e-01 + (1e-05, 19970500.0] 1.50000e+00 + + >>> samples = pd.DataFrame([[0.5, 1.5], [0.75, 1.25]] , index=[0, 1]) + >>> samples.columns = pd.MultiIndex.from_arrays([[125, 125], [2, 2], [1, 1], [1e-05, 19970500.0]], names=['MAT', 'MT', 'L', 'E']) + >>> sandy.Samples(samples).to_pert() + MAT 125 + MT 2 + L 1 + SMP E + 0 1.00000e-05 5.00000e-01 + 1.99705e+07 1.50000e+00 + 1 1.00000e-05 7.50000e-01 + 1.99705e+07 1.25000e+00 + + >>> sandy.Samples(samples).to_pert(smp=0) + MAT 125 + MT 2 + L 1 + ENERGY + (0.0, 1e-05] 5.00000e-01 + (1e-05, 19970500.0] 1.50000e+00 + """ + data = self.data.copy() + levels = list(np.arange(data.columns.nlevels)) + if smp is not None: + if isinstance(smp, int) and smp in data.index: + pert = data.loc[smp].to_frame().unstack(level=levels[0:-1]) + pert.columns = pert.columns.droplevel(level=None) + pert = sandy.Pert(pert) + else: + print(f"{smp} is not correct or do not exist") + else: + def foo(df): + df = df.stack() + df.index = df.index.droplevel(level='SMP') + return df + pert = data.groupby('SMP').apply(foo).fillna(0) + return pert diff --git a/sandy/core/xs.py b/sandy/core/xs.py index 0373268a..2f979c5c 100644 --- a/sandy/core/xs.py +++ b/sandy/core/xs.py @@ -156,38 +156,124 @@ def reshape(self, eg): df = pd.DataFrame(xsnew, index=enew, columns=df.columns) return self.__class__(df) - def custom_perturbation(self, mat, mt, pert): + def custom_perturbation(self, pert, **kwargs): """ - Apply a custom perturbation to a given cross section identified by - a MAT and MT number. + Function to apply perturbations to individual XS or a group of XS. Parameters ---------- - mat : `int` - MAT material number of the xs to which perturbations are to be - applied - mt : `int` - MT reaction number of the xs to which perturbations are to be - applied - pert : `sandy.Pert` - tabulated perturbations + pert : `sandy.Pert` or `pd.Dataframe` + tabulated perturbations with the index representing the energy grid + and the columns representing MT. + **kwargs : `dict` + keyword argument to pass to `sandy.Xs._recontruct_sums`. + + Parameters for _recontruct_sums + --------------------- + drop : `bool`, optional + Keep only mts present in the original file. The default is True. + inplace : `bool`, optional + Argument to define whether the output is a new object or + overwrite the original object.. The default is False. Returns ------- - `Xs` - cross section instance with given series MAT/MT perturbed + `sandy.Xs` + Xs disturbed and reconstructed. + + Examples + -------- + Test single perturbation: + >>> endf6 = sandy.get_endf6_file('jeff_33','xs', 10010) + >>> xs = sandy.Xs.from_endf6(endf6) + >>> col = pd.MultiIndex.from_arrays([[125], [2]], names=('MAT', 'MT')) + >>> pert = pd.DataFrame([1, 1.05], index =pd.IntervalIndex.from_breaks(pd.Index([10, 100]).insert(0, 0)), columns=col) + >>> pert_xs = xs.custom_perturbation(sandy.Pert(pert)) + >>> (pert_xs.data.loc[:, (125, 2)] / xs.data.loc[:, (125, 2)]).round(2).unique() + array([1. , 1.05]) + + >>> (pert_xs.data.loc[[1.00000e-05, 10, 20.0, 100, 10000.0], (125, 2)] / xs.data.loc[[1.00000e-05, 10, 20.0, 100, 10000.0], (125, 2)]).round(2) + E + 1.00000e-05 1.00000e+00 + 1.00000e+01 1.00000e+00 + 2.00000e+01 1.05000e+00 + 1.00000e+02 1.05000e+00 + 1.00000e+04 1.00000e+00 + Name: (125, 2), dtype: float64 + + Multiple perturbation: + Test single perturbation (Redundant): + >>> endf6 = sandy.get_endf6_file('jeff_33','xs', 260560) + >>> xs = sandy.Xs.from_endf6(endf6) + >>> col = pd.MultiIndex.from_arrays([[2631], [1]], names=('MAT', 'MT')) + >>> pert = pd.DataFrame([1, 1.05], index =pd.IntervalIndex.from_breaks(pd.Index([1.94000e+08, 1.96000e+08]).insert(0, 0)), columns=col) + >>> pert_xs = xs.custom_perturbation(sandy.Pert(pert)) + >>> (pert_xs.data.loc[[1.92000e+08, 1.94000e+08, 1.96000e+08, 1.98000e+08], [(2631, 1), (2631, 2), (2631, 3)]] / xs.data.loc[[1.92000e+08, 1.94000e+08, 1.96000e+08, 1.98000e+08], [(2631, 1), (2631, 2), (2631, 3)]]).round(2) + MAT 2631 + MT 1 2 3 + E + 1.92000e+08 1.00000e+00 1.00000e+00 1.00000e+00 + 1.94000e+08 1.00000e+00 1.00000e+00 1.00000e+00 + 1.96000e+08 1.03000e+00 1.05000e+00 1.00000e+00 + 1.98000e+08 1.00000e+00 1.00000e+00 1.00000e+00 + + Multiple perturbation(redundant + non redundant, mt perturb max = 3) + >>> col = pd.MultiIndex.from_arrays([[2631, 2631], [3, 2]], names=('MAT', 'MT')) + >>> pert = pd.DataFrame([[1, 1.05], [1.05, 1]], index =pd.IntervalIndex.from_breaks(pd.Index([1.94000e+08, 1.96000e+08+1]).insert(0, 0)), columns=col) + >>> pert_xs = xs.custom_perturbation(sandy.Pert(pert)) + >>> pert_xs.data.loc[[1.92000e+08, 1.94000e+08, 1.96000e+08, 1.98000e+08], [(2631, 1), (2631, 2), (2631, 3), (2631, 5)]] / xs.data.loc[[1.92000e+08, 1.94000e+08, 1.96000e+08, 1.98000e+08], [(2631, 1), (2631, 2), (2631, 3), (2631, 5)]] + MAT 2631 + MT 1 2 3 5 + E + 1.92000e+08 1.03353e+00 1.05000e+00 1.00121e+00 1.00000e+00 + 1.94000e+08 1.03360e+00 1.05000e+00 1.00106e+00 1.00000e+00 + 1.96000e+08 1.01702e+00 1.00000e+00 1.05085e+00 1.05000e+00 + 1.98000e+08 1.00016e+00 1.00000e+00 1.00044e+00 1.00000e+00 + + Multiple perturbation(non redundant): + >>> col = pd.MultiIndex.from_arrays([[2631, 2631], [5, 2]], names=('MAT', 'MT')) + >>> pert = pd.DataFrame([[1, 1.05], [1.05, 1]], index =pd.IntervalIndex.from_breaks(pd.Index([1.94000e+08, 1.96000e+08+1]).insert(0, 0)), columns=col) + >>> pert_xs = xs.custom_perturbation(sandy.Pert(pert)) + >>> pert_xs.data.loc[[1.92000e+08, 1.94000e+08, 1.96000e+08, 1.98000e+08], [(2631, 1), (2631, 2), (2631, 3), (2631, 5)]] / xs.data.loc[[1.92000e+08, 1.94000e+08, 1.96000e+08, 1.98000e+08], [(2631, 1), (2631, 2), (2631, 3), (2631, 5)]] + MAT 2631 + MT 1 2 3 5 + E + 1.92000e+08 1.03353e+00 1.05000e+00 1.00121e+00 1.00000e+00 + 1.94000e+08 1.03360e+00 1.05000e+00 1.00106e+00 1.00000e+00 + 1.96000e+08 1.01702e+00 1.00000e+00 1.05085e+00 1.05000e+00 + 1.98000e+08 1.00016e+00 1.00000e+00 1.00044e+00 1.00000e+00 """ - if (mat, mt) not in self.data: - msg = f"could not find MAT{mat}/MT{mt}, " +\ - "perturbation will not be applied" - logging.warning(msg) - u_xs = self - else: - enew = np.union1d(self.data.index.values, pert.right.index.values) - u_xs = self.reshape(enew) - u_pert = pert.reshape(enew) - u_xs.data[(mat, mt)] = u_xs.data[(mat, mt)] * u_pert.right.values - return self.__class__(u_xs.data) + pert_ = sandy.Pert(pert) if not isinstance(pert, sandy.Pert) else pert + # Reshape (all to the right): + index = self.data.index + enew = index.union(pert.right.index.values) + enew = enew[(enew <= index.max()) & (enew >= index.min())] + u_xs = self.reshape(enew).data + u_pert = pert_.reorder(self.data.columns).reshape(enew).right + # Redundant xs: + for mat in u_pert.columns.get_level_values('MAT'): + mt = u_xs[(mat)].columns.get_level_values('MT') + mt_pert_o = u_pert[(mat)].columns.get_level_values('MT') + parent = pd.Index(self.__class__.redundant_xs.keys()) + mask = mt_pert_o.isin(parent) + if mt_pert_o.max() == 3: + u_pert = u_pert.join(pd.concat([u_pert[(mat, 3)]]*len(mt[mt > 3]), + axis=1, + keys=zip([mat]*len(mt[mt > 3]), + mt[mt > 3]))) + elif len(mt_pert_o[mask]) != 0: + for mt_parent in mt_pert_o[mask].sort_values(ascending=False): + mt_pert = mt[(mt.isin(redundant_xs[mt_parent]))] + mt_pert = mt_pert[~(mt_pert.isin(mt_pert_o))] + if len(mt_pert) != 0: + u_pert = u_pert.join(pd.concat([u_pert[(mat, mt_parent)]]*len(mt_pert), + axis=1, + keys=zip([mat]*len(mt_pert), + mt_pert))) + # Apply pertubation: + u_xs.loc[:, u_pert.columns] = u_xs.loc[:, u_pert.columns]\ + .multiply(u_pert, axis='columns') + return self.__class__(u_xs)._reconstruct_sums(**kwargs) def filter_energies(self, energies): mask = self.data.index.isin(energies) @@ -362,6 +448,20 @@ def foo(l, r): def _reconstruct_sums(self, drop=True, inplace=False): """ Reconstruct redundant xs. + + Parameters + ---------- + drop : `bool`, optional + Keep only mts present in the original file. The default is True. + inplace : `bool`, optional + Argument to define whether the output is a new object or + overwrite the original object.. The default is False. + + Returns + ------- + `sandy.Xs` + Sandy object with the redundant xs calculated. + """ df = self.data.copy() for mat in self.data.columns.get_level_values("MAT").unique(): diff --git a/sandy/pert.py b/sandy/pert.py index 9dde8232..914f0249 100644 --- a/sandy/pert.py +++ b/sandy/pert.py @@ -12,6 +12,17 @@ ] +col = pd.MultiIndex.from_arrays([[2631, 2631], [5, 2]], names=('MAT', 'MT')) +pert = pd.DataFrame([[1, 1.05], [1.05, 1]], + index=pd.IntervalIndex.from_breaks(pd.Index([10, 100]) + .insert(0, 0)), + columns=col) + +pert_ind = pd.Series([1, 1.05], index=pd.IntervalIndex + .from_breaks(pd.Index([10, 100]) + .insert(0, 0))) + + class Pert(): """ Container for tabulated multigroup perturbation coefficients. @@ -42,8 +53,8 @@ class Pert(): def __repr__(self): return self.data.__repr__() - def __init__(self, series, **kwargs): - self.data = pd.Series(series, dtype=float, **kwargs) + def __init__(self, df, **kwargs): + self.data = pd.DataFrame(df, dtype=float, **kwargs) @property def data(self): @@ -93,8 +104,23 @@ def right(self): ------- `pandas.Series` perturbations + + Examples + -------- + >>> sandy.Pert(pert_ind).right + 0 + 1.00000e+01 1.00000e+00 + 1.00000e+02 1.05000e+00 + + >>> sandy.Pert(pert).right + MAT 2631 + MT 5 2 + 1.00000e+01 1.00000e+00 1.05000e+00 + 1.00000e+02 1.05000e+00 1.00000e+00 """ - return pd.Series(self.data.values, index=self.data.index.right) + return pd.DataFrame(self.data.values, + index=self.data.index.right, + columns=self.data.columns) @property def left(self): @@ -107,8 +133,23 @@ def left(self): ------- `pandas.Series` perturbations + + Examples + -------- + >>> sandy.Pert(pert_ind).left + 0 + 0.00000e+00 1.00000e+00 + 1.00000e+01 1.05000e+00 + + >>> sandy.Pert(pert).left + MAT 2631 + MT 5 2 + 0.00000e+00 1.00000e+00 1.05000e+00 + 1.00000e+01 1.05000e+00 1.00000e+00 """ - return pd.Series(self.data.values, index=self.data.index.left) + return pd.DataFrame(self.data.values, + index=self.data.index.left, + columns=self.data.columns) @property def mid(self): @@ -121,10 +162,60 @@ def mid(self): ------- `pandas.Series` perturbations + + Examples + -------- + >>> sandy.Pert(pert_ind).mid + 0 + 5.00000e+00 1.00000e+00 + 5.50000e+01 1.05000e+00 + + >>> sandy.Pert(pert).mid + MAT 2631 + MT 5 2 + 5.00000e+00 1.00000e+00 1.05000e+00 + 5.50000e+01 1.05000e+00 1.00000e+00 """ - return pd.Series(self.data.values, index=self.data.index.mid) + return pd.DataFrame(self.data.values, + index=self.data.index.mid, + columns=self.data.columns) - def reshape(self, eg, inplace=False): + def reorder(self, columns): + """ + Sort the levels in the columns of `Sandy.Pert` + + Parameters + ---------- + columns : `pd.Multiindex` + New columns order. + + Returns + ------- + `sandy.Pert` + Pert order with sort levels. + + Examples + -------- + >>> col = pd.MultiIndex.from_arrays([[5, 2], [2631, 2631]], names=('MT', 'MAT')) + >>> pert = pd.DataFrame([[1, 1.05], [1.05, 1]], index =pd.IntervalIndex.from_breaks(pd.Index([1.94000e+08, 1.96000e+08+1]).insert(0, 0)), columns=col) + >>> col_new = pd.MultiIndex.from_arrays([[2631, 2631], [2, 5]], names=('MAT', 'MT')) + >>> sandy.Pert(pert).reorder(col_new) + MAT 2631 + MT 2 5 + ENERGY + (0.0, 194000000.0] 1.05000e+00 1.00000e+00 + (194000000.0, 196000001.0] 1.00000e+00 1.05000e+00 + """ + data = self.data + for name in columns.names: + if name not in self.data.columns.names: + print("Levels do not match") + return + data_reorder = data.reorder_levels(columns.names, axis=1) + col = columns.intersection(data_reorder.columns) + return sandy.Pert(data_reorder[col]) + + def reshape(self, eg, inplace=False, right_values=1): """ Interpolate perturbation over new energy grid structure using `bfill` method. @@ -150,33 +241,71 @@ def reshape(self, eg, inplace=False): if the given energy grid is not monotonic `value.Error` if negative values are found in the given energy grid + + Examples + -------- + >>> sandy.Pert(pert).reshape([0, 1, 5, 25, 50, 100, 1.95000e+08]).data + MAT 2631 + MT 5 2 + ENERGY + (0.0, 1.0] 1.00000e+00 1.05000e+00 + (1.0, 5.0] 1.00000e+00 1.05000e+00 + (5.0, 10.0] 1.00000e+00 1.05000e+00 + (10.0, 25.0] 1.05000e+00 1.00000e+00 + (25.0, 50.0] 1.05000e+00 1.00000e+00 + (50.0, 100.0] 1.05000e+00 1.00000e+00 + (100.0, 195000000.0] 1.00000e+00 1.00000e+00 + + >>> sandy.Pert(pert_ind).reshape([0, 1, 5, 25, 50, 100, 1.95000e+08]).data + 0 + ENERGY + (0.0, 1.0] 1.00000e+00 + (1.0, 5.0] 1.00000e+00 + (5.0, 10.0] 1.00000e+00 + (10.0, 25.0] 1.05000e+00 + (25.0, 50.0] 1.05000e+00 + (50.0, 100.0] 1.05000e+00 + (100.0, 195000000.0] 1.00000e+00 """ - index = pd.Index(eg) - if not index.is_monotonic_increasing: + eg_ = pd.Index(eg) + df = self.right + if not eg_.is_monotonic_increasing: raise sandy.Error("energy grid is not monotonic increasing") - if (index < 0).any(): + if (eg_ < 0).any(): raise ValueError("found negative values in the energy grid") - enew = self.right.index.union(index).unique().astype(float).values + enew = df.index.union(eg_).unique().astype(float).values # remove zero if any, it will be automatically added by `Pert` enew = enew[enew != 0] # this part prevents errors in "scipy.interp1d" when x.size == 1 - x = self.right.index.values - y = self.right.values - if y.size == 1: - # this must be done after that enew is created - x = np.insert(x, 0, 0) - y = np.insert(y, 0, 0) - pertnew = sandy.shared.reshape_bfill( - x, - y, + name = df.index.name + if df.shape[1] == 1: + index = df.index.values + values = df.values + if values.size == 1: + # this must be done after that enew is created + index = np.insert(index, 0, 0) + values = np.insert(values, 0, 0) + pertnew = sandy.shared.reshape_bfill( + index, + values, enew, - left_values=1, - right_values=1, + left_values="first", + right_values=right_values, ) - series = pd.Series(pertnew, index=enew) + data = pd.DataFrame(pertnew, index=enew, columns=df.columns) + else: + data = df.apply(lambda x: sandy.shared.reshape_bfill( + x.index.values, + x.values, + enew, + left_values="first", + right_values=right_values, + )).set_index(enew).rename_axis(name) + data = data.loc[(data.index >= eg_.min()) & + (data.index <= eg_.max()), :] if not inplace: - return Pert(series) - self.data = series + return self.__class__(data) + self.data = data def truncate(self, low=0.0, high=2.0): """ @@ -196,6 +325,23 @@ def truncate(self, low=0.0, high=2.0): ------- `sandy.Pert` perturbation instance with truncated values. + + Examples + -------- + >>> pert = pd.DataFrame([[-1, 2.55], [2.55, -1]], index =pd.IntervalIndex.from_breaks(pd.Index([10, 100]).insert(0, 0)), columns=col) + >>> sandy.Pert(pert).truncate().data + MAT 2631 + MT 5 2 + ENERGY + (0.0, 10.0] 0.00000e+00 2.00000e+00 + (10.0, 100.0] 2.00000e+00 0.00000e+00 + + >>> pert_ind = pd.Series([-1, 2.05], index=pd.IntervalIndex.from_breaks(pd.Index([10, 100]).insert(0, 0))) + >>> sandy.Pert(pert_ind).truncate().data + 0 + ENERGY + (0.0, 10.0] 0.00000e+00 + (10.0, 100.0] 2.00000e+00 """ data = self.data.copy() data = data.where(data <= high, high).where(data >= low, low) @@ -222,6 +368,23 @@ def recenter(self, low=0.0, high=2.0, value=1.): ------- `sandy.Pert` perturbation instance with truncated values. + + Examples + -------- + >>> pert = pd.DataFrame([[-1, 2.55], [2.55, -1]], index =pd.IntervalIndex.from_breaks(pd.Index([10, 100]).insert(0, 0)), columns=col) + >>> sandy.Pert(pert).recenter().data + MAT 2631 + MT 5 2 + ENERGY + (0.0, 10.0] 1.00000e+00 1.00000e+00 + (10.0, 100.0] 1.00000e+00 1.00000e+00 + + >>> pert_ind = pd.Series([-1, 2.05], index=pd.IntervalIndex.from_breaks(pd.Index([10, 100]).insert(0, 0))) + >>> sandy.Pert(pert_ind).recenter().data + 0 + ENERGY + (0.0, 10.0] 1.00000e+00 + (10.0, 100.0] 1.00000e+00 """ data = self.data.copy() data = data.where(data <= high, value).where(data >= low, value) @@ -259,7 +422,7 @@ def from_bin(cls, elow, ehigh, coeff): Parameters ---------- - elow : TYPE + elow : `float` lower energy boundary in eV. ehigh : `float` upper energy boundary in eV. @@ -274,10 +437,10 @@ def from_bin(cls, elow, ehigh, coeff): Examples -------- >>> Pert.from_bin(1e-5, 1e-4, 0.05) - ENERGY - (0.0, 1e-05] 1.00000e+00 - (1e-05, 0.0001] 1.05000e+00 - dtype: float64 + 0 + ENERGY + (0.0, 1e-05] 1.00000e+00 + (1e-05, 0.0001] 1.05000e+00 """ pert = cls([1, 1 + coeff], index=[elow, ehigh]) return pert diff --git a/sandy/pfns.py b/sandy/pfns.py index 5b85195b..f8e22366 100644 --- a/sandy/pfns.py +++ b/sandy/pfns.py @@ -186,6 +186,64 @@ def get_table(self, mat, mt, k): .sort_index(axis="index") \ .sort_index(axis="columns") + def reshape(self, enew): + """ + Reshape the outgoing energy. + + Parameters + ---------- + enew : 1D iterable or `float` + new energy grid. + + Returns + ------- + `sandy.Edistr` + Energy distribution instance over new grid. + + Examples + -------- + >>> orig = Edistr(minimal_edistrtest) + >>> orig.reshape(1.5) + MAT MT K EIN EOUT VALUE + 0 9437 18 0 1.00000e+00 1.00000e-05 4.00000e-01 + 1 9437 18 0 1.00000e+00 1.50000e+00 4.00000e-01 + 2 9437 18 0 1.00000e+00 2.00000e+07 6.00000e-01 + 3 9437 18 0 2.00000e+00 1.00000e-04 2.00000e-01 + 4 9437 18 0 2.00000e+00 1.00000e+00 7.00000e-01 + 5 9437 18 0 2.00000e+00 1.50000e+00 7.00000e-01 + 6 9437 18 0 2.00000e+00 1.00000e+07 1.00000e-01 + + >>> orig.reshape([1.5e-5, 1.5, 15e6]) + MAT MT K EIN EOUT VALUE + 0 9437 18 0 1.00000e+00 1.00000e-05 4.00000e-01 + 1 9437 18 0 1.00000e+00 1.50000e-05 4.00000e-01 + 2 9437 18 0 1.00000e+00 1.50000e+00 4.00000e-01 + 3 9437 18 0 1.00000e+00 1.50000e+07 5.50000e-01 + 4 9437 18 0 1.00000e+00 2.00000e+07 6.00000e-01 + 5 9437 18 0 2.00000e+00 1.00000e-04 2.00000e-01 + 6 9437 18 0 2.00000e+00 1.00000e+00 7.00000e-01 + 7 9437 18 0 2.00000e+00 1.50000e+00 7.00000e-01 + 8 9437 18 0 2.00000e+00 1.00000e+07 1.00000e-01 + """ + enew_ = np.array(enew) if isinstance(enew, np.ndarray) else enew + + def foo(df, enew): + df_ = df.loc[:, ["EOUT", "VALUE"]].set_index("EOUT") + enew_ = np.union1d(df_.index.values, enew) + enew_ = enew_[(enew_ >= df_.index.min()) & (enew_ <= df_.index.max())] + new_edistr = sandy.shared.reshape_differential( + df_.index.values, + df_.values, + enew_, + ) + new_edistr = pd.DataFrame(new_edistr, index=enew_, + columns=df_.columns) + new_edistr.index.name = 'EOUT' + return new_edistr + edistr_reshape = self.data.groupby(['MAT', 'MT', 'K', 'EIN'])\ + .apply(foo, enew_).reset_index() + return self.__class__(edistr_reshape) + def add_energy_point(self, mat, mt, k, enew): """ Add outgoing energy distribution at one additional incident energy by @@ -378,7 +436,7 @@ def normalize(self): df = pd.concat(out) return self.__class__(df) - def custom_perturbation(self, pert, mat, mt, k, ein_low, ein_high): + def custom_perturbation(self, pert): """ Given a peruration object (fractions), a MAT number, a MT number, a subsection number, a lower and an upper incoming energy bound, @@ -387,18 +445,8 @@ def custom_perturbation(self, pert, mat, mt, k, ein_low, ein_high): Parameters ---------- - pert : `sandy.Pert` + pert : `sandy.Pert` or `pd.Series` perturbation object. - mat : `int` - MAT number. - mt : `int` - MT number. - k : `int` - subsection. - ein_low : TYPE - lower energy boundary in eV. - ein_high : TYPE - upper energy boundary in eV. Returns ------- @@ -417,31 +465,48 @@ def custom_perturbation(self, pert, mat, mt, k, ein_low, ein_high): Examples -------- >>> orig = Edistr(minimal_edistrtest) - >>> pert = sandy.Pert([1.3], index=[1e-3]) - >>> orig.custom_perturbation(pert, 9437, 18, 0, 1.5, 2.5) + >>> pert = pd.Series([1.3], index=[1e-3]) + >>> columns = pd.MultiIndex.from_product([[9437], [18], [0], [1.5], [2.5]], names=['MAT', 'MT', 'K', 'ELO', 'EHI']) + >>> df = pd.DataFrame(pert.values, index=pert.index, columns=columns) + >>> pert_ = sandy.Pert(df) + >>> orig.custom_perturbation(pert_) MAT MT K EIN EOUT VALUE - 0 9437 18 0 1.00000e+00 1.00000e-05 4.00000e-01 - 1 9437 18 0 1.00000e+00 2.00000e+07 6.00000e-01 - 2 9437 18 0 2.00000e+00 1.00000e-04 2.60000e-01 - 3 9437 18 0 2.00000e+00 1.00000e+00 7.00000e-01 - 4 9437 18 0 2.00000e+00 1.00000e+07 1.00000e-01 + 0 9437 18 0 1.00000e+00 1.00000e-05 4.00000e-08 + 1 9437 18 0 1.00000e+00 2.00000e+07 6.00000e-08 + 2 9437 18 0 2.00000e+00 1.00000e-04 3.75000e-07 + 3 9437 18 0 2.00000e+00 1.00000e+00 1.75000e-07 + 4 9437 18 0 2.00000e+00 1.00000e+07 2.50000e-08 """ - data = self.data.copy() - condition = (data.MT == mt) &\ - (data.MAT == mat) &\ - (data.K == k) &\ - (data.EIN < ein_high) &\ - (data.EIN >= ein_low) - dfs = [] - dfs.append(data[~condition]) - for ein, df in data[condition].groupby("EIN"): - series = pert.reshape(df.EOUT).data.loc[df.EOUT] - # truncate extremes and replace them with boundaries - px = sandy.Pert(series).truncate() - df.VALUE *= px.data.values - dfs.append(df) - out = pd.concat(dfs) - return self.__class__(out) + pert_ = sandy.Pert(pert) if not isinstance(pert, sandy.Pert) else pert + + # New energy grid: + energy_grid = self.data.loc[:, 'EOUT'].unique() + enew = np.union1d(energy_grid, pert.right.index) + enew = enew[(enew <= energy_grid.max()) & (enew >= energy_grid.min())] + + # Reshape to new energy grid and columns estructure: + u_pert = pert.reshape(enew, right_values=0).right + + # Apply the perturbation: + def foo(df, pert): + ein = df.loc[:, 'EIN'].unique()[0] + mat = df.loc[:, 'MAT'].unique()[0] + mt = df.loc[:, 'MT'].unique()[0] + col = pert.columns + mask = (mat == col.get_level_values('MAT')) & \ + (mt == col.get_level_values('MT')) & \ + (ein >= col.get_level_values('ELO')) & \ + (ein <= col.get_level_values('EHI')) + pert_ = pert.loc[:, mask] + if not pert_.empty: + pert_ = pert_.iloc[:, [0]]\ + .reindex(index=df.loc[:, "EOUT"].values)\ + .values.flatten() + df["VALUE"] = df['VALUE'].values + pert_ + return df + pert_edistr = self.data.groupby(['MAT', 'MT', 'K', 'EIN'])\ + .apply(foo, u_pert) + return self.__class__(pert_edistr).normalize() def _perturb(self, pert, method=2, normalize=True, **kwargs): """Perturb energy distributions given a set of perturbations.