Local adjoint source for Random Ray

[j-fletcher]

Description

Adds ability to specify fixed, localized adjoint sources for the random ray solver. This capability will ultimately be useful for the addition of automated CADIS weight-windowing and source biasing using the adjoint random ray solver, and in this implementation, enables FW-CADIS weight windows to be generated for local variance reduction.

There are two possible methods of constructing local adjoint sources. The first is by the straightforward specification of a geometric region of interest and energy discretization desired for the adjoint source; the initial forward calculation proceeds as previously to set $Q^{\dagger} = 1/\phi$, but only within flat source regions in the region of interest. Elsewhere, the adjoint source simply remains 0. This approach, although still requiring an initial forward calculation, ensures that variance reduction is performed in all regions of phase space counted as part of the response (e.g., each individual detector and each energy group therein) uniformly. The user can designate a target openmc.Tallies object containing each tally on which variance reduction should be performed, then the weight windower will use geometric and energy constraints inferred from CellFilter, CellInstanceFilter, DistribcellFilter, UniverseFilter, MaterialFilter, and EnergyFilter instances on each tally to populate adjoint sources in flat source regions, as described in the method above.

The second is by user specification of the adjoint source term directly, using the distribution objects already available for constructing forward openmc.SourceBase instances, for example when information is available about a certain detector's response function. This approach can potentially result in unequal variance across different detector responses if multiple such sources are used to generate weight window and source biasing parameters in CADIS, but removes the need to perform an initial forward solve when pursuing local variance reduction.

This PR includes both of these capabilities, with the original implementation of the forward-weighted method by @jtramm. It includes an additional model in openmc.examples, based on the existing random_ray_three_region_cube with the addition of two "detector" volumes, which can be used to demonstrate both options as shown below.

Specifying Forward-Weighted Local Adjoint Sources for FW-CADIS

First, we'll import the test geometry and define a weight window mesh over the entire geometry.

import openmc
from openmc.examples import random_ray_three_region_cube_with_detectors

model = random_ray_three_region_cube_with_detectors()

ww_mesh = openmc.RegularMesh()
n = 7
width = 35.0
ww_mesh.dimension = (n, n, n)
ww_mesh.lower_left = (0.0, 0.0, 0.0)
ww_mesh.upper_right = (width, width, width)

Next, although there are several tallies included in the model, we'll pick out the tallies which cover our phase space regions of interest (detectors 1 and 2). We'll then add these to an FW-CADIS WeightWindowGenerator on the model and generate our weight windows for local variance reduction:

target_tallies = openmc.Tallies()
for tally in list(model.tallies):
    if tally.name in {"Detector 1 Tally", "Detector 2 Tally"}:
        target_tallies.append(tally)

wwg = openmc.WeightWindowGenerator(
    method="fw_cadis", 
    targets=target_tallies,
    mesh=ww_mesh, 
    max_realizations=model.settings.batches
)
model.settings.weight_window_generators = wwg
model.settings.random_ray['volume_estimator'] = 'naive'

model.run()

The resulting weight window lower bounds, depicted below, show a decreasing trend in the direction from the forward source (in the "bottom left" corner; $(X, Y, Z) = [0.0, 0.0, 0.0]\times [5.0, 5.0, 5.0]$) towards Detector 1 (centered at (2.5, 32.5, 2.5)) and especially towards Detector 2 (centered further from the source, at (32.5, 32.5, 32.5)), as desired for simultaneous local variance reduction on both tallies.

Specifying Localized Sources for Non-Forward-Weighted Adjoint Simulations

Outside of FW-CADIS weight window generation, this PR also adds the ability to specify localized sources for pure adjoint simulations. If supplied by the user, this removes the requirement to run an initial forward solve before performing adjoint transport, but in general the forward-weighted method is still recommended for applications like weight window generation with multiple detectors because the uniform reduction of variance in tally regions of interest is guaranteed.

At present, these sources are subject to the same constraints as all forward fixed sources for random ray: they must be isotropic in angle, use a discrete (i.e. multigroup) energy distribution, and the spatial distribution must either be a point, or constrained by Cell, Universe or Material. In this demonstration, we use a Box distribution identical to that specified on model.settings.random_ray['ray_source'], but constrain it to the flat source regions which make up Detector 1. As in the previous example, the problem is run with a single energy group.

model = random_ray_three_region_cube_with_detectors()

detector1_cells = model.geometry.get_cells_by_name("detector 1")

strengths = [1.0]
midpoints = [100.0]
energy_distribution = openmc.stats.Discrete(x=midpoints, p=strengths)

adj_space = openmc.stats.Box([0.0, 0.0, 0.0], [35.0, 35.0, 35.0], only_fissionable=False)
adj_source = openmc.IndependentSource(energy=energy_distribution, space=adj_space, constraints={'domains': detector1_cells})

model.settings.random_ray['adjoint'] = True
model.settings.random_ray['adjoint_source'] = adj_source
model.settings.random_ray['volume_estimator'] = 'naive'

model.run()

As shown below, the addition of the 'adjoint_source' entry has allowed us to map the adjoint flux with respect to Detector 1 throughout the geometry without running an initial forward calculation.

Fixes #3710

[jtramm]

PR is looking fantastic -- nice work on this! Just a few minor questions/comments to consider.

[jtramm]

Looks good -- nice work @j-fletcher!