4.1.4. Absorption Reaction Rate

This tutorial uses VolumePostprocessor with the sigma_a cross-section multiplier to integrate the absorption reaction rate over a selected volume.

4.1.4.1. Create a reference transport problem

The \(2 \times 1\) cm domain contains a purely absorbing material with \(\Sigma_a=2\) \(\text{cm}^{-1}\) and a uniform source of strength 2. Reflecting boundaries produce the uniform solution \(\phi=1\).

from mpi4py import MPI
from pyopensn.aquad import GLCProductQuadrature2DXY
from pyopensn.context import Finalize
from pyopensn.logvol import RPPLogicalVolume
from pyopensn.mesh import OrthogonalMeshGenerator
from pyopensn.post import VolumePostprocessor
from pyopensn.solver import DiscreteOrdinatesProblem, SteadyStateSourceSolver
from pyopensn.source import VolumetricSource
from pyopensn.xs import MultiGroupXS

rank = MPI.COMM_WORLD.rank
x_nodes = [0.1 * i for i in range(21)]
y_nodes = [0.1 * i for i in range(11)]
mesh = OrthogonalMeshGenerator(node_sets=[x_nodes, y_nodes]).Execute()
mesh.SetUniformBlockID(0)

xs = MultiGroupXS()
xs.CreateSimpleOneGroup(sigma_t=2.0, c=0.0)
source = VolumetricSource(block_ids=[0], group_strength=[2.0])
quadrature = GLCProductQuadrature2DXY(n_polar=2, n_azimuthal=8, scattering_order=0)
problem = DiscreteOrdinatesProblem(
    mesh=mesh,
    num_groups=1,
    groupsets=[{"groups_from_to": (0, 0), "angular_quadrature": quadrature}],
    xs_map=[{"block_ids": [0], "xs": xs}],
    volumetric_sources=[source],
    boundary_conditions=[
        {"name": name, "type": "reflecting"}
        for name in ("xmin", "xmax", "ymin", "ymax")
    ],
)
solver = SteadyStateSourceSolver(problem=problem)
solver.Initialize()
solver.Execute()

4.1.4.2. Integrate \(\Sigma_a\phi\) over a subvolume

The selected rectangle is \(0.5 \le x \le 1.5\) cm and \(0.2 \le y \le 0.8\) cm. The same postprocessor can instead use block_ids when the desired region follows material blocks.

sample_volume = RPPLogicalVolume(
    xmin=0.5, xmax=1.5, ymin=0.2, ymax=0.8, infz=True
)
absorption = VolumePostprocessor(
    problem=problem,
    value_type="integral",
    logical_volumes=[sample_volume],
    xs_multiplier="sigma_a",
)
absorption.Execute()
absorption_rate = float(absorption.GetValue()[0][0])

if rank == 0:
    print(f"Absorption reaction rate={absorption_rate:.8e}")
assert abs(absorption_rate - 1.2) < 1.0e-6

if "opensn_console" not in globals():
    from IPython import get_ipython
    if get_ipython() is not None:
        Finalize()
        MPI.Finalize()