Provably Convergent and Robust Newton-Raphson Method: A New Dawn in Primitive Variable Recovery for Relativistic MHD

Abstract

A long-standing and formidable challenge faced by all conservative schemes for relativistic magnetohydrodynamics (RMHD) is the recovery of primitive variables from conservative ones. This process involves solving highly nonlinear equations subject to physical constraints. An ideal solver should be "robust, accurate, and fast -- it is at the heart of all conservative RMHD schemes," as emphasized in [S.C. Noble et al., ApJ, 641:626-637, 2006]. Despite over three decades of research, seeking efficient solvers that can provably guarantee stability and convergence remains an open problem. This paper presents the first theoretical analysis for designing a robust, physical-constraint-preserving (PCP), and provably (quadratically) convergent Newton-Raphson (NR) method for primitive variable recovery in RMHD. Our key innovation is a unified approach for the initial guess, devised based on sophisticated analysis. It ensures that the NR iteration consistently converges and adheres to physical constraints. Given the extreme nonlinearity and complexity of the iterative function, the theoretical analysis is highly nontrivial and technical. We discover a pivotal inequality for delineating the convexity and concavity of the iterative function and establish theories to guarantee the PCP property and convergence. We also develop theories to determine a computable initial guess within a theoretical "safe" interval. Intriguingly, we find that the unique positive root of a cubic polynomial always falls within this interval. Our PCP NR method is versatile and can be seamlessly integrated into any RMHD scheme that requires the recovery of primitive variables, potentially leading to a broad impact in this field. As an application, we incorporate it into a discontinuous Galerkin method, resulting in fully PCP schemes. Several numerical experiments demonstrate the efficiency and robustness of the PCP NR method.

Figures (7)

• Figure 1: A potential failure scenario in solving \ref{['equ:fu_exa']} using the NR method. When $\xi_1 \le \xi_a$, $f_a(\xi_1)\le 0$ so that the values of $\mathcal{W}(\xi_1)$ and ${\mathcal{F}}(\xi_1)$ become complex or ill-defined, resulting in failure. Here $\xi_a$ is the largest non-negative root of $f_a(\xi)$.
• Figure 1: Three cases of NR iterations.
• Figure 1: Two patterns of ${\mathcal{F}}'(\xi)$ observed in extensive random experiments for three EOSs.
• Figure 1: The profiles of ${\mathcal{F}}(\xi)$ with $\mathcal{W}(\xi)$ computed by \ref{['eq:DefWxi']} (left) and \ref{['fu:express2']} (right), respectively. Here we set $\gamma=5/3$, $D=1$, $E=10^{8}$, $m^2=9.999999999\times 10^{15}$, $B=10^4$, and $\tau =1$.
• Figure 2: Four possible patterns of $f(\xi)$ described in Lemma \ref{['lem:NRconverges']}.

Theorems & Definitions (32)

• Definition 2.1
• Theorem 2.2: wu2017admissiblewu2018physical
• Remark 2.3
• Remark 2.4
• Theorem 3.1: wu2017admissiblewu2018physical
• Definition 3.2
• Definition 3.3
• Theorem 3.4
• Proof 1
• Theorem 3.5