Realizability-preserving finite element discretizations of the $M_1$ model for dose calculation in proton therapy

Abstract

We present a deterministic framework for proton therapy dose calculation based on finite element discretizations of the energy-dependent $M_1$ moment model. The nonlinear $M_1$ system is derived from the Fokker--Planck equation for charged particles and closed using an entropy-based approximation of the second moment. Energy is treated as a pseudo-time coordinate. The zeroth and first moments of the proton fluence are evolved backward in energy. To ensure hyperbolicity and physical admissibility, we employ a monolithic convex limiting (MCL) strategy. Representing the standard continuous Galerkin discretization in terms of auxiliary `bar' states, we construct a nonlinear scheme that is provably invariant domain preserving (IDP) w.r.t. convex realizable sets consisting of all admissible states. The realizability of the bar states is enforced using the MCL technology for homogeneous hyperbolic systems. The forcing induced by stiff scattering is incorporated using Strang-type operator splitting. We use an explicit strong-stability-preserving Runge--Kutta method for the radiation transport subproblem and exact integration in the forcing steps, which guarantees the IDP property. The deposited dose is defined as the integral of a weighted zeroth moment over a bounded energy range. It is accumulated during the backward-in-energy evolution. Numerical experiments demonstrate that the proposed Strang-MCL method produces accurate and physically consistent dose distributions.

Figures (5)

• Figure 1: One-dimensional dose distribution of a $62\,\mathrm{MeV}$ proton beam in a water phantom computed with the proposed MCL scheme on meshes with $N_h\in\{257, 513, 1025, 2049\}$ and $\mathrm{CFL}=0.5$. Scattering is neglected. Results are compared with the reference solution \ref{['eq:refSol']}.
• Figure 2: Three-dimensional dose distribution of a $62\,\mathrm{MeV}$ proton beam in a water phantom computed with the proposed realizability-preserving MCL scheme on a uniform hexahedral mesh with $N_h = 257 \times 97 \times 97$ nodes and $\mathrm{CFL}=0.5$. Shown are slices at $z = 0.75\,\mathrm{cm}$ without scattering (top) and with scattering (bottom).
• Figure 3: Three-dimensional dose distributions of a $62\,\mathrm{MeV}$ proton beam in a water phantom computed with the proposed MCL scheme on a mesh with $N_h = 257\times 97 \times 97$ nodes and $\mathrm{CFL}=0.5$ with and without scattering effects. Results are integrated over $y$-$z$ planes and compared with the reference solution \ref{['eq:refSol']}.
• Figure 4: Three-dimensional dose distribution of a $65\,\mathrm{MeV}$ proton beam in a heterogeneous multi-material geometry computed with the proposed MCL scheme on a uniform hexahedral mesh with $N_h = 257 \times 97 \times 97$ nodes and $\mathrm{CFL}=0.5$. Shown are slices at $z = 0.75\,\mathrm{cm}$ with dose (top) and deposited energy density (bottom).
• Figure 5: Dose distribution for two perpendicular proton beams in a water phantom computed with the proposed MCL scheme on a uniform rectangular mesh with $N_h = 257\times 257$ nodes and $\mathrm{CFL}= 0.5$.

Theorems & Definitions (2)

• remark 1
• remark 2