Diffusion Synthetic Acceleration for polytopic discretisations of Boltzmann transport

This paper presents a computational study demonstrating that a modified interior penalty (MIP) formulation for Diffusion Synthetic Acceleration (DSA) applied to polytopic discontinuous Galerkin discretizations of $S_N$ transport equations maintains robust convergence across various optical and scattering regimes, outperforming the classical symmetric interior penalty (SIP) approach which can lose stability in intermediate conditions.

The Gist

Imagine you are trying to navigate a very crowded, foggy room (this is the Boltzmann Transport Equation). You want to get from one side to the other, but the room is filled with people (particles) bumping into each other constantly.

In the world of nuclear reactors, medical radiation therapy, and space shielding, scientists use computers to simulate this "foggy room" to predict how radiation moves. The computer tries to solve this by taking small steps, checking where the particles are, and adjusting. This process is called Source Iteration.

The Problem: The "Stuck in Mud" Effect

When the room is very foggy and people are bumping into each other constantly (a highly scattering environment), the computer's simple step-by-step method gets stuck. It's like trying to walk through waist-deep mud; you take a step, slip back a bit, take another step, and slip back again. It takes forever to get anywhere, or the computer might just give up and say, "I'm done," even though it hasn't found the right answer.

The Solution: The "GPS Shortcut" (DSA)

To fix this, scientists use a trick called Diffusion Synthetic Acceleration (DSA).

Think of the computer's main calculation as a hiker trying to find a path through a dense forest. The hiker is slow and gets lost easily.

• The Trick: The hiker pulls out a GPS (the Diffusion model). The GPS doesn't know every single tree, but it knows the general shape of the forest and the direction of the exit. It gives the hiker a "correction" or a shortcut to get back on track.
• The Result: Instead of wandering aimlessly, the hiker takes the GPS's advice, corrects their path, and reaches the destination much faster.

The Twist: Two Types of GPS (SIP vs. MIP)

The paper investigates two different ways to build this "GPS" (mathematically called Discontinuous Galerkin methods) on complex, irregular maps (called Polytopic Meshes—think of a map made of irregularly shaped puzzle pieces rather than perfect squares).

1. The Standard GPS (SIP - Symmetric Interior Penalty):
   This is the traditional way of building the GPS. It works great in clear weather or light fog. However, the paper found that in thick fog (optically thick regimes), this GPS sometimes gets confused. It might tell the hiker to go in the wrong direction, causing the whole process to crash or stall. It's like a GPS that works fine on a highway but fails completely in a dense jungle.

2. The "Toughened" GPS (MIP - Modified Interior Penalty):
   This is the new, improved version tested in the paper. The scientists realized that in the thick fog, the standard GPS wasn't "punishing" the hiker enough for taking wrong turns. The MIP method adds a "safety floor." It ensures that even in the thickest fog, the GPS remains strong and doesn't lose its grip.

   • The Analogy: If the Standard GPS is a regular umbrella, the MIP is a heavy-duty storm umbrella. In light rain, they are the same. In a hurricane, the regular one flips inside out, but the heavy-duty one keeps you dry.

What the Experiments Showed

The researchers ran thousands of simulations (like running the hiker through different types of forests) to see which GPS worked best. They tested:

• Fog density: How thick the fog was.
• Map complexity: How weirdly shaped the puzzle pieces were.
• Step size: How detailed the map was.

The Findings:

• The Standard GPS (SIP) often failed when the fog got too thick or the map got too weird. It would stop working entirely.
• The Toughened GPS (MIP) worked reliably in almost every situation. Even in the most challenging, thick-fog scenarios, it kept the hiker moving forward, often cutting the time needed to solve the problem by more than half.
• The Boundary Conditions: They also tested different rules for the edges of the room (like "bounce back" vs. "walk out"). While these mattered a little bit, the choice between the Standard and Toughened GPS was the most important factor.

Why This Matters

This research is a big deal for:

• Nuclear Engineers: Designing safer reactors and better shielding.
• Doctors: Planning radiation therapy to kill cancer cells without hurting healthy tissue.
• Space Agencies: Protecting astronauts from cosmic radiation.

By proving that the MIP method is robust, the paper gives engineers a reliable tool to simulate these complex, dangerous environments faster and more accurately, without the computer getting stuck in the "mud."

In short: The paper says, "If you are simulating radiation in a complex, foggy world, don't use the old, fragile GPS. Use the new, reinforced MIP GPS, and you'll get your answer fast and without crashing."

Technical Summary

1. Problem Statement

The paper addresses the computational challenges associated with solving the monoenergetic, isotropically scattering Boltzmann Transport Equation (BTE) using Discrete Ordinates ($S_N$) methods discretized in space via Polytopic Discontinuous Galerkin (DG) methods.

• The Core Issue: In optically thick and highly scattering regimes (common in nuclear shielding, reactor physics, and radiotherapy), standard Source Iteration (SI) methods suffer from severe slow convergence or pseudo-convergence. The contraction factor approaches unity, making the solution process computationally prohibitive.
• The Acceleration Need: Diffusion Synthetic Acceleration (DSA) is the standard remedy, approximating the error dynamics with a diffusion equation to provide an additive correction. However, for DSA to be robust, the diffusion correction must be discretized in a manner compatible with the transport discretization.
• Specific Challenge: On general polytopic meshes (e.g., Voronoi tessellations) and with high-order DG methods, standard Symmetric Interior Penalty (SIP) diffusion formulations can fail. Specifically, in thick regimes, the SIP penalty may become too weak relative to the upwind dissipation of the transport solver, leading to instability or divergence of the accelerated iteration.

2. Methodology

The authors conduct a comprehensive computational study comparing different DSA formulations on polytopic meshes.

• Transport Discretization: The $S_N$ equations are solved using a Polytopic Discontinuous Galerkin (DGFEM) method with upwind numerical fluxes. The spatial domain is discretized using bounded Voronoi tessellations.
• Acceleration Scheme: The study implements DSA where the diffusion correction is solved on the same mesh and in the same polynomial space as the transport problem.
• Key Comparison: The paper contrasts two interior-penalty diffusion formulations:
  1. Classical Symmetric Interior Penalty (SIP): Uses standard arithmetic averages and penalty parameters scaled by element size and polynomial degree.
  2. Modified Interior Penalty (MIP): A variant where the penalty parameter is enforced to have a "transport-scale floor." Specifically, the penalty on each face is set to $\max(\sigma_{SIP}, C_{face})$, where $C_{face}$ is derived from the angular quadrature weights and face orientation. This ensures the diffusion penalty never under-penalizes relative to the transport dissipation.
• Boundary Conditions: The study evaluates two boundary treatments for the diffusion equation:
  • Homogeneous Dirichlet: $\phi = 0$ on the boundary.
  • Marshak (Robin): An asymptotically correct condition derived from the diffusion limit ($n \cdot (D\nabla \phi) + \kappa \phi = 0$).
• Experimental Design: The authors perform numerical experiments varying:
  • Total macroscopic cross-section ($\sigma_t$) and scattering ratio ($c$).
  • Angular quadrature order ($N_Q$).
  • Mesh refinement ($h$-refinement) and polynomial degree ($p$-refinement).
  • Mesh anisotropy (using Lloyd-smoothed vs. highly distorted Voronoi meshes).
  • Metric: The primary metric is the empirical convergence factor ($\rho$) of the outer iteration loop.

3. Key Contributions

1. Systematic Evaluation on Polytopal Meshes: This is one of the first dedicated empirical studies quantifying DSA performance on general polytopic (Voronoi) meshes, addressing a gap in the literature where most studies focus on simplicial or Cartesian grids.
2. Validation of MIP Robustness: The study provides strong numerical evidence that MIP-based DSA maintains robustness across a wide range of optical thicknesses and scattering ratios, whereas SIP-based DSA loses stability in the intermediate-to-thick regimes.
3. Penalty Scaling Analysis: The paper demonstrates that the failure of SIP in thick regimes is due to the penalty parameter becoming too small relative to the upwind transport dissipation. The MIP formulation explicitly corrects this by matching the transport dissipation scale face-by-face.
4. Boundary Condition Impact: The study quantifies the trade-offs between Dirichlet and Marshak boundary conditions, showing that while Marshak is asymptotically correct, Dirichlet can sometimes yield better convergence in intermediate regimes, though both converge to similar behavior in the deep diffusion limit.
5. Cost-Benefit Analysis: The authors analyze the computational overhead of DSA, showing that while the diffusion solve adds cost per iteration, the reduction in iteration count leads to significant speed-ups, particularly as the number of angular ordinates ($N_Q$) increases.

4. Key Results

• Convergence Stability:
  • SIP: Exhibits a clear "loss of robustness." As optical thickness ($\sigma_t$) increases, SIP-based DSA often diverges or stalls in the intermediate regime (e.g., $\sigma_t \approx 17$ in their tests).
  • MIP: Remains robust and convergent for all tested parameters. In challenging optically thick, highly scattering settings ($c=0.999$), MIP-based schemes achieve convergence factors typically below 0.6.
• Parameter Sensitivity:
  • Scattering Ratio ($c$): As $c \to 1$, unaccelerated SI fails completely. MIP-DSA continues to converge, though iteration counts increase.
  • Angular Quadrature ($N_Q$): The cost of the diffusion correction is largely independent of $N_Q$. As $N_Q$ increases, the transport sweep dominates the cost, and the DSA overhead is effectively amortized, leading to greater total speed-ups.
  • Mesh Refinement ($h$) and Polynomial Degree ($p$): Convergence factors collapse onto a curve when plotted against optical thickness ($h\sigma_t$). MIP remains effective under $p$-refinement (up to $p=3$ tested), whereas SIP diverges earlier for higher $p$.
  • Anisotropy: Highly anisotropic meshes degrade performance for both methods, but MIP significantly mitigates the instability region compared to SIP.
• Boundary Conditions: The choice between Dirichlet and Marshak has a secondary effect. Dirichlet often performs slightly better in the intermediate regime, while Marshak is theoretically preferred for the diffusion limit. However, MIP is robust to both choices.

5. Significance

This work is significant for the advancement of deterministic transport solvers in complex geometries:

• Enabling Polytopic Meshes: It validates the use of polytopic DG methods (which offer superior geometric flexibility and reduced element counts for complex domains) in conjunction with DSA, a critical requirement for large-scale nuclear and medical physics simulations.
• Algorithmic Reliability: By identifying the specific failure mode of standard SIP formulations in thick regimes and confirming MIP as a robust alternative, the paper provides a practical guideline for implementing stable accelerators in modern transport codes.
• Efficiency: The results confirm that DSA is essential for high-fidelity simulations in scattering-dominated regimes, offering substantial reductions in wall-clock time compared to unaccelerated source iteration, provided the correct penalty formulation (MIP) is used.

In summary, the paper establishes that Modified Interior Penalty (MIP) is the necessary formulation for Diffusion Synthetic Acceleration when using Polytopic Discontinuous Galerkin methods, ensuring robust convergence across the full spectrum of transport regimes, including highly scattering and optically thick environments.