A Method for Thermal Radiation Transport Using Backward Characteristic Tracing

TL;DR

The authors present a backward-in-time long-characteristics method (MOC) for thermal radiation transport and couple it to a gray low-order (LO) closure to enforce energy conservation in radiation-hydrodynamics. By tracing MOC rays backward to the end of each timestep, terminating rays at an optical-depth threshold, and deriving timestep-dependent multigroup LO coefficients, they achieve accurate, diffusion-limit-consistent solutions with potential speedups in several test problems. Across plane-parallel, Marshak-wave, crooked pipe, Mach 45 shock, and Coax configurations, the method demonstrates good agreement with analytic and semi-analytic solutions and competitive performance in regimes where optically thick regions dominate, though Coax shows current slower performance relative to IMC-DDMC and S$_n$ in the present implementation. The work highlights the viability of backward MOC with LO closures for radiation transport and hydrodynamics, and points to future enhancements such as angular adaptivity and tensor-train compression to further improve efficiency and accuracy.

Abstract

Thermal radiation transport is a challenging problem in computational physics that has long been approached primarily in one of a few standard ways: approximate moment methods (for instance P$_1$ or M$_1$), implicit Monte Carlo, discrete ordinates, and long characteristics. In this work we consider the efficacy of the Method of (Long) Characteristics (MOC) applied to thermal radiation transport. Along the way we develop three major ideas: transporting MOC particles backwards in time from quadrature grids at the end of the timestep, limiting the computational cost of these backward characteristics by terminating transport once optical depths along rays become sufficiently large, and timestep-dependent closures with multigroup MOC solutions for a gray low-order system. We apply this method to a suite of standard radiation transport and radiation hydrodynamics test problems. We compare the method to several standard analytic and semi-analytic solutions, as well as implicit Monte Carlo, P$_1$, and discrete ordinates (S$_n$). We see that the method: gives excellent agreement with known results, has stability for large time steps, has the diffusion limit for large spatial cells, and achieves $\sim$20-70\% performance improvement when terminating optical depths at O(10-100) in the grey Marshak and crooked pipe problems. However, for the Coax radiation-hydrodynamics problem, we see that MOC is approximately two to three times slower than IMC-DDMC and S$_n$ in its current implementation.

Figures (17)

• Figure 1: A 2D mesh with MOC samples shown on their quadrature grid at the end of the cycle ($t_{\rm f}$; left) and beginning of cycle ($t_0$; right) after transport backwards in time. Directions of arrows indicate the movement of the particle going forwards in time, the opposite direction to which particles are actually moved during numerical updates. At the end of the cycle, the samples are colocated at a vertex and lie on a prescribed quadrature grid, allowing for fast and accurate integrations over solid angle.
• Figure 2: Illustration of MOC particle intensities being interpolated from previous end-of-timestep support point particles (cyan) with the same directions at vertices to current beginning-of-timestep position to which current MOC particle (red) was transported backwards in time to. Note that this algorithm generalizes to adaptive mesh refinement.
• Figure 3: Diagram of AMR and DLB for HO and LO systems in a spatial cell undergoing one level of refinement (dashed gray lines) and the upper right cell being sent to another MPI rank (hatched area). The blue solid nodes and arrows represent the original vertex-based MOC/HO intensity and face-based LO normal flux, respectively. The red nodes and arrows represent the new MOC intensity particles and interior normal fluxes, respectively, obtained from interpolation. The blue dotted arrows represent the partitioned LO normal flux from the original cell, determined from old-new face area ratios (in a first pass) and continuity corrections (in a second pass). To evaluate refinement in the hatched area, the original blue arrow and node data is needed on the rank containing the hatched area.
• Figure 4: Radiation energy density versus position for the plane-parallel vacuum problem at $t=0.003$ sh. The analytic solution (black dashed) is compared to the MOC (blue) and LO (orange) solutions. Left panel: triangular Chebyshev-Legendre discrete ordinate quadrature of order 6. Right panel: the same quadrature but of order 24. The LO solution is not strictly causal. We also observe that the LO solution is strongly influenced by ray effects from the MOC sector.
• Figure 5: Eddington tensor components versus position for the plane-parallel vacuum problem at $t=0.003$ sh. The analytic solution (black dashed) is compared to MOC solutions with triangular Chebyshev-Legendre quadratures of order 6 (blue) and order 24 (orange). Left panel: the diagonal $xx$-component. Right panel: the diagonal $yy$-component.