From f9418ad8eb75ad6ea9dbb23d16de30a9f40fa752 Mon Sep 17 00:00:00 2001 From: Ethan Peterson Date: Tue, 28 Feb 2023 14:20:06 -0500 Subject: [PATCH 1/4] reworking Be assembly example --- fns_be/fns_be.py | 388 +++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 388 insertions(+) create mode 100644 fns_be/fns_be.py diff --git a/fns_be/fns_be.py b/fns_be/fns_be.py new file mode 100644 index 00000000..daf525eb --- /dev/null +++ b/fns_be/fns_be.py @@ -0,0 +1,388 @@ +#!/usr/bin/env/ python3 +import numpy as np +import argparse + + +import openmc + + +def _parse_args(): + # Create argument parser + parser = argparse.ArgumentParser() + parser.add_argument('-b', '--batches', type=int, default=100) + parser.add_argument('-p', '--particles', type=int, default=10000) + group = parser.add_mutually_exclusive_group(required=True) + group.add_argument('-n', '--neutron_spectrum', action='store_true') + group.add_argument('-g', '--gamma_heating', action='store_true') + group.add_argument('-r', '--reaction_rates', action='store_true') + parser.add_argument('-s', '--threads', type=int) + + # Parse and return commandline arguments. + return parser.parse_args() + + +def main(): + ############################################################################ + # Analysis of Beryllium Experiment + # + ############################################################################ + args = vars(_parse_args()) + + # Instantiate Model object + model = openmc.Model() + + # Build materials + air = openmc.Material(name='air') + air.add_nuclide('N14', 0.7886608412924202, 'ao') + air.add_nuclide('O16', 0.21133915870757974, 'ao') + + beryllium = openmc.Material(name='Beryllium') + beryllium.add_nuclide('Be9', 0.9927289710133025, 'ao') + beryllium.add_nuclide('C12', 0.0006299237839521456, 'ao') + beryllium.add_nuclide('O16', 0.0040693555162183695, 'ao') + beryllium.add_nuclide('Al27', 0.0023701485875583395, 'ao') + beryllium.add_nuclide('Fe54', 1.1783584234715633e-05, 'ao') + beryllium.add_nuclide('Fe56', 0.00018497707234766441, 'ao') + beryllium.add_nuclide('Fe57', 4.271927287144984e-06, 'ao') + beryllium.add_nuclide('Fe58', 5.685150990914985e-07, 'ao') + + tube = openmc.Material(name='drawer-tube') + tube.add_nuclide('Mn55', 0.00997121925432881, 'ao') + tube.add_nuclide('Cr50', 0.008470987544935695, 'ao') + tube.add_nuclide('Cr52', 0.1633545628084274, 'ao') + tube.add_nuclide('Cr53', 0.018523096125301272, 'ao') + tube.add_nuclide('Cr54', 0.004610790689015631, 'ao') + tube.add_nuclide('Fe54', 0.041430336359409875, 'ao') + tube.add_nuclide('Fe56', 0.6503676787547126, 'ao') + tube.add_nuclide('Fe57', 0.015019825961606423, 'ao') + tube.add_nuclide('Fe58', 0.001998863105791888, 'ao') + tube.add_nuclide('Ni58', 0.05871820932193519, 'ao') + tube.add_nuclide('Ni60', 0.022618029628936442, 'ao') + tube.add_nuclide('Ni61', 0.0009831938364803666, 'ao') + tube.add_nuclide('Ni62', 0.0031349384315041144, 'ao') + tube.add_nuclide('Ni64', 0.0007982681776143339, 'ao') + + if args['reaction_rates']: + # Add reaction rate elements to materials in small quantities + air.add_element('Al', 1e-8, 'ao') + air.add_element('Nb', 1e-8, 'ao') + air.add_element('In', 1e-8, 'ao') + air.add_element('Au', 1e-8, 'ao') + air.add_element('Fe', 1e-8, 'ao') + air.add_element('Zr', 1e-8, 'ao') + air.add_element('Ni', 1e-8, 'ao') + air.add_element('Ti', 1e-8, 'ao') + air.add_nuclide('U235', 1e-8, 'ao') + air.add_nuclide('Li6', 1e-8, 'ao') + air.add_nuclide('Be9', 1e-8, 'ao') + beryllium.add_element('Nb', 1e-8, 'ao') + beryllium.add_element('In', 1e-8, 'ao') + beryllium.add_element('Au', 1e-8, 'ao') + beryllium.add_element('Zr', 1e-8, 'ao') + beryllium.add_element('Ni', 1e-8, 'ao') + beryllium.add_element('Ti', 1e-8, 'ao') + beryllium.add_nuclide('U235', 1e-8, 'ao') + beryllium.add_nuclide('Li6', 1e-8, 'ao') + tube.add_element('Al', 1e-8, 'ao') + tube.add_element('Nb', 1e-8, 'ao') + tube.add_element('In', 1e-8, 'ao') + tube.add_element('Au', 1e-8, 'ao') + tube.add_element('Zr', 1e-8, 'ao') + tube.add_element('Ti', 1e-8, 'ao') + tube.add_nuclide('U235', 1e-8, 'ao') + tube.add_nuclide('Li6', 1e-8, 'ao') + + air.set_density('atom/b-cm', 4.921e-05) + beryllium.set_density('atom/b-cm', 0.122149) + tube.set_density('atom/b-cm', 0.0200898) + + ############################################################################ + # Build Geometry + + # Surfaces + s1 = openmc.ZPlane(-10.0, boundary_type='vacuum') + s2 = openmc.ZPlane(20.0) + s3 = openmc.ZPlane(65.54) + s4 = openmc.ZPlane(70.00) + s5 = openmc.ZCylinder(r=31.5) + s6 = openmc.ZCylinder(r=35.0, boundary_type='vacuum') + s7 = openmc.XPlane(-2.46) + s8 = openmc.XPlane(2.46) + s9 = openmc.XPlane(-2.54) + s10 = openmc.XPlane(2.54) + s11 = openmc.YPlane(-2.46) + s12 = openmc.YPlane(2.46) + s13 = openmc.YPlane(-2.54) + s14 = openmc.YPlane(2.54) + s15 = openmc.ZPlane(19.9) + s16 = openmc.ZPlane(21.0) + s17 = openmc.ZPlane(22.0) + s18 = openmc.ZPlane(23.4) + s19 = openmc.ZPlane(24.0) + s20 = openmc.ZPlane(25.0) + s21 = openmc.ZPlane(27.0) + s22 = openmc.ZPlane(29.0) + s23 = openmc.ZPlane(31.0) + s24 = openmc.ZPlane(33.0) + s25 = openmc.ZPlane(35.0) + s26 = openmc.ZPlane(37.0) + s27 = openmc.ZPlane(39.0) + s28 = openmc.ZPlane(41.0) + s29 = openmc.ZPlane(43.0) + s30 = openmc.ZPlane(45.0) + s31 = openmc.ZPlane(48.0) + s32 = openmc.ZPlane(51.0) + s33 = openmc.ZPlane(54.0) + s34 = openmc.ZPlane(57.0) + s35 = openmc.ZPlane(60.0) + s36 = openmc.ZPlane(63.0) + s37 = openmc.ZPlane(66.0) + s38 = openmc.Sphere(r=0.5) + s39 = openmc.ZPlane(-15.0) + s40 = openmc.ZPlane(75.0, boundary_type='vacuum') + + # Helper Regions + xybox = +s7 & -s8 & +s11 & -s12 + zcyl = +s2 & -s3 & -s5 + region1 = +s1 & -s2 & -s5 & +s38 & ~(+s15 & -s2 & xybox) + region2 = +s3 & -s4 & -s5 & ~(+s3 & -s37 & xybox) + region3 = (zcyl & -s9) | (zcyl & +s10) | (zcyl & -s13) | (zcyl & +s14) + region4 = (+s2 & -s3 & +s9 & -s7 & +s13 & -s14) | \ + (+s2 & -s3 & +s8 & -s10 & +s13 & -s14) | \ + (+s2 & -s3 & +s7 & -s8 & +s13 & -s11) | \ + (+s2 & -s3 & +s7 & -s8 & +s12 & -s14) + region30 = (+s39 & -s1 & -s6) | (+s1 & -s4 & -s6 & +s5) | (+s4 & -s40 & -s6) + + # Cells + cell1 = openmc.Cell(region=region1, fill=air, name='cell 1') + cell2 = openmc.Cell(region=region2, fill=air, name='cell 2') + cell3 = openmc.Cell(region=region3, fill=beryllium, name='cell 3') + cell4 = openmc.Cell(region=region4, fill=tube, name='cell 4') + cell5 = openmc.Cell(region=+s15 & -s2 & xybox, fill=air, name='cell 5') + cell6 = openmc.Cell(region=+s2 & -s16 & xybox, fill=beryllium, name='cell 6') + cell7 = openmc.Cell(region=+s16 & -s17 & xybox, fill=beryllium, name='cell 7') + cell8 = openmc.Cell(region=+s17 & -s18 & xybox, fill=beryllium, name='cell 8') + cell9 = openmc.Cell(region=+s18 & -s19 & xybox, fill=beryllium, name='cell 9') + cell10 = openmc.Cell(region=+s19 & -s20 & xybox, fill=beryllium, name='cell 10') + cell11 = openmc.Cell(region=+s20 & -s21 & xybox, fill=beryllium, name='cell 11') + cell12 = openmc.Cell(region=+s21 & -s22 & xybox, fill=beryllium, name='cell 12') + cell13 = openmc.Cell(region=+s22 & -s23 & xybox, fill=beryllium, name='cell 13') + cell14 = openmc.Cell(region=+s23 & -s24 & xybox, fill=beryllium, name='cell 14') + cell15 = openmc.Cell(region=+s24 & -s25 & xybox, fill=beryllium, name='cell 15') + cell16 = openmc.Cell(region=+s25 & -s26 & xybox, fill=beryllium, name='cell 16') + cell17 = openmc.Cell(region=+s26 & -s27 & xybox, fill=beryllium, name='cell 17') + cell18 = openmc.Cell(region=+s27 & -s28 & xybox, fill=beryllium, name='cell 18') + cell19 = openmc.Cell(region=+s28 & -s29 & xybox, fill=beryllium, name='cell 19') + cell20 = openmc.Cell(region=+s29 & -s30 & xybox, fill=beryllium, name='cell 20') + cell21 = openmc.Cell(region=+s30 & -s31 & xybox, fill=beryllium, name='cell 21') + cell22 = openmc.Cell(region=+s31 & -s32 & xybox, fill=beryllium, name='cell 22') + cell23 = openmc.Cell(region=+s32 & -s33 & xybox, fill=beryllium, name='cell 23') + cell24 = openmc.Cell(region=+s33 & -s34 & xybox, fill=beryllium, name='cell 24') + cell25 = openmc.Cell(region=+s34 & -s35 & xybox, fill=beryllium, name='cell 25') + cell26 = openmc.Cell(region=+s35 & -s36 & xybox, fill=beryllium, name='cell 26') + cell27 = openmc.Cell(region=+s36 & -s3 & xybox, fill=beryllium, name='cell 27') + cell28 = openmc.Cell(region=+s3 & -s37 & xybox, fill=air, name='cell 28') + cell29 = openmc.Cell(region=-s38, fill=air, name='cell 29') + cell30 = openmc.Cell(region=region30, name='cell 30') + + cells = [cell1, cell2, cell3, cell4, cell5, cell6, cell7, cell8, cell9, cell10, + cell11, cell12, cell13, cell14, cell15, cell16, cell17, cell18, + cell19, cell20, cell21, cell22, cell23, cell24, cell25, cell26, + cell27, cell28, cell29, cell30 + ] + + # Add Geometry to Model + model.geometry = openmc.Geometry(cells) + + ############################################################################### + # Define problem settings + + # Indicate how many particles to run + settings = openmc.Settings(run_mode='fixed source') + settings.batches = args['batches'] + settings.particles = args['particles'] + + # Create source distribution + # upper neutron energy bins + x = 1e6*np.array([1.0010e-11, 3.2241e-07, + 5.3156e-07, 8.7640e-07, 1.4449e-06, 2.3823e-06, 3.9278e-06, + 6.4758e-06, 1.0677e-05, 1.7603e-05, 2.9023e-05, 4.7850e-05, + 7.8891e-05, 1.3007e-04, 2.1445e-04, 3.5357e-04, 5.8293e-04, + 9.6110e-04, 1.2341e-03, 1.5846e-03, 2.0346e-03, 2.6125e-03, + 3.3546e-03, 4.3073e-03, 5.5307e-03, 7.1016e-03, 9.1186e-03, + 1.1709e-02, 1.5034e-02, 1.9304e-02, 2.1874e-02, 2.4787e-02, + 2.8087e-02, 3.1827e-02, 3.6065e-02, 4.0867e-02, 4.6308e-02, + 5.2474e-02, 5.9461e-02, 6.7378e-02, 7.6349e-02, 8.6515e-02, + 9.8035e-02, 1.1109e-01, 1.2588e-01, 1.4264e-01, 1.6163e-01, + 1.8315e-01, 2.0754e-01, 2.3517e-01, 2.6649e-01, 3.0197e-01, + 3.4217e-01, 3.8774e-01, 4.3936e-01, 4.9786e-01, 5.6415e-01, + 6.3927e-01, 7.2438e-01, 8.2084e-01, 9.3013e-01, 1.0540e+00, + 1.1943e+00, 1.3533e+00, 1.5335e+00, 1.7377e+00, 1.8498e+00, + 1.9691e+00, 2.0961e+00, 2.2313e+00, 2.3752e+00, 2.5284e+00, + 2.6914e+00, 2.8650e+00, 3.0498e+00, 3.2465e+00, 3.4559e+00, + 3.6787e+00, 3.9160e+00, 4.1686e+00, 4.4374e+00, 4.7236e+00, + 5.0282e+00, 5.3525e+00, 5.6978e+00, 6.0652e+00, 6.4564e+00, + 6.8728e+00, 7.3161e+00, 7.7879e+00, 8.2902e+00, 8.8249e+00, + 9.3940e+00, 9.9999e+00, 1.0157e+01, 1.0317e+01, 1.0480e+01, + 1.0645e+01, 1.0812e+01, 1.0983e+01, 1.1156e+01, 1.1331e+01, + 1.1510e+01, 1.1691e+01, 1.1875e+01, 1.2062e+01, 1.2252e+01, + 1.2445e+01, 1.2641e+01, 1.2840e+01, 1.3042e+01, 1.3248e+01, + 1.3456e+01, 1.3668e+01, 1.3883e+01, 1.4102e+01, 1.4324e+01, + 1.4550e+01, 1.4779e+01, 1.5012e+01, 1.5248e+01, 1.5488e+01]) + + # bin probabilities + p = np.array([0.0, 1.5142e-07, + 2.2732e-09, 4.2225e-09, 7.4848e-09, 1.4264e-08, 8.3975e-08, + 1.8398e-07, 2.2450e-07, 1.3922e-07, 1.6817e-07, 2.9754e-07, + 3.8068e-06, 3.0541e-06, 2.2612e-06, 6.9372e-06, 7.2049e-06, + 8.7622e-06, 7.8013e-06, 1.4320e-05, 1.1820e-05, 1.6544e-05, + 1.4791e-05, 1.7624e-05, 2.8404e-05, 2.4899e-05, 3.7633e-05, + 4.4237e-05, 4.6320e-05, 6.1572e-05, 3.7185e-05, 5.3362e-05, + 4.8831e-05, 5.0292e-05, 5.7202e-05, 6.9230e-05, 8.0602e-05, + 8.3190e-05, 9.7450e-05, 1.0531e-04, 1.2632e-04, 1.4874e-04, + 1.7906e-04, 3.7225e-04, 4.9933e-04, 5.3824e-04, 6.0762e-04, + 7.0593e-04, 8.0965e-04, 9.5392e-04, 1.0785e-03, 1.2232e-03, + 1.3867e-03, 1.5803e-03, 1.6473e-03, 1.8238e-03, 2.0605e-03, + 2.2042e-03, 2.3040e-03, 2.5211e-03, 2.5709e-03, 2.5872e-03, + 2.5765e-03, 2.7699e-03, 2.8528e-03, 2.5945e-03, 1.3898e-03, + 1.4298e-03, 1.3270e-03, 1.3489e-03, 1.3820e-03, 1.4312e-03, + 1.3760e-03, 1.4329e-03, 1.4558e-03, 1.3518e-03, 1.4053e-03, + 1.2861e-03, 1.2741e-03, 1.1711e-03, 1.1937e-03, 1.0563e-03, + 1.0018e-03, 8.8451e-04, 7.9827e-04, 7.9293e-04, 7.5872e-04, + 6.9228e-04, 6.2956e-04, 5.1710e-04, 5.0750e-04, 5.1007e-04, + 4.1280e-04, 3.5649e-04, 9.0768e-05, 8.2287e-05, 9.2862e-05, + 9.1407e-05, 9.3708e-05, 7.9567e-05, 8.8737e-05, 8.7841e-05, + 1.1227e-04, 1.6798e-04, 1.5985e-04, 1.6563e-04, 2.1025e-04, + 4.1363e-04, 7.4899e-04, 7.8183e-04, 5.1771e-04, 4.5938e-04, + 4.6458e-04, 9.1020e-04, 2.6083e-03, 9.5007e-04, 5.1474e-03, + 3.0897e-02, 2.3565e-01, 4.0901e-01, 2.2296e-01, 1.4419e-01]) + + source = openmc.Source() + source.particle = 'neutron' + source.space = openmc.stats.Point() + source.angle = openmc.stats.Isotropic() + source.energy = openmc.stats.Tabular(x, p, interpolation='histogram') + settings.source = source + settings.output = {'tallies': False} + + ############################################################################### + # Build tallies based on the desired run + cell_filter = openmc.CellFilter(cells[4:26]) + + if args['neutron_spectrum']: + # neutron energy spectrum tally + energy_n = openmc.EnergyFilter(x) + particle_n = openmc.ParticleFilter(['neutron']) + tally_n = openmc.Tally(name="Neutron Spectrum") + tally_n.filters = [cell_filter, energy_n, particle_n] + tally_n.scores = ['flux'] + model.tallies = openmc.Tallies([tally_n]) + cwd = 'neutron_spectrum' + + elif args['gamma_heating']: + # Tallies gamma heating + settings.survival_biasing = True + settings.photon_transport = True + settings.electron_treatment = 'ttb' + particle_p = openmc.ParticleFilter(['neutron','photon','electron','positron']) + tally_p = openmc.Tally(name="Gamma heating") + tally_p.filters = [cell_filter, particle_p] + tally_p.scores = ['heating'] + model.tallies = openmc.Tallies([tally_p]) + cwd = 'gamma_heating' + + elif args['reaction_rates']: + particle_n = openmc.ParticleFilter(['neutron']) + + tally_1 = openmc.Tally(name="Au-197 (n,gamma)") + tally_1.nuclides = ['Au197'] + tally_1.filters = [cell_filter, particle_n] + tally_1.scores = ['(n,gamma)'] + + tally_2 = openmc.Tally(name="In115 (n,n')") + tally_2.nuclides = ['In115'] + tally_2.filters = [cell_filter, particle_n] + tally_2.scores = ['(n,nc)'] + + tally_3 = openmc.Tally(name="Al-27 (n,alpha)") + tally_3.nuclides = ['Al27'] + tally_3.filters = [cell_filter, particle_n] + tally_3.scores = ['(n,a)'] + + tally_4 = openmc.Tally(name="Nb-93 (n,2n)") + tally_4.nuclides = ['Nb93'] + tally_4.filters = [cell_filter, particle_n] + tally_4.scores = ['(n,2n)'] + + tally_5 = openmc.Tally(name="Ni-58 (n,2n)") + tally_5.nuclides = ['Ni58'] + tally_5.filters = [cell_filter, particle_n] + tally_5.scores = ['(n,2n)'] + + tally_6 = openmc.Tally(name="Ti-48(n,p)") + tally_6.nuclides = ['Ti48'] + tally_6.filters = [cell_filter, particle_n] + tally_6.scores = ['(n,p)'] + + tally_7 = openmc.Tally(name="Ti-49 (n,np)") + tally_7.nuclides = ['Ti49'] + tally_7.filters = [cell_filter, particle_n] + tally_7.scores = ['(n,np)'] + + tally_8 = openmc.Tally(name="Ti-47 (n,p)") + tally_8.nuclides = ['Ti47'] + tally_8.filters = [cell_filter, particle_n] + tally_8.scores = ['(n,p)'] + + tally_9 = openmc.Tally(name="Ti-48 (n,np)") + tally_9.nuclides = ['Ti48'] + tally_9.filters = [cell_filter, particle_n] + tally_9.scores = ['(n,np)'] + + tally_10 = openmc.Tally(name="Fe-56 (n,p)") + tally_10.nuclides = ['Fe56'] + tally_10.filters = [cell_filter, particle_n] + tally_10.scores = ['(n,p)'] + + tally_11 = openmc.Tally(name="Zr-90 (n,2n)") + tally_11.nuclides = ['Zr90'] + tally_11.filters = [cell_filter, particle_n] + tally_11.scores = ['(n,2n)'] + + tally_12 = openmc.Tally(name="Ni-58 (n,p)") + tally_12.nuclides = ['Ni58'] + tally_12.filters = [cell_filter, particle_n] + tally_12.scores = ['(n,p)'] + + tally_13 = openmc.Tally(name="Li-6 (n,tritium)") + tally_13.nuclides = ['Li6'] + tally_13.filters = [cell_filter, particle_n] + tally_13.scores = ['(n,Xt)'] + + tally_14 = openmc.Tally(name="Tritium production") + tally_14.filters = [cell_filter, particle_n] + tally_14.scores = ['H3-production'] + + tally_15 = openmc.Tally(name="U-235 (n,f)") + tally_15.nuclides = ['U235'] + tally_15.filters = [cell_filter, particle_n] + tally_15.scores = ['fission'] + + tally_16 = openmc.Tally(name="Be-9 (n,tritium)") + tally_16.nuclides = ['Be9'] + tally_16.filters = [cell_filter, particle_n] + tally_16.scores = ['(n,Xt)'] + + model.tallies = openmc.Tallies([tally_1, tally_2, tally_3, tally_4, + tally_5, tally_6, tally_7, tally_8, + tally_9, tally_10, tally_11, tally_12, + tally_13, tally_14, tally_15, tally_16]) + cwd = 'reaction_rates' + + model.settings = settings + + model.run(cwd=cwd, threads=args['threads']) + + +if __name__ == "__main__": + main() From 81a11349805a79704416fed77f3c2766d785f51f Mon Sep 17 00:00:00 2001 From: Ethan Peterson Date: Wed, 1 Mar 2023 11:39:26 -0500 Subject: [PATCH 2/4] Apply suggestions from @pshriwise code review Co-authored-by: Patrick Shriwise --- fns_be/fns_be.py | 17 ++++++----------- 1 file changed, 6 insertions(+), 11 deletions(-) diff --git a/fns_be/fns_be.py b/fns_be/fns_be.py index daf525eb..70ede82c 100644 --- a/fns_be/fns_be.py +++ b/fns_be/fns_be.py @@ -2,7 +2,6 @@ import numpy as np import argparse - import openmc @@ -22,10 +21,7 @@ def _parse_args(): def main(): - ############################################################################ - # Analysis of Beryllium Experiment - # - ############################################################################ + """Analysis of Beryllium Experiment""" args = vars(_parse_args()) # Instantiate Model object @@ -147,10 +143,10 @@ def main(): region1 = +s1 & -s2 & -s5 & +s38 & ~(+s15 & -s2 & xybox) region2 = +s3 & -s4 & -s5 & ~(+s3 & -s37 & xybox) region3 = (zcyl & -s9) | (zcyl & +s10) | (zcyl & -s13) | (zcyl & +s14) - region4 = (+s2 & -s3 & +s9 & -s7 & +s13 & -s14) | \ - (+s2 & -s3 & +s8 & -s10 & +s13 & -s14) | \ - (+s2 & -s3 & +s7 & -s8 & +s13 & -s11) | \ - (+s2 & -s3 & +s7 & -s8 & +s12 & -s14) + region4 = (+s2 & -s3 & +s9 & -s7 & +s13 & -s14) + region4 |= (+s2 & -s3 & +s8 & -s10 & +s13 & -s14) + region4 |= (+s2 & -s3 & +s7 & -s8 & +s13 & -s11) + region4 |= (+s2 & -s3 & +s7 & -s8 & +s12 & -s14) region30 = (+s39 & -s1 & -s6) | (+s1 & -s4 & -s6 & +s5) | (+s4 & -s40 & -s6) # Cells @@ -278,7 +274,6 @@ def main(): tally_n.scores = ['flux'] model.tallies = openmc.Tallies([tally_n]) cwd = 'neutron_spectrum' - elif args['gamma_heating']: # Tallies gamma heating settings.survival_biasing = True @@ -381,7 +376,7 @@ def main(): model.settings = settings - model.run(cwd=cwd, threads=args['threads']) + return model.run(cwd=cwd, threads=args['threads']) if __name__ == "__main__": From 0b1e58b50e3278489cc629726da352e64153aec2 Mon Sep 17 00:00:00 2001 From: Ethan Peterson Date: Wed, 1 Mar 2023 15:55:08 -0500 Subject: [PATCH 3/4] addressing more review comments --- fns_be/fns_be.py | 163 ++++++++++++++++++++++++++--------------------- 1 file changed, 90 insertions(+), 73 deletions(-) mode change 100644 => 100755 fns_be/fns_be.py diff --git a/fns_be/fns_be.py b/fns_be/fns_be.py old mode 100644 new mode 100755 index 70ede82c..00adf35d --- a/fns_be/fns_be.py +++ b/fns_be/fns_be.py @@ -1,4 +1,4 @@ -#!/usr/bin/env/ python3 +#!/usr/bin/env python3 import numpy as np import argparse @@ -6,23 +6,32 @@ def _parse_args(): - # Create argument parser + """Parse and return commandline arguments""" parser = argparse.ArgumentParser() parser.add_argument('-b', '--batches', type=int, default=100) parser.add_argument('-p', '--particles', type=int, default=10000) - group = parser.add_mutually_exclusive_group(required=True) - group.add_argument('-n', '--neutron_spectrum', action='store_true') - group.add_argument('-g', '--gamma_heating', action='store_true') - group.add_argument('-r', '--reaction_rates', action='store_true') parser.add_argument('-s', '--threads', type=int) + parser.add_argument('-c', '--cwd', type=str) + group = parser.add_argument_group('tallies') + group.add_argument('-n', '--neutron_spectrum', action='store_true', + default=False) + group.add_argument('-g', '--gamma_heating', action='store_true', + default=False) + group.add_argument('-r', '--reaction_rates', action='store_true', + default=False) - # Parse and return commandline arguments. return parser.parse_args() def main(): """Analysis of Beryllium Experiment""" - args = vars(_parse_args()) + + # Parse commandline arguments + args = _parse_args() + # If neither gamma heating or rxn rate tallies are desired, default to only + # the neutron spectrum tallies + if not (args.gamma_heating or args.reaction_rates): + args.neutron_spectrum = True # Instantiate Model object model = openmc.Model() @@ -58,7 +67,7 @@ def main(): tube.add_nuclide('Ni62', 0.0031349384315041144, 'ao') tube.add_nuclide('Ni64', 0.0007982681776143339, 'ao') - if args['reaction_rates']: + if args.reaction_rates: # Add reaction rate elements to materials in small quantities air.add_element('Al', 1e-8, 'ao') air.add_element('Nb', 1e-8, 'ao') @@ -150,36 +159,36 @@ def main(): region30 = (+s39 & -s1 & -s6) | (+s1 & -s4 & -s6 & +s5) | (+s4 & -s40 & -s6) # Cells - cell1 = openmc.Cell(region=region1, fill=air, name='cell 1') - cell2 = openmc.Cell(region=region2, fill=air, name='cell 2') - cell3 = openmc.Cell(region=region3, fill=beryllium, name='cell 3') - cell4 = openmc.Cell(region=region4, fill=tube, name='cell 4') - cell5 = openmc.Cell(region=+s15 & -s2 & xybox, fill=air, name='cell 5') - cell6 = openmc.Cell(region=+s2 & -s16 & xybox, fill=beryllium, name='cell 6') - cell7 = openmc.Cell(region=+s16 & -s17 & xybox, fill=beryllium, name='cell 7') - cell8 = openmc.Cell(region=+s17 & -s18 & xybox, fill=beryllium, name='cell 8') - cell9 = openmc.Cell(region=+s18 & -s19 & xybox, fill=beryllium, name='cell 9') - cell10 = openmc.Cell(region=+s19 & -s20 & xybox, fill=beryllium, name='cell 10') - cell11 = openmc.Cell(region=+s20 & -s21 & xybox, fill=beryllium, name='cell 11') - cell12 = openmc.Cell(region=+s21 & -s22 & xybox, fill=beryllium, name='cell 12') - cell13 = openmc.Cell(region=+s22 & -s23 & xybox, fill=beryllium, name='cell 13') - cell14 = openmc.Cell(region=+s23 & -s24 & xybox, fill=beryllium, name='cell 14') - cell15 = openmc.Cell(region=+s24 & -s25 & xybox, fill=beryllium, name='cell 15') - cell16 = openmc.Cell(region=+s25 & -s26 & xybox, fill=beryllium, name='cell 16') - cell17 = openmc.Cell(region=+s26 & -s27 & xybox, fill=beryllium, name='cell 17') - cell18 = openmc.Cell(region=+s27 & -s28 & xybox, fill=beryllium, name='cell 18') - cell19 = openmc.Cell(region=+s28 & -s29 & xybox, fill=beryllium, name='cell 19') - cell20 = openmc.Cell(region=+s29 & -s30 & xybox, fill=beryllium, name='cell 20') - cell21 = openmc.Cell(region=+s30 & -s31 & xybox, fill=beryllium, name='cell 21') - cell22 = openmc.Cell(region=+s31 & -s32 & xybox, fill=beryllium, name='cell 22') - cell23 = openmc.Cell(region=+s32 & -s33 & xybox, fill=beryllium, name='cell 23') - cell24 = openmc.Cell(region=+s33 & -s34 & xybox, fill=beryllium, name='cell 24') - cell25 = openmc.Cell(region=+s34 & -s35 & xybox, fill=beryllium, name='cell 25') - cell26 = openmc.Cell(region=+s35 & -s36 & xybox, fill=beryllium, name='cell 26') - cell27 = openmc.Cell(region=+s36 & -s3 & xybox, fill=beryllium, name='cell 27') - cell28 = openmc.Cell(region=+s3 & -s37 & xybox, fill=air, name='cell 28') - cell29 = openmc.Cell(region=-s38, fill=air, name='cell 29') - cell30 = openmc.Cell(region=region30, name='cell 30') + cell1 = openmc.Cell(region=region1, fill=air) + cell2 = openmc.Cell(region=region2, fill=air) + cell3 = openmc.Cell(region=region3, fill=beryllium) + cell4 = openmc.Cell(region=region4, fill=tube) + cell5 = openmc.Cell(region=+s15 & -s2 & xybox, fill=air) + cell6 = openmc.Cell(region=+s2 & -s16 & xybox, fill=beryllium) + cell7 = openmc.Cell(region=+s16 & -s17 & xybox, fill=beryllium) + cell8 = openmc.Cell(region=+s17 & -s18 & xybox, fill=beryllium) + cell9 = openmc.Cell(region=+s18 & -s19 & xybox, fill=beryllium) + cell10 = openmc.Cell(region=+s19 & -s20 & xybox, fill=beryllium) + cell11 = openmc.Cell(region=+s20 & -s21 & xybox, fill=beryllium) + cell12 = openmc.Cell(region=+s21 & -s22 & xybox, fill=beryllium) + cell13 = openmc.Cell(region=+s22 & -s23 & xybox, fill=beryllium) + cell14 = openmc.Cell(region=+s23 & -s24 & xybox, fill=beryllium) + cell15 = openmc.Cell(region=+s24 & -s25 & xybox, fill=beryllium) + cell16 = openmc.Cell(region=+s25 & -s26 & xybox, fill=beryllium) + cell17 = openmc.Cell(region=+s26 & -s27 & xybox, fill=beryllium) + cell18 = openmc.Cell(region=+s27 & -s28 & xybox, fill=beryllium) + cell19 = openmc.Cell(region=+s28 & -s29 & xybox, fill=beryllium) + cell20 = openmc.Cell(region=+s29 & -s30 & xybox, fill=beryllium) + cell21 = openmc.Cell(region=+s30 & -s31 & xybox, fill=beryllium) + cell22 = openmc.Cell(region=+s31 & -s32 & xybox, fill=beryllium) + cell23 = openmc.Cell(region=+s32 & -s33 & xybox, fill=beryllium) + cell24 = openmc.Cell(region=+s33 & -s34 & xybox, fill=beryllium) + cell25 = openmc.Cell(region=+s34 & -s35 & xybox, fill=beryllium) + cell26 = openmc.Cell(region=+s35 & -s36 & xybox, fill=beryllium) + cell27 = openmc.Cell(region=+s36 & -s3 & xybox, fill=beryllium) + cell28 = openmc.Cell(region=+s3 & -s37 & xybox, fill=air) + cell29 = openmc.Cell(region=-s38, fill=air) + cell30 = openmc.Cell(region=region30) cells = [cell1, cell2, cell3, cell4, cell5, cell6, cell7, cell8, cell9, cell10, cell11, cell12, cell13, cell14, cell15, cell16, cell17, cell18, @@ -187,6 +196,9 @@ def main(): cell27, cell28, cell29, cell30 ] + for c in cells: + c.name = f'cell {c.id}' + # Add Geometry to Model model.geometry = openmc.Geometry(cells) @@ -195,8 +207,8 @@ def main(): # Indicate how many particles to run settings = openmc.Settings(run_mode='fixed source') - settings.batches = args['batches'] - settings.particles = args['particles'] + settings.batches = args.batches + settings.particles = args.particles # Create source distribution # upper neutron energy bins @@ -263,19 +275,23 @@ def main(): ############################################################################### # Build tallies based on the desired run + model.tallies = openmc.Tallies() cell_filter = openmc.CellFilter(cells[4:26]) + usrcwd = args.cwd + cwd = '' - if args['neutron_spectrum']: + if args.neutron_spectrum: # neutron energy spectrum tally energy_n = openmc.EnergyFilter(x) particle_n = openmc.ParticleFilter(['neutron']) tally_n = openmc.Tally(name="Neutron Spectrum") tally_n.filters = [cell_filter, energy_n, particle_n] tally_n.scores = ['flux'] - model.tallies = openmc.Tallies([tally_n]) - cwd = 'neutron_spectrum' - elif args['gamma_heating']: - # Tallies gamma heating + model.tallies.append(tally_n) + cwd += 'neutron_spectrum' + + if args.gamma_heating: + # gamma heating tally settings.survival_biasing = True settings.photon_transport = True settings.electron_treatment = 'ttb' @@ -283,100 +299,101 @@ def main(): tally_p = openmc.Tally(name="Gamma heating") tally_p.filters = [cell_filter, particle_p] tally_p.scores = ['heating'] - model.tallies = openmc.Tallies([tally_p]) - cwd = 'gamma_heating' - - elif args['reaction_rates']: - particle_n = openmc.ParticleFilter(['neutron']) + model.tallies.append(tally_p) + cwd = cwd + '-gamma_heating' if cwd else 'gamma_heating' + if args.reaction_rates: + # reaction rate tallies tally_1 = openmc.Tally(name="Au-197 (n,gamma)") tally_1.nuclides = ['Au197'] - tally_1.filters = [cell_filter, particle_n] + tally_1.filters = [cell_filter] tally_1.scores = ['(n,gamma)'] tally_2 = openmc.Tally(name="In115 (n,n')") tally_2.nuclides = ['In115'] - tally_2.filters = [cell_filter, particle_n] + tally_2.filters = [cell_filter] tally_2.scores = ['(n,nc)'] tally_3 = openmc.Tally(name="Al-27 (n,alpha)") tally_3.nuclides = ['Al27'] - tally_3.filters = [cell_filter, particle_n] + tally_3.filters = [cell_filter] tally_3.scores = ['(n,a)'] tally_4 = openmc.Tally(name="Nb-93 (n,2n)") tally_4.nuclides = ['Nb93'] - tally_4.filters = [cell_filter, particle_n] + tally_4.filters = [cell_filter] tally_4.scores = ['(n,2n)'] tally_5 = openmc.Tally(name="Ni-58 (n,2n)") tally_5.nuclides = ['Ni58'] - tally_5.filters = [cell_filter, particle_n] + tally_5.filters = [cell_filter] tally_5.scores = ['(n,2n)'] tally_6 = openmc.Tally(name="Ti-48(n,p)") tally_6.nuclides = ['Ti48'] - tally_6.filters = [cell_filter, particle_n] + tally_6.filters = [cell_filter] tally_6.scores = ['(n,p)'] tally_7 = openmc.Tally(name="Ti-49 (n,np)") tally_7.nuclides = ['Ti49'] - tally_7.filters = [cell_filter, particle_n] + tally_7.filters = [cell_filter] tally_7.scores = ['(n,np)'] tally_8 = openmc.Tally(name="Ti-47 (n,p)") tally_8.nuclides = ['Ti47'] - tally_8.filters = [cell_filter, particle_n] + tally_8.filters = [cell_filter] tally_8.scores = ['(n,p)'] tally_9 = openmc.Tally(name="Ti-48 (n,np)") tally_9.nuclides = ['Ti48'] - tally_9.filters = [cell_filter, particle_n] + tally_9.filters = [cell_filter] tally_9.scores = ['(n,np)'] tally_10 = openmc.Tally(name="Fe-56 (n,p)") tally_10.nuclides = ['Fe56'] - tally_10.filters = [cell_filter, particle_n] + tally_10.filters = [cell_filter] tally_10.scores = ['(n,p)'] tally_11 = openmc.Tally(name="Zr-90 (n,2n)") tally_11.nuclides = ['Zr90'] - tally_11.filters = [cell_filter, particle_n] + tally_11.filters = [cell_filter] tally_11.scores = ['(n,2n)'] tally_12 = openmc.Tally(name="Ni-58 (n,p)") tally_12.nuclides = ['Ni58'] - tally_12.filters = [cell_filter, particle_n] + tally_12.filters = [cell_filter] tally_12.scores = ['(n,p)'] tally_13 = openmc.Tally(name="Li-6 (n,tritium)") tally_13.nuclides = ['Li6'] - tally_13.filters = [cell_filter, particle_n] + tally_13.filters = [cell_filter] tally_13.scores = ['(n,Xt)'] tally_14 = openmc.Tally(name="Tritium production") - tally_14.filters = [cell_filter, particle_n] + tally_14.filters = [cell_filter] tally_14.scores = ['H3-production'] tally_15 = openmc.Tally(name="U-235 (n,f)") tally_15.nuclides = ['U235'] - tally_15.filters = [cell_filter, particle_n] + tally_15.filters = [cell_filter] tally_15.scores = ['fission'] tally_16 = openmc.Tally(name="Be-9 (n,tritium)") tally_16.nuclides = ['Be9'] - tally_16.filters = [cell_filter, particle_n] + tally_16.filters = [cell_filter] tally_16.scores = ['(n,Xt)'] - model.tallies = openmc.Tallies([tally_1, tally_2, tally_3, tally_4, - tally_5, tally_6, tally_7, tally_8, - tally_9, tally_10, tally_11, tally_12, - tally_13, tally_14, tally_15, tally_16]) - cwd = 'reaction_rates' + model.tallies.extend([tally_1, tally_2, tally_3, tally_4, tally_5, + tally_6, tally_7, tally_8, tally_9, tally_10, + tally_11, tally_12, tally_13, tally_14, tally_15, + tally_16]) + + cwd = cwd + '-rxn_rates' if cwd else 'rxn_rates' + cwd = usrcwd if usrcwd is not None else cwd model.settings = settings - return model.run(cwd=cwd, threads=args['threads']) + return model.run(cwd=cwd, threads=args.threads) if __name__ == "__main__": From 9fbcb487479789837ff976a6052c99dad6fa8261 Mon Sep 17 00:00:00 2001 From: Ethan Peterson Date: Wed, 1 Mar 2023 22:51:58 -0500 Subject: [PATCH 4/4] put all arg manipulation in _parse_args --- fns_be/fns_be.py | 13 ++++++++----- 1 file changed, 8 insertions(+), 5 deletions(-) diff --git a/fns_be/fns_be.py b/fns_be/fns_be.py index 00adf35d..0c7d473e 100755 --- a/fns_be/fns_be.py +++ b/fns_be/fns_be.py @@ -20,7 +20,14 @@ def _parse_args(): group.add_argument('-r', '--reaction_rates', action='store_true', default=False) - return parser.parse_args() + args = parser.parse_args() + + # If neither gamma heating or rxn rate tallies are desired, default to only + # the neutron spectrum tallies + if not (args.gamma_heating or args.reaction_rates): + args.neutron_spectrum = True + + return args def main(): @@ -28,10 +35,6 @@ def main(): # Parse commandline arguments args = _parse_args() - # If neither gamma heating or rxn rate tallies are desired, default to only - # the neutron spectrum tallies - if not (args.gamma_heating or args.reaction_rates): - args.neutron_spectrum = True # Instantiate Model object model = openmc.Model()