feature request: dilute_initial to add inital isotopes based on material

[shimwell]

The openmc.deplete.Operator accepts a dilute_initial argument that is a boolean.

When set to True a small quantity of burnable nuclide are added to the material which helps numerically solve the inventory calculation.

However dilute_initial can add nuclides that would not be generated by the irradiation of the material.

For example using dilute_initial=True when irradiating Iron with neutrons would add nuclides that would not be created such as U, Pu, Cs to the material inventory.

minimal example to reproduce

import openmc
import openmc.deplete

# MATERIALS

model = openmc.Model()

material = openmc.Material()
material.add_nuclide('Fe56', 1)
material.set_density('g/cm3',1)
material.volume = 100000
material.depletable = True

materials = openmc.Materials([material])
model.materials = materials

# GEOMETRY

# surfaces
sph1 = openmc.Sphere(r=75, boundary_type='vacuum')

# cells, makes a simple sphere cell
shield_cell = openmc.Cell(region=-sph1)
shield_cell.fill = material
shield_cell.volume = 100000

# sets the geometry to the universe that contains just the one cell
universe = openmc.Universe(cells=[shield_cell])
geometry = openmc.Geometry(universe)
model.geometry = geometry

# creates a simple isotropic neutron source in the center with 14MeV neutrons
my_source = openmc.Source()
my_source.space = openmc.stats.Point((0, 0, 0))
my_source.angle = openmc.stats.Isotropic()
my_source.energy = openmc.stats.Discrete([14e6], [1])

# specifies the simulation computational intensity
settings = openmc.Settings()
settings.batches = 2
settings.particles = 1_000_00
settings.inactive = 0
settings.run_mode = "fixed source"
settings.source = my_source

model.settings = settings

chain_filename = "chain-nndc-b7.1.xml"
openmc.deplete.Chain.from_xml(chain_filename)

geometry.export_to_xml()
settings.export_to_xml()
materials.export_to_xml()
operator = openmc.deplete.Operator(
    model=model,
    chain_file=chain_filename,
    normalization_mode="source-rate"
)

time_steps = [1] + [60] * (60*24) # 1 second followed by 60 second blocks for 24 hours
source_rates = [1e21] + [0.] * (60*24) # 1 pulse followed by cooling

integrator = openmc.deplete.PredictorIntegrator(
    operator=operator,
    timesteps=time_steps,
    source_rates=source_rates,
)

integrator.integrate()

results = openmc.deplete.Results.from_hdf5("depletion_results.h5")

all_materials = []
for step in range(len(results.get_times())):
    material = results.export_to_materials(step)[0]  # zero index as one mat in problem
    all_materials.append(material)
    print(material)

[paulromano]

@shimwell I just realized that if you pass reduce_chain=True to the Operator class, you should get exactly the behavior that you want: 1) the chain will be reduced to only the isotopes present in depletable materials and their possible progeny (determined by walking through the chain), 2) the "burnable" nuclides are determined based on the chain, and 3) dilute_initial is applied based on the burnable nuclides.

[shimwell]

oh nice,

So for a "typical" fixed source activation simulation in fusion I guess that I would set reduce_chain, dilute_initial and normalization_mode in the various places. Here is a complete example in case it is handy for anyone else.

operator = openmc.deplete.Operator(
    model=model,
    chain_file='chain-nndc-b7.1.xml',
    normalization_mode="source-rate",  # set for fixed source simulation, otherwise defaults to fission simulation
    dilute_initial=0  # set to zero to avoid adding small amounts of isotopes, defaults to adding small amounts of fissionable isotopes
    reduce_chain=True  # reduced to only the isotopes present in depletable materials and their possible progeny
)

time_steps = [365*24*60*60] * 5 + [365*24*60*60] * 5
source_rates = [1e21]*5 + [0] * 5

integrator = openmc.deplete.PredictorIntegrator(
    operator=operator, timesteps=time_steps,source_rates=source_rates
)

integrator.integrate()

I think this issue can now be closed as we don't need to change the dilute_initial method. Thanks @paulromano