<jupyter_start><jupyter_text>Openmc Overview The purpose of this notebook is to provide an introduction to the capabilities of OpenMC.  The coverage is intended to provide a fairly deep look at a geometrically simple problem of an infinite pin-cell. The basic topics to be covered will be:1. material specification2. geometry description3. tally specification4. control parameter specificationThe analysis will include a simple criticality calculation and neutron flux energy spectrum measurement.  This will be followed by a somewhat more elaborate treatment of fuel material specification and some exploration into the consequences of one or two important details of fuel material composition.  A burn-up analysis will also be conducted to demonstrate the code capabilities and get a first look at the evolution of important fission products.<jupyter_code>import matplotlib.pyplot as plt
import numpy as np
%matplotlib inline
import math
import openmc<jupyter_output><empty_output><jupyter_text>Material Specification First we will define the fuel<jupyter_code>fuel = openmc.Material(name="uo2")
fuel.add_element("U",1,percent_type="ao", enrichment=4.25,enrichment_type='wo')
fuel.add_element("O",2)
fuel.set_density("g/cc",10.4)<jupyter_output><empty_output><jupyter_text>Documentation on the OpenMC Material object, along with all other OpenMC features can be found [here](https://docs.openmc.org/en/stable/) Next we define the cladding<jupyter_code>clad = openmc.Material(name='clad')
clad.add_element("Zr",1)
clad.set_density('g/cc',6.0)<jupyter_output><empty_output><jupyter_text>and water moderator.<jupyter_code>water = openmc.Material(name='water')
water.add_element('O',1)
water.add_element('H',2)
water.set_density('g/cc',1.0)
water.add_s_alpha_beta('c_H_in_H2O')<jupyter_output><empty_output><jupyter_text>The last line adds some additional physics.  The kinematics of neutron scattering interactions depends on whether the target atom is *free* or *bound in a molecule*.  This dependence is stronger for light atoms.  The required modification to scattering interactions is captured in, so called, $S(\alpha,\beta)$ relation.  In the code block above we "add" the $S(\alpha,\beta)$ relation for hydrogen bound in water.Any time you have a significant atom fraction of light atoms bound in molecules you should be looking to add an appropriate $S(\alpha,\beta)$. These are all of the materials that we are going to include for this simple model.  In order for the OpenMC neutron transport code to use the materials, they need to be written to an XML file called 'materials.xml'. The OpenMC Python interface provides tools for doing this in the `Materials` object.<jupyter_code>materials = openmc.Materials()
materials += [fuel, clad, water]
materials.export_to_xml()<jupyter_output><empty_output><jupyter_text>GeometryNext we will define the geometry.  As a first example we will create an infinite array of infinitely long pin cells.The geometry definition process we will follow for this example will comprise these steps:1. Definition of surfaces (in some cases including boundary conditions)2. Definition of regions which are collections of half-spaces3. Definition of cells which combines half-spaces with materials; and4. Organization of cells into a universe.  Surface Definitions<jupyter_code>R_fo = 0.39 # cm, fuel outer radius
R_ci = 0.40 # cm, clad inner radius
R_co = 0.46 # cm, clad outer radius
Pitch = 1.26 # cm, distance between centerline of adjacent cells

fuel_outer_radius = openmc.ZCylinder(r=R_fo)
clad_inner_radius = openmc.ZCylinder(r=R_ci)
clad_outer_radius = openmc.ZCylinder(r=R_co)

left = openmc.XPlane(x0=-Pitch/2., boundary_type='reflective')
right = openmc.XPlane(x0=Pitch/2., boundary_type='reflective')
bottom = openmc.YPlane(y0=-Pitch/2., boundary_type='reflective')
top = openmc.YPlane(y0=Pitch/2., boundary_type='reflective')<jupyter_output><empty_output><jupyter_text>Definition of Regions<jupyter_code>fuel_region = -fuel_outer_radius
gap_region = +fuel_outer_radius & -clad_inner_radius
clad_region = +clad_inner_radius & -clad_outer_radius

water_region = +left & -right & +bottom & -top & +clad_outer_radius<jupyter_output><empty_output><jupyter_text>Definition of Cells<jupyter_code>fuel_cell = openmc.Cell(name='fuel')
fuel_cell.fill = fuel
fuel_cell.region = fuel_region

gap_cell = openmc.Cell(name='gap')
gap_cell.region = gap_region
gap_cell.fill = None # no material = void

clad_cell = openmc.Cell(name='clad')
clad_cell.region = clad_region
clad_cell.fill = clad

water_cell = openmc.Cell(name='moderator')
water_cell.region = water_region
water_cell.fill = water<jupyter_output><empty_output><jupyter_text>Definition of a Universe<jupyter_code>root_universe = openmc.Universe(cells=(fuel_cell,gap_cell,clad_cell,water_cell))

geometry = openmc.Geometry()
geometry.root_universe = root_universe
geometry.export_to_xml()<jupyter_output><empty_output><jupyter_text>Plot the Geometry<jupyter_code>colors = {}
colors[fuel] = 'yellow'
colors[clad] = 'grey'
colors[water] = 'blue'
p = openmc.Plot()
p.width = [Pitch,Pitch]
p.pixels = [400,400]
p.origin = [0.0, 0.0, 0.0]
p.color_by = 'material'
p.colors = colors

openmc.plot_inline(p)<jupyter_output><empty_output><jupyter_text>Now that the geometry and material composition of the reactor is specified, we need to tell OpenMC what information we hope to get out of the model.  This is done through *tallies*.  Re-do Geometry with a Finite Pincell<jupyter_code>h_core = 300.; 
h_fuel = 200.;

r_fuel = 0.42;
r_pin = 0.45;

P_D = 1.4;

pitch = P_D*2*r_pin;

fuel_temp = 900; # K, guess at fuel temperature
mod_temp = 600; # K, moderator temperature

# fuel cylinder:
fuel_cyl = openmc.model.RightCircularCylinder([0.,0.,-h_fuel/2.],
                                             h_fuel, r_fuel);

fuel.volume = np.pi*(r_fuel**2)*h_fuel;

# pin cylinder
pin_cyl = openmc.model.RightCircularCylinder([0.,0.,-(h_fuel+(r_pin-r_fuel))/2.],
                                            h_fuel+(r_pin-r_fuel)*2.,r_pin);

# pin cell container
core_cell = openmc.model.RectangularParallelepiped(-pitch/2.,pitch/2.,
                                                  -pitch/2.,pitch/2.,
                                                  -h_core/2.,h_core/2.,
                                                  boundary_type='reflective');

fuel_cell = openmc.Cell();
fuel_cell.region = -fuel_cyl
fuel_cell.fill = fuel;
fuel_cell.temperature = fuel_temp;


clad_cell = openmc.Cell();
clad_cell.region = +fuel_cyl & -pin_cyl;
clad_cell.fill = clad;

mod_cell = openmc.Cell();
mod_cell.region = +pin_cyl & -core_cell;
mod_cell.fill = water

root_univ = openmc.Universe();
root_univ.add_cells([fuel_cell,clad_cell,mod_cell]);

geometry = openmc.Geometry();
geometry.root_universe = root_univ;

materials.export_to_xml();
geometry.export_to_xml();<jupyter_output><empty_output><jupyter_text>Tally Specification<jupyter_code># flux energy spectrum in fuel cell
cell_filter = openmc.CellFilter(fuel_cell)

energy_bins = np.logspace(-3,7,num=1001)
energy_filter = openmc.EnergyFilter(values=energy_bins)

flux_tally = openmc.Tally(name='flux')
flux_tally.scores = ['flux']
flux_tally.filters = [cell_filter,energy_filter]


tallies = openmc.Tallies([flux_tally])

tallies.export_to_xml()<jupyter_output><empty_output><jupyter_text>Settings<jupyter_code>settings = openmc.Settings()
settings.batches = 300
settings.inactive = 50
settings.particles = 20000

point = openmc.stats.Point((0,0,0))
source = openmc.Source(space=point)

settings.source = source

settings.export_to_xml()
openmc.run()<jupyter_output>%%%%%%%%%%%%%%%
                           %%%%%%%%%%%%%%%%%%%%%%%%
                        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                      %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                                    %%%%%%%%%%%%%%%%%%%%%%%%
                                     %%%%%%%%%%%%%%%%%%%%%%%%
                 ###############      %%%%%%%%%%%%%%%%%%%%%%%%
                ##################     %%%%%%%%%%%%%%%%%%%%%%%
                ###################     %%%%%%%%%%%%%%%%%%%%%%%
                ####################     %%%%%%%%%%%%%%%%%%%%%%
                #####################     %%%%%%%%%%%%%%%%%%%%%
                ######################     %%%%%%%%%%%%%%%%%%%%
                #######################     %%%%%%%%%%%%%%%%%%
                 #######################     %%%%%%%%%%%%%%%%%
                 #####################[...]<jupyter_text>Access Tally and Plot Flux Spectrum<jupyter_code>sp = openmc.StatePoint('statepoint.300.h5')
sp.tallies
flux = sp.get_tally(name='flux')
flux_df = flux.get_pandas_dataframe()
flux_vals = flux_df['mean'].to_numpy()

energy_x = 0.5*(energy_bins[0:-1]+energy_bins[1:])

plt.loglog(energy_x,flux_vals)
plt.grid()
plt.xlabel('Energy [eV]')
plt.ylabel('flux [n/cm**2-s]')<jupyter_output><empty_output><jupyter_text>Notice the deep "pits" in the flux energy distribution.  It is worthwhile to investigate why those features exist.  To do so, more tallies could be added; specifically to measure capture - in total, as well as by specific isotopes.   Bonu Section: Depletion<jupyter_code>import openmc.deplete

chain = openmc.deplete.Chain.from_xml("chain_endfb71_pwr.xml")

chain.nuclide_dict
operator = openmc.deplete.Operator(geometry,settings,"chain_endfb71_pwr.xml")

power = 1e4;
days = 24*3600; # number of seconds in a day

time_steps = [0.1*days, 0.1*days, 0.1*days, 1.*days, 30*days, 30*days, 180*days]

integrator = openmc.deplete.PredictorIntegrator(operator,time_steps,power)

integrator.integrate()<jupyter_output>%%%%%%%%%%%%%%%
                           %%%%%%%%%%%%%%%%%%%%%%%%
                        %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                      %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                   %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
                                    %%%%%%%%%%%%%%%%%%%%%%%%
                                     %%%%%%%%%%%%%%%%%%%%%%%%
                 ###############      %%%%%%%%%%%%%%%%%%%%%%%%
                ##################     %%%%%%%%%%%%%%%%%%%%%%%
                ###################     %%%%%%%%%%%%%%%%%%%%%%%
                ####################     %%%%%%%%%%%%%%%%%%%%%%
                #####################     %%%%%%%%%%%%%%%%%%%%%
                ######################     %%%%%%%%%%%%%%%%%%%%
                #######################     %%%%%%%%%%%%%%%%%%
                 #######################     %%%%%%%%%%%%%%%%%
                 #####################[...]
{'fraction_non_alphanumeric': 0.19580288551620761, 'fraction_numerical': 0.013771781899943788, 'mean_word_length': 3.921622867680959, 'pattern_counts': {'":': 0, '<': 59, '<?xml version=': 0, '>': 59, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 10}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 14, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '15', 'dir': 'jupyter-structured-clean-dedup', 'id': '18994', 'license': 'MIT', 'n_tokens_mistral': 3389, 'n_tokens_neox': 3067, 'n_words': 739, 'path': 'examples/overview/.ipynb_checkpoints/OpenMC Overview-checkpoint.ipynb', 'repo_name': 'stu314159/er362_omc', 'size': '10674'}