Obtaining initial source for fix source mode

Hi all

I’m trying to obtain an initial_source.h5 file for a fixed source simulation. I have been using a bit of a hack to get this to work in version 0.12.2 but it appears to no longer work on the develop branch.

I’m just wondering if anyone knows what has changed to stop this working.

import xml.etree.ElementTree as ET
import numpy as np
import openmc


def create_initial_particles(source, number_of_particles=2000):
    """Accepts an openmc source and creates an initial_source.h5 that can be
    used to find initial xyz, direction and energy of the partice source
    """

    # no real materials are needed for finding the source
    mats = openmc.Materials([])

    # just a minimal geometry
    outer_surface = openmc.Sphere(r=100, boundary_type="vacuum")
    cell = openmc.Cell(region=-outer_surface)
    universe = openmc.Universe(cells=[cell])
    geom = openmc.Geometry(universe)

    # Instantiate a Settings object
    sett = openmc.Settings()
    sett.run_mode = (
        "eigenvalue"  # this will fail but it will write the initial_source.h5 file first
        # "fixed source"  # write_initial_source is not supported by this run mode
    )
    sett.particles = number_of_particles
    sett.batches = 1
    sett.inactive = 0
    sett.write_initial_source = True

    sett.source = source

    model = openmc.model.Model(geom, mats, sett)

    model.export_to_xml()

    # this just adds write_initial_source == True to the settings.xml
    tree = ET.parse("settings.xml")
    root = tree.getroot()
    elem = ET.SubElement(root, "write_initial_source")
    elem.text = "true"
    tree.write("settings.xml")

    # This will crash but it should write the initial_source.h5
    openmc.run()


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])

output from openmc 0.12.2 with initial source .h5 file produced

output from openmc 0.13 dev branch with initial source .h5 not file produced

Have you got OpenMC built in debug mode? You can attach GDB and probe from there, guessing in the case it bugs out before writing the file. Is there not a ‘less broken’ way you can create this file?

@Shimwell I tried this on the develop branch and I was able to produce an initial_source.h5 file without running into a segfault. Looking at the traceback you have there, it looks like something went wrong with parallel HDF5. Your parallel HDF5 installation is linked against MPICH, but the fact that libmpi.so.40 shows up in the backtrace makes me suspect you compiled with OpenMPI.

That aside, I love your ingenuity in figuring out how to get an initial source file, but we (OpenMC developers) ought to make it easier! At the very least, it looks like write_initial_source is not connected to the Python API which should be fixed.