Preconditioned iterative solvers for constrained high-order implicit shock tracking methods

TL;DR

The paper tackles the challenge of solving linearized, saddle-point systems arising in constrained high-order implicit shock tracking (HOIST) methods. It develops a family of matrix-based preconditioners that blend constrained preconditioning with standard DG preconditioners (block Jacobi and block ILU MDF) and avoids forming dense Hessian blocks, while also introducing a two-level $p$-multigrid strategy. Extensive numerical experiments on inviscid Euler problems show that ILU-based, BILU/ilu preconditioners offer the best practical performance, with the Hessian regularization parameter $\\gamma$ dominating iteration counts and $p$-multigrid providing limited gains. The results support applying these preconditioners to large-scale HOIST problems and motivate future work on parallel scalability and viscous extensions.

Abstract

High-order implicit shock tracking (fitting) is a class of high-order numerical methods that use numerical optimization to simultaneously compute a high-order approximation to a conservation law solution and align elements of the computational mesh with non-smooth features. This alignment ensures that non-smooth features are perfectly represented by inter-element jumps and high-order basis functions approximate smooth regions of the solution without nonlinear stabilization, which leads to accurate approximations on traditionally coarse meshes. In this work, we devise a family of preconditioners for the saddle point linear system that defines the step toward optimality at each iteration of the optimization solver so Krylov solvers can be effectively used. Our preconditioners integrate standard preconditioners from constrained optimization with popular preconditioners for discontinuous Galerkin discretizations such as block Jacobi, block incomplete LU factorizations with minimum discarded fill reordering, and p-multigrid. Thorough studies are performed using two inviscid compressible flow problems to evaluate the effectivity of each preconditioner in this family and their sensitivity to critical shock tracking parameters such as the mesh and Hessian regularization, linearization state, and resolution of the solution space.

Figures (16)

• Figure 1: Example two-dimensional mesh (left) ($10$ nodes and $9$ elements) and corresponding sparsity structure of ${\bm{J}}_{\bm{u}}$ (right) for a polynomial degree of $p=1$ and a single conservation law ($m=1$). This choice leads to $9$ blocks of size $3 \times 3$ for ${\bm{J}}_{\bm{u}}$.
• Figure 2: Sparsity structure of $(\partial{\bm{R}}/ \partial{\bm{u}})^T$ for the mesh in Figure \ref{['fig:Ju_sprs']}, polynomial degrees of $p=1,~p'=2$ and a single conservation law ($m=1$). This choice results in $9$ blocks of size $6 \times 3$ for $\partial{\bm{R}}/\partial{\bm{u}}$.
• Figure 3: Sparsity structure of ${\bm{B}}_{{\bm{u}}{\bm{u}}}$ (left) and ${\bm{B}}_{{\bm{y}}{\bm{y}}}$ (right) (assuming no boundary constraints, i.e. ${\boldsymbol{\phi}}({\bm{y}})={\bm{y}}$) for mesh depicted in Figure \ref{['fig:Ju_sprs']} with polynomial degrees of $p=1,~p'=2$ and a single conservation law ($m=1$).
• Figure 4: Sparsity structure of $(\partial{\bm{R}} /{\partial{\bm{x}}})^T$ for mesh depicted in Figure \ref{['fig:Ju_sprs']} with polynomial degrees $p=1,~p'=2$ and a single conservation law ($m=1$).
• Figure 5: Example of mesh restriction/prolongation for a second order mesh ($q=2$) with one element (left). The original element is restricted to $q=1$ (middle) removing the high order nodes $4,5,6$. Prolongation (right) is performed by inserting high order nodes interpolating the low order element.

Theorems & Definitions (6)

• Remark 1
• Remark 2
• Remark 3
• Remark 4
• Remark 5
• Remark 6