A high-order accurate unstructured finite volume Newton–Krylov algorithm for inviscid compressible flows

Abstract

A fast implicit Newton–Krylov finite volume algorithm has been developed for high-order unstructured steady-state computation of inviscid compressible flows. The matrix-free generalized minimal residual (GMRES) algorithm is used for solving the linear system arising from implicit discretization of the governing equations, avoiding expensive and complex explicit computation of the high-order Jacobian matrix. The solution process has been divided into two phases: start-up and Newton iterations. In the start-up phase an approximate solution with the general characteristics of the steady-state flow is computed by using a defect correction procedure. At the end of the start-up phase, the linearization of the flow field is accurate enough for steady-state solution, and a quasi-Newton method is used, with an infinite time step and very rapid convergence. A proper limiter implementation for efficient convergence of the high-order discretization is discussed and a new formula for limiting the high-order terms of the reconstruction polynomial is introduced. The accuracy, fast convergence and robustness of the proposed high-order unstructured Newton–Krylov solver for different speed regimes is demonstrated for the second, third and fourth-order discretization. The possibility of reducing computational cost required for a given level of accuracy by using high-order discretization is examined.

Introduction

In computational aerodynamics, as in scientific computing generally, we wish to obtain an accurate solution as quickly as possible. When accuracy requirements are strict, high-order methods are preferable to second-order methods in the sense of providing a better solution on a given mesh, and therefore a comparable quality solution on a coarser mesh, when compared with second-order schemes. However, because of the higher computational cost per control volume and the possibility of degradation in convergence rates for high-order methods, the question of whether one can obtain an accurate solution more quickly with a high-order method hinges on the ability to achieve rapid convergence with such schemes. In this paper, we will present evidence that a high-order finite-volume scheme can in fact reach a given level of accuracy with less computational effort than a second-order scheme. Specifically, we will present a high-order accurate unstructured mesh finite-volume solver for compressible, inviscid flow as well as the techniques we use to ensure rapid convergence of an implicit Newton-GMRES scheme for both third- and fourth-order discretizations.

High-order finite-volume methods – by which we mean anything beyond second-order accurate – were first applied to aerodynamic flows by Rogers et al. [30], with a third-order discretization implemented in INS3D, and by Godfrey et al. [14], who applied ENO schemes in computational aerodynamics. More recently, De Rango and Zingg [8] dramatically reduced numerical error in drag using relatively coarse structured grids for steady turbulent flow over a 2D airfoil by applying a globally third-order accurate algorithm. Their results provide a convincing demonstration of the accuracy benefits of high-order methods compared with a second-order method for practical flows. Zingg et al. [36] compared different high-order accurate flux discretization techniques for laminar and turbulent flows (including transition) in both subsonic and transonic speed regimes, showing that high-order discretization produces comparably accurate solutions much more efficiently than a second-order method.

Our research in high-order unstructured solvers is motivated by the desire to combine the accuracy and efficiency benefits seen in the application of high-order methods on structured meshes with the geometric and adaptive flexibility of unstructured meshes. Also, to provide a natural upgrade path for the substantial base of existing unstructured finite-volume flow solvers, we focus on high-order finite-volume schemes rather than high-order discontinuous Galerkin methods.

The roots of high-order unstructured mesh finite-volume methods lie in the work of Barth and various co-workers, who published a series of papers describing techniques for using least-squares reconstruction to obtain high-order accuracy [5], [3], [4], [6]. While the techniques described in these papers are applicable in principle to all orders of accuracy, the results presented in them do not extend beyond third-order accuracy (quadratic reconstruction). Ollivier-Gooch [24], [25] experimented with solution-dependent weightings as a way to compute oscillation-free least-squares reconstructions for schemes of up to fourth-order accuracy. Delanaye and Essers [10] and Geuzaine et al. [13] proposed a quadratic reconstruction finite volume scheme, including a new approach to monotonicity enforcement. They computed the inviscid flux directly from their quadratic reconstruction; however, viscous terms were obtained through a linear interpolation and were therefore only second order. For monotonicity enforcement, they used a discontinuity detector to eliminate the high-order terms from the reconstruction in the vicinity of discontinuities.

Achieving steady-state convergence by time marching, or pseudo-time marching, requires that all transient phenomena in the domain be damped or convected out of the domain. Implicit methods use large time steps to both smear and rapidly propagate errors; convergence rates are typically dictated by non-linear stability issues or by inexact linearization of a non-linear residual. In the limit of very large time steps, implicit time advance schemes reduce to Newton’s method [34], [32], [2], [6]. In principle, Newton’s method gives quadratic convergence rates for steady flow solution when starting from a good enough initial guess. In practice, however, forming and solving the linear system at each Newton iteration accurately enough to achieve quadratic convergence is difficult enough that most researchers choose a quasi-Newton method. Quasi-Newton methods sacrifice quadratic convergence in favor of lower memory usage and lower cost per iteration, with the net result of reducing overall CPU time [27]. Quasi-Newton methods are generally categorized as approximate Newton (in which an approximate linearization is used) and inexact Newton (in which the linear system is not solved exactly).

Regardless of the details of the space and time discretization techniques, implicit methods inevitably require the solution of a large linear system resulting from the linearization of the fluid flow equations in time. While approximate factorization techniques are very effective linear solvers for structured meshes, modern unstructured mesh implicit schemes in CFD generally rely on iterative linear solvers, especially the generalized minimal residual (GMRES) method [31], which was developed particularly for non-symmetric systems such as those arising from implicit time discretization on unstructured meshes. In the matrix-free variant of GMRES, the matrix-vector products required by the GMRES algorithm are computed by using directional derivatives, eliminating the problem of explicitly forming the Jacobian matrix, and thereby considerably reducing both memory usage and programming effort. This is especially helpful for high-order unstructured mesh solvers where the full (analytic) Jacobian calculation is costly and difficult.

The efficiency of GMRES in solving a linear system depends strongly on the conditioning of the linear system, and compressible flow problems on unstructured meshes do not typically have well-conditioned Jacobian matrices. This is especially true for high-order methods, as we shall see in Section 5. Effective preconditioning is therefore essential to rapid convergence of Newton-GMRES schemes.

To retain the low memory usage of a matrix-free GMRES method, we could use a matrix-free preconditioner, such as the lower–upper symmetric Gauss–Seidel (LU-SGS) preconditioner introduced by Luo et al. for 3D compressible flows [16]. They completely eliminated the storage of the preconditioning Jacobian matrix by approximating the Jacobian with numerical fluxes. However, they did not perform full Newton iterations, probably because of the stability considerations for their LU-SGS preconditioner, and their convergence rates remained nearly linear. We have experimented with LU-SGS preconditioning for a high-order matrix-free Newton-GMRES algorithm for compressible flow [20]. For supersonic flows, LU-SGS-GMRES was about as efficient for the third-order discretization as for the second-order one. However, for flows that are more difficult to converge, for fourth-order discretization, and for full Newton iterations (infinite time step), the LU-SGS preconditioner proved to be inadequate for our needs.

Currently, the most prevalent preconditioning approach for compressible flows is to use a low-order Jacobian as the preconditioner matrix and incomplete lower–upper (ILU) factorization to apply the preconditioner. For structured meshes, Pueyo and Zingg [28] presented an efficient matrix-free Newton–GMRES solver for steady-state aerodynamic flow computations using this approach. After a thorough parametric study for inviscid, laminar and turbulent 2D flows, they concluded that the approximate Newton method using matrix-free GMRES preconditioned with a first-order Jacobian and ILU(2) provides the best efficiency in terms of CPU time for most cases. Subsequently, Nichols and Zingg [23] developed a 3D multi-block Newton–Krylov solver for the Euler equations using the same approach, concluding that for this case, ILU(1) gives an appropriate balance between good preconditioning and low computational time per iteration.

For unstructured meshes, Venkatakrishnan and Mavriplis [35] developed an approximate Newton-GMRES implicit solver for computing compressible inviscid and turbulent flows around a multi-element airfoil. They compared different preconditioning strategies and found out that GMRES with ILU(0) preconditioning and a first-order Jacobian had the best performance. Manzano et al. [17] presented an efficient ILU preconditioned matrix-free Newton-GMRES algorithm for 3D unstructured meshes. They used different levels of fill (ILU(1–3)) depending on the case and the flux residual to achieve optimum performance. Delanaye et al. [12] presented the first-ever third-order accurate ILU preconditioned matrix-free Newton–GMRES solver for the Euler and Navier–Stokes equations. They showed that convergence can stall for stiff problems when using a high-order Jacobian and ILU(0) decomposition. Full convergence was achieved by allowing additional fill (ILU(1)) in the decomposition.

In this paper, we will present our recent work on third- and fourth-order accurate schemes for the Euler equations. In addition to verifying full-order accuracy in our solutions, we have also attained excellent efficiency in convergence to steady state, with the third- and fourth-order schemes comparable in efficiency to the second-order scheme. Our discretization scheme uses Barth and Frederickson’s [5] reconstruction algorithm as a foundation for a fourth-order accurate solution of the compressible Euler equations. The details of our spatial discretization scheme are given in Section 3, including high-order boundary condition enforcement and a better-converging variant of Delanaye and Essers’ [10] monotonicity enforcement scheme. Section 4 gives a comprehensive description of our time advance scheme, including formation of the preconditioning matrix and our strategies for both start-up and Newton iterations; the complex issue of preconditioning, which is only summarized in this section, is discussed fully in a companion paper [22]. The accuracy, efficiency and robustness of the current algorithm are demonstrated in Section 5 through several test cases.