A geometric multigrid solution technique for the incompressible Navier-Stokes equations in three dimensions is presented, utilizing the concept of discretely divergence-free finite elements without requiring the explicit construction of a basis on each mesh level. For this purpose, functions are constructed in an a priori manner spanning the subspace of discretely divergence-free functions for the Rannacher-Turek finite element pair under consideration. Compared to mixed formulations, this approach yields smaller system matrices with no saddle point structure. This prevents the use of complex Schur complement solution techniques, and more general preconditioners can be employed. While constructing a basis for discretely divergence-free finite elements may pose significant challenges and prevent the use of a structured assembly routine, a basis is utilized only on the coarsest mesh level of the multigrid algorithm. On finer grids, this information is extrapolated to prescribe boundary conditions efficiently. Here, special attention is required for geometries introducing bifurcations in the flow. In such cases, so-called "global" functions with an extended support are defined, which can be used to prescribe the net flux through different branches. Various numerical examples for meshes with different shapes and boundary conditions illustrate the strengths, limitations, and future challenges of this solution concept.

We present the application of a time-simultaneous multigrid algorithm closely related to multigrid waveform relaxation for stabilized convection-diffusion equations in the regime of small diffusion coefficients. We use Galerkin finite elements and the Crank-Nicolson scheme for discretization in space and time. The multigrid method blocks all time steps for each spatial unknown, enhancing parallelization in space. While the number of iterations of the solver is bounded above for the 1D heat equation, convergence issues arise in convection-dominated cases. In singularly perturbed advection-diffusion scenarios, Galerkin FE discretizations are known to show instabilities in the numerical solution. We explore a higher-order variational multiscale stabilization, aiming to enhance solution smoothness and improve convergence without compromising accuracy.

This work focuses on the solution of the convection–diffusion equation, especially for small diffusion coefficients, employing a time-simultaneous multigrid algorithm, which is closely related to multigrid waveform relaxation. For discretization purposes, linear finite elements are used while the Crank–Nicolson scheme acts as the time integrator. By combining all time steps into a global linear system of equations and rearranging the degrees of freedom, a space-only problem is formed with vector-valued unknowns for each spatial node. The generalized minimal residual method with block Jacobi preconditioning can be used to numerically solve the (spatial) problem, allowing a higher degree of parallelization in space. A time-simultaneous multigrid approach is applied, utilizing space-only coarsening and the aforementioned solution techniques for smoothing purposes. Numerical studies analyze the iterative solution technique for 1D test problems. For the heat equation, the number of iterations stays bounded independently of the number of time steps, the time increment, and the spatial resolution. However, convergence issues arise in situations where the diffusion coefficient is small compared to the grid size and the magnitude of the velocity field. Therefore, a higher-order variational multiscale stabilization is used to improve the convergence behavior and solution smoothness without compromising its accuracy in convection-dominated scenarios.

This work is focused on an accelerated global-in-time solution strategy for the Oseen equations, which highly exploits the augmented Lagrangian methodology to improve the convergence behavior of the Schur complement iteration. The main idea of the solution strategy is to block the individual linear systems of equations at each time step into a single all-at-once saddle point problem. By elimination of all velocity unknowns, the resulting implicitly defined equation can then be solved using a global-in-time pressure Schur complement (PSC) iteration. To accelerate the convergence behavior of this iterative scheme, the augmented Lagrangian approach is exploited by modifying the momentum equation for all time steps in a strongly consistent manner. While the introduced discrete grad-div stabilization does not modify the solution of the discretized Oseen equations, the quality of customized PSC preconditioners drastically improves and, hence, guarantees a rapid convergence. This strategy comes at the cost that the involved auxiliary problem for the velocity field becomes ill conditioned so that standard iterative solution strategies are no longer efficient. Therefore, a highly specialized multigrid solver based on modified intergrid transfer operators and an additive block preconditioner is extended to solution of the all-at-once problem. The potential of the proposed overall solution strategy is discussed in several numerical studies as they occur in commonly used linearization techniques for the incompressible Navier–Stokes equations.

In this work, a new global-in-time solution strategy for incompressible flow problems is presented, which highly exploits the pressure Schur complement (PSC) approach for the construction of a space–time multigrid algorithm. For linear problems like the incompressible Stokes equations discretized in space using an inf-sup-stable finite element pair, the fundamental idea is to block the linear systems of equations associated with individual time steps into a single all-at-once saddle point problem for all velocity and pressure unknowns. Then the pressure Schur complement can be used to eliminate the velocity fields and set up an implicitly defined linear system for all pressure variables only. This algebraic manipulation allows the construction of parallel-in-time preconditioners for the corresponding all-at-once Picard iteration by extending frequently used sequential PSC preconditioners in a straightforward manner. For the construction of efficient solution strategies, the so defined preconditioners are employed in a GMRES method and then embedded as a smoother into a space–time multigrid algorithm, where the computational complexity of the coarse grid problem highly depends on the coarsening strategy in space and/or time. While commonly used finite element intergrid transfer operators are used in space, tailor-made prolongation and restriction matrices in time are required due to a special treatment of the pressure variable in the underlying time discretization. The so defined all-at-once multigrid solver is extended to the solution of the nonlinear Navier–Stokes equations by using Newton’s method for linearization of the global-in-time problem. In summary, the presented multigrid solution strategy only requires the efficient solution of time-dependent linear convection–diffusion–reaction equations and several independent Poisson-like problems. In order to demonstrate the potential of the proposed solution strategy for viscous fluid simulations with large time intervals, the convergence behavior is examined for various linear and nonlinear test cases.

Abstract In this work, the convergence behavior of a time-simultaneous two-grid algorithm for the one-dimensional heat equation is studied using Fourier arguments in space. The underlying linear system of equations is obtained by a finite element or finite difference approximation in space while the semi-discrete problem is discretized in time using the ϑ-scheme. The simultaneous treatment of all time instances leads to a global system of linear equations which provides the potential for a higher degree of parallelization of multigrid solvers due to the increased number of degrees of freedom per spatial unknown. It is shown that the all-at-once system based on an equidistant discretization in space and time stays well conditioned even if the number of blocked time-steps grows arbitrarily. Furthermore, mesh-independent convergence rates of the considered two-grid algorithm are proved by adopting classical Fourier arguments in space without assuming periodic boundary conditions. The rate of convergence with respect to the Euclidean norm does not deteriorate arbitrarily if the number of blocked time steps increases and, hence, underlines the potential of the solution algorithm under investigation. Numerical studies demonstrate why minimizing the spectral norm of the iteration matrix may be practically more relevant than improving the asymptotic rate of convergence.

Building on recent advances in the analysis and design of algebraic flux correction (AFC) schemes, new differentiable limiter functions are constructed for efficient nonlinear solution strategies. The proposed scaling parameters are used to limit artificial diffusion operators incorporated into the residual of a high order target scheme to produce accurate and bound-preserving finite element approximations to hyperbolic problems. Due to this stabilization procedure, the occurring system becomes highly nonlinear and the efficient computation of corresponding solutions is a challenging task. The presented regularization approach makes the AFC residual twice continuously differentiable so that Newton’s method converges quadratically for sufficiently good initial guesses. Furthermore, the performance of each nonlinear iteration is improved by expressing the Jacobian as the sum and product of matrices having the same sparsity pattern as the Galerkin system matrix. Eventually, the AFC methodology constructed is validated numerically by applying it to several numerical benchmarks.

We discuss how “parallel-in-space & simultaneous-in-time” Newton-multigrid approaches can be designed which improve the scaling behavior of the spatial parallelism by reducing the latency costs. The idea is to solve many time steps at once and therefore solving fewer but larger systems. These large systems are reordered and interpreted as a space-only problem leading to multigrid algorithm with semi-coarsening in space and line smoothing in time direction. The smoother is further improved by embedding it as a preconditioner in a Krylov subspace method. As a prototypical application, we concentrate on scalar partial differential equations (PDEs) with up to many thousands of time steps which are discretized in time, resp., space by finite difference, resp., finite element methods. For linear PDEs, the resulting method is closely related to multigrid waveform relaxation and its theoretical framework. In our parabolic test problems the numerical behavior of this multigrid approach is robust w.r.t. the spatial and temporal grid size and the number of simultaneously treated time steps. Moreover, we illustrate how corresponding time-simultaneous fixed-point and Newton-type solvers can be derived for nonlinear nonstationary problems that require the described solution of linearized problems in each outer nonlinear step. As the main result, we are able to generate much larger problem sizes to be treated by a large number of cores so that the combination of the robustly scaling multigrid solvers together with a larger degree of parallelism allows a faster solution procedure for nonstationary problems.

This work extends the algebraic flux correction (AFC) paradigm to finite element discretizations of conservation laws for symmetric tensor fields. The proposed algorithms are designed to enforce discrete maximum principles and preserve the eigenvalue range of evolving tensors. To that end, a continuous Galerkin approximation is modified by adding a linear artificial diffusion operator and a nonlinear antidiffusive correction. The latter is decomposed into edge-based fluxes and constrained to prevent violations of local bounds for the minimal and maximal eigenvalues. In contrast to the flux-corrected transport (FCT) algorithm developed previously by the author and existing slope limiting techniques for stress tensors, the admissible eigenvalue range is defined implicitly and the limited antidiffusive terms are incorporated into the residual of the nonlinear system. In addition to scalar limiters that use a common correction factor for all components of a tensor-valued antidiffusive flux, tensor limiters are designed using spectral decompositions. The new limiter functions are analyzed using tensorial extensions of the existing AFC theory for scalar convection-diffusion equations. The proposed methodology is backed by rigorous proofs of eigenvalue range preservation and Lipschitz continuity. Convergence of pseudo time-stepping methods to stationary solutions is demonstrated in numerical studies.

This book introduces a very general framework for the analysis and design of bound-preserving finite element methods. The main focus is on new limiting techniques designed to control the range of solution values for advected scalar quantities or the eigenvalue range of symmetric tensors.

This work extends the flux-corrected transport (FCT) methodology to arbitrary order continuous finite element discretizations of scalar conservation laws on simplex meshes. Using Bernstein polynomials as local basis functions, we constrain the total variation of the numerical solution by imposing local discrete maximum principles on the Bézier net. The design of accuracy-preserving FCT schemes for high order Bernstein–Bézier finite elements requires the development of new algorithms and/or generalization of limiting techniques tailored for linear and multilinear Lagrange elements. In this paper, we propose (i) a new discrete upwinding strategy leading to local extremum bounded low order approximations with compact stencils, (ii) high order variational stabilization based on the difference between two gradient approximations, and (iii) new localized limiting techniques for antidiffusive element contributions. The optional use of a smoothness indicator, based on a second derivative test, makes it possible to potentially avoid unnecessary limiting at smooth extrema and achieve optimal convergence rates for problems with smooth solutions. The accuracy of the proposed schemes is assessed in numerical studies for the linear transport equation in 1D and 2D.

This paper presents a new approach to constraining the eigenvalue range of symmetric tensors in numerical advection schemes based on the flux-corrected transport (FCT) algorithm and a continuous finite element discretization. In the context of element-based FEM-FCT schemes for scalar conservation laws, the numerical solution is evolved using local extremum diminishing (LED) antidiffusive corrections of a low order approximation which is assumed to satisfy the relevant inequality constraints. The application of a limiter to antidiffusive element contributions guarantees that the corrected solution remains bounded by the local maxima and minima of the low order predictor. The FCT algorithm to be presented in this paper guarantees the LED property for the maximal and minimal eigenvalues of the transported tensor at the low order evolution step. At the antidiffusive correction step, this property is preserved by limiting the antidiffusive element contributions to all components of the tensor in a synchronized manner. The definition of the element-based correction factors for FCT is based on perturbation bounds for auxiliary tensors which are constrained to be positive semidefinite to enforce the generalized LED condition. The derivation of sharp bounds involves calculating the roots of polynomials of degree up to 3. As inexpensive and numerically stable alternatives, limiting techniques based on appropriate estimates are considered. The ability of the new limiters to enforce local bounds for the eigenvalue range is confirmed by numerical results for 2D advection problems. (C) 2017 Elsevier Inc. All rights reserved.

This paper deals with the problem of positivity preservation in numerical algorithms for simulating fiber suspension flows. In contrast to fiber orientation models based on the Advani–Tucker evolution equations for even-order orientation tensors, the probability distribution function of fiber orientation is approximated using the Galerkin discretization of the Fokker–Planck equation with Fourier basis functions or spherical harmonics. This procedure leads to a natural generalization of orientation tensor models replacing ad hoc closure approximations by Galerkin equations for the fine-scale components. After introducing an operator splitting approach to solving the discretized Fokker–Planck equation, we present conditions and correction techniques that guarantee physically correct distribution functions. As the reader will see, the derivation of these conditions is independent of the space dimension and their applicability is not limited to the simulation of fiber suspensions.

This paper presents a new approach to synchronized limiting of density, energy, and pressure jumps in discontinuous Galerkin (DG) methods for the Euler equations of gas dynamics. A vertex-based version of the Barth-Jespersen limiter for scalar quantities is generalized to constrain the gradients of the conservative variables in a way which guarantees that all quantities of interest remain in the range of admissible values. The bounds for the corresponding inequality constraints are designed to enforce local maximum principles in regions of strong density variations and become less restrictive in smooth regions. The proposed limiting strategy guarantees positivity preservation and leads to closed-form expressions for the synchronized gradient correction factors without the need to solve inequality-constrained optimization problems. A numerical study is performed for two-dimensional test problems.

This work addresses the design of failsafe flux limiters for systems of conserved quantities and derived variables in numerical schemes for the equations of gas dynamics. Building on Zalesak's multidimensional flux-corrected transport (FCT) algorithm, we construct a new positivity-preserving limiter for the density, total energy, and pressure. The bounds for the underlying inequality constraints are designed to enforce local maximum principles in regions of strong density variations and become less restrictive in smooth regions. The proposed approach leads to closed-form expressions for the synchronized correction factors without the need to solve inequality-constrained optimization problems. A numerical study is performed for the compressible Euler equations discretized using a finite element based FCT scheme.