First-order System Least-squares Research Articles

Improved rates for a space–time FOSLS of parabolic PDEs

We consider the first-order system space–time formulation of the heat equation introduced by Bochev and Gunzburger (in: Bochev and Gunzburger (eds) Applied mathematical sciences, vol 166, Springer, New York, 2009), and analyzed by Führer and Karkulik (Comput Math Appl 92:27–36, 2021) and Gantner and Stevenson (ESAIM Math Model Numer Anal 55(1):283–299 2021), with solution components (u1,u2)=(u,-∇xu)\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$(u_1,\ extbf{u}_2)=(u,-\ abla _\ extbf{x} u)$$\\end{document}. The corresponding operator is boundedly invertible between a Hilbert space U and a Cartesian product of L2\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$L_2$$\\end{document}-type spaces, which facilitates easy first-order system least-squares (FOSLS) discretizations. Besides L2\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$L_2$$\\end{document}-norms of ∇xu1\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\ abla _\ extbf{x} u_1$$\\end{document} and u2\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\ extbf{u}_2$$\\end{document}, the (graph) norm of U contains the L2\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$L_2$$\\end{document}-norm of ∂tu1+divxu2\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\\partial _t u_1 +{{\\,\ extrm{div}\\,}}_\ extbf{x} \ extbf{u}_2$$\\end{document}. When applying standard finite elements w.r.t. simplicial partitions of the space–time cylinder, estimates of the approximation error w.r.t. the latter norm require higher-order smoothness of u2\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\ extbf{u}_2$$\\end{document}. In experiments for both uniform and adaptively refined partitions, this manifested itself in disappointingly low convergence rates for non-smooth solutions u. In this paper, we construct finite element spaces w.r.t. prismatic partitions. They come with a quasi-interpolant that satisfies a near commuting diagram in the sense that, apart from some harmless term, the aforementioned error depends exclusively on the smoothness of ∂tu1+divxu2\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$\\partial _t u_1 +{{\\,\ extrm{div}\\,}}_\ extbf{x} \ extbf{u}_2$$\\end{document}, i.e., of the forcing term f=(∂t-Δx)u\\documentclass[12pt]{minimal} \\usepackage{amsmath} \\usepackage{wasysym} \\usepackage{amsfonts} \\usepackage{amssymb} \\usepackage{amsbsy} \\usepackage{mathrsfs} \\usepackage{upgreek} \\setlength{\\oddsidemargin}{-69pt} \\begin{document}$$f=(\\partial _t-\\Delta _x)u$$\\end{document}. Numerical results show significantly improved convergence rates.

APFOS-Net: Asymptotic preserving scheme for anisotropic elliptic equations with deep neural network

In this paper, a new asymptotic preserving (AP) scheme is proposed for the anisotropic elliptic equations. Different from previous AP schemes, the actual one is based on first-order system least-squares for second-order partial differential equations, and it is uniformly well-posed with respect to anisotropic strength. The numerical computation is realized by a deep neural network (DNN), where least-squares functionals are employed as loss functions to determine parameters of DNN. Numerical results show that the current AP scheme is easy for implementation and is robust to approximate solutions or to identify the anisotropic parameter in various 2D and 3D tests.

Three-dimensional inversion for sparse potential data using first-order system least squares with application to gravity anomalies in Western Queensland

SUMMARY We present an inversion algorithm tailored for point gravity data. As the data are from multiple surveys, it is inconsistent with regards to spacing and accuracy. An algorithm design objective is the exact placement of gravity observations to ensure no interpolation of the data is needed prior to any inversion. This is accommodated by discretization using an unstructured tetrahedral finite-element mesh for both gravity and density with mesh nodes located at all observation points and a first-order system least-squares (FOSLS) formulation for the gravity modelling equations. Regularization follows the Bayesian framework where we use a differential operator approximation of an exponential covariance kernel, avoiding the usual requirement of inverting large dense covariance matrices. Rather than using higher order basis functions with continuous derivatives across element faces, regularization is also implemented with a FOSLS formulation using vector-valued property function (density and its gradient). Minimization of the cost function, comprised of data misfit and regularization, is achieved via a Lagrange multiplier method with the minimum of the gravity FOSLS functional as a constraint. The Lagrange variations are combined into a single equation for the property function and solved using an integral form of the pre-conditioned conjugate gradient method (I-PCG). The diagonal entries of the regularization operator are used as the pre-conditioner to minimize computational costs and memory requirements. Discretization of the differential operators with the finite-element method (FEM) results in matrix systems that are solved with smoothed aggregation algebraic multigrid pre-conditioned conjugate gradient (AMG-PCG). After their initial setup, the AMG-PCG operators and coarse grid solvers are reused in each iteration step, further reducing computation time. The algorithm is tested on data from 23 surveys with a total of 6519 observation points in the Mt Isa–Cloncurry region in north–west Queensland, Australia. The mesh had about 2.5 million vertices and 16.5 million cells. A synthetic case was also tested using the same mesh and error measures for localized concentrations of high and low densities. The inversion results for different parameters are compared to each other as well as to lower order smoothing. Final inversion results are shown with and without depth weighting and compared to previous geological studies for the Mt Isa–Cloncurry region.

Further results on a space-time FOSLS formulation of parabolic PDEs

In [2019, Space-time least-squares finite elements for parabolic equations, arXiv:1911.01942] by Führer and Karkulik, well-posedness of a space-time First-Order System Least-Squares formulation of the heat equation was proven. In the present work, this result is generalized to general second order parabolic PDEs with possibly inhomogenoeus boundary conditions, and plain convergence of a standard adaptive finite element method driven by the least-squares estimator is demonstrated. The proof of the latter easily extends to a large class of least-squares formulations.

Numerical results for adaptive (negative norm) constrained first order system least squares formulations

We perform a followup computational study of the recently proposed space–time first order system least squares ( FOSLS ) method subject to constraints referred to as CFOSLS where we now combine it with the new capability we have developed, namely, parallel adaptive mesh refinement (AMR) in 4D. The AMR is needed to alleviate the high memory demand in the combined space time domain and also allows general (4D) meshes that better follow the physics in space–time. With an extensive set of computational experiments, performed in parallel, we demonstrate the feasibility of the combined space–time AMR approach in both two space plus time and three space plus time dimensions.

Deep least-squares methods: An unsupervised learning-based numerical method for solving elliptic PDEs

This paper studies an unsupervised deep learning-based numerical approach for solving partial differential equations (PDEs). The approach makes use of the deep neural network to approximate solutions of PDEs through the compositional construction and employs least-squares functionals as loss functions to determine parameters of the deep neural network. There are various least-squares functionals for a partial differential equation. This paper focuses on the so-called first-order system least-squares (FOSLS) functional studied in [3], which is based on a first-order system of scalar second-order elliptic PDEs. Numerical results for second-order elliptic PDEs in one dimension are presented.

Adaptive First-Order System Least-Squares Finite Element Methods for Second-Order Elliptic Equations in Nondivergence Form

This paper studies adaptive first-order system least-squares finite element methods (LSFEMs) for second-order elliptic partial differential equations in nondivergence form. Unlike the classical finite element methods which use weak formulations of PDEs that are not applicable for the nondivergence equation, the first-order least-squares formulations naturally have stable weak forms without using integration by parts, allow simple finite element approximation spaces, and have built-in a posteriori error estimators for adaptive mesh refinements. The nondivergence equation is first written as a system of first-order equations by introducing the gradient as a new variable. Then, two versions of least-squares finite element methods using simple $C^0$ finite elements are developed in the paper. One is the $L^2$-LSFEM which uses linear elements, and the other is the W-LSFEM with a mesh-dependent weight to ensure the optimal convergence. Under the assumption that the PDE has a unique solution, a priori and a posteriori error estimates are presented. With an extra assumption on the operator regularity, convergence in standard norms for the W-LSFEM is also discussed. $L^2$-error estimates are derived for both formulations. Extensive numerical experiments for continuous, discontinuous, and even degenerate coefficients on smooth and singular solutions are performed to test the accuracy and efficiency of the proposed methods.

An optimal adaptive tensor product wavelet solver of a space-time FOSLS formulation of parabolic evolution problems

In this work, we construct a well-posed first-order system least squares (FOSLS) simultaneously space-time formulation of parabolic PDEs. Using an adaptive wavelet solver, this problem is solved with the best possible rate in linear complexity. Thanks to the use of a basis that consists of tensor products of wavelets in space and time, this rate is equal to that when solving the corresponding stationary problem. Our findings are illustrated by numerical results.

A Fluidity-Based First-Order System Least-Squares Method for Ice Sheets

This paper describes two first-order system least-squares (FOSLS) formulations of a nonlinear Stokes flow model for glaciers and ice sheets. In Glen's law, the most commonly used constitutive equation for ice rheology, the ice viscosity becomes infinite as the velocity gradients (strain rates) approach zero, which typically occurs near the ice surface where deformation rates are low, or when the basal slip velocities are high. The computational difficulties associated with the infinite viscosity are often overcome by an arbitrary modification of Glen's law that bounds the maximum viscosity. In this paper, two FOSLS formulations (the viscosity formulation and fluidity formulation) are presented. The viscosity formulation is a FOSLS representation of the standard nonlinear Stokes problem. The new fluidity formulation exploits the fact that only the product of the viscosity and strain rate appears in the nonlinear Stokes problem, a quantity that, in fact, approaches zero as the strain rate goes to zero. The fluidity formulation is expressed in terms of a new set of variables and overcomes the problem of infinite viscosity. The new formulation is well posed, and linearizations are essentially $H^1$-elliptic around approximations for which the viscosity is bounded. A nested iteration (NI) Newton-FOSLS-AMG (algebraic multigrid) approach is used to solve both nonlinear Stokes problems, in which most of the iterations are performed on the coarsest grid. Both formulations demonstrate optimal finite element convergence. However, the fluidity formulation is more accurate. The fluidity formulation results in linear systems that are more efficiently solved by AMG.

A Comparison of Finite Element Spaces for $H$(div) Conforming First-Order System Least Squares

First-order system least squares (FOSLS) is a commonly used technique in a wide range of physical applications. FOSLS discretizations are straightforward to implement and offer many advantages over traditional Galerkin or saddle-point formulations. Often, these problems are formulated in $H(div)$ spaces and implemented using $H(div)$-conforming elements. These elements have fewer regularity assumptions than the more commonly used $H^1$-conforming elements and are therefore often believed to be more suited for less regular problems. This paper compares the approximation properties of the $H(div)$-conforming Raviart--Thomas and Brezzi--Douglas--Marini elements to $H^1$-conforming piecewise polynomials for $H(div)$ problem formulations. Furthermore, for each problem, a corresponding $H^1$-formulation is derived and compared. The traditional Poisson and Stokes problems are examined with less regular solutions addressed by adaptive refinement strategies. The results imply that for smooth problems or problems exhibiting a point singularity $H^1$-conforming piecewise polynomials are a viable, and sometimes more efficient, choice even if an $H(div)$-formulation is used. Smooth problems do additionally benefit by using an $H^1$-formulation instead of an $H(div)$-formulation.

A first order system least squares method for the Helmholtz equation

We present a first order system least squares (FOSLS) method for the Helmholtz equation at high wave number k, which always leads to a Hermitian positive definite algebraic system. By utilizing a non-trivial solution decomposition to the dual FOSLS problem which is quite different from that of the standard finite element methods, we give an error analysis to the hp-version of the FOSLS method where the dependence on the mesh size h, the approximation order p, and the wave number k is given explicitly. In particular, under some assumption of the boundary of the domain, the L2 norm error estimate of the scalar solution from the FOSLS method is shown to be quasi optimal under the condition that kh/p is sufficiently small and the polynomial degree p is at least O(logk). Numerical experiments are given to verify the theoretical results.

Div First-Order System LL* (FOSLL*) for Second-Order Elliptic Partial Differential Equations

The first-order system LL* (FOSLL*) approach for general second-order elliptic partial differential equations was proposed and analyzed in [Z. Cai et al., SIAM J. Numer. Anal., 39 (2001), pp. 1418--1445], in order to retain the full efficiency of the $L^2$ norm first-order system least-squares (FOSLS) approach while exhibiting the generality of the inverse-norm FOSLS approach. The FOSLL* approach of Cai et al. was applied to the div-curl system with added slack variables, and hence it is quite complicated. In this paper, we apply the FOSLL* approach to the div system and establish its well-posedness. For the corresponding finite element approximation, we obtain a quasi-optimal a priori error bound under the same regularity assumption as the standard Galerkin method, but without the restriction to sufficiently small mesh size. Unlike the FOSLS approach, the FOSLL* approach does not have a free a posteriori error estimator. We then propose an explicit residual error estimator and establish its reliability and efficiency bounds.

Newton-[formula omitted] method for the second-order semi-linear elliptic partial differential equations

Newton’s method with first-order system least squares (FOSLS) finite element method has been widely used to approximately solve a system of nonlinear partial differential equations (Adler et al., 2010 [9], Codd et al., 2003 [10], Manteuffel et al., 2006 [11]). In this paper, we propose to use the first order system LL∗ method to find a correction in each Newton’s iteration which provides an L2-approximation of the second-order semi-linear elliptic partial differential equations while the typical Newton-FOSLS method providesH1-approximations. The numerical tests have been conducted to validate the theory.

Error Analysis for Constrained First-Order System Least-Squares Finite-Element Methods

In this paper, a general error analysis is provided for finite-element discretizations of partial differential equations in a saddle-point form with divergence constraint. In particular, this extends upon the work of [J. H. Adler and P. S. Vassilevski, Springer Proc. Math. Statist. 45, Springer, New York, 2013, pp. 1--19], giving a general error estimate for finite-element problems augmented with a divergence constraint and showing that these estimates are obtained for problems such as diffusion and Stokes' using the first-order system least-squares (FOSLS) finite-element method. The main result is that by enforcing the constraint on a ${\mathbf H}^1$-equivalent FOSLS formulation one maintains optimal convergence of the FOSLS functional (i.e., the energy norm of the error) while guaranteeing the conservation of the divergence constraint (i.e., mass conservation in some examples). The error estimates and results depend on using finite elements for the constraint space that are inf-sup stable when paired with the spaces used for the original unknowns. This includes using discontinuous spaces on coarse meshes and pairing with standard bilinear or biquadratic elements in order to confirm the results.

First-order system least squares with inhomogeneous boundary conditions

We derive well-posed first-order system least-squares formulations of second-order elliptic boundary value problems, in which homogeneous essential boundary conditions are appended by means of additional equations, rather than by their incorporation in the trial space. This approach has the advantage that it applies equally well to inhomogeneous boundary conditions.

The least-squares pseudo-spectral method for Navier–Stokes equations

A spectral collocation approximation of first-order system least squares for incompressible Stokes equations was analyzed in Kim et al. (2004) [12], and finite element approximations for incompressible Navier–Stokes equations were developed in Bochev et al. (1998,1999) [9,10]. The aim of this paper is to analyze the first-order system least-squares pseudo-spectral method for incompressible Navier–Stokes equations. The paper will be an extension of the result in Kim et al. (2004) [12] to the Navier–Stokes equations. Our least-squares functional is defined by the sum of discrete spectral norms of a first-order system of equations corresponding to the Navier–Stokes equations based on Legendre–Gauss–Lobatto points. We show its ellipticity and continuity over an appropriate product space, and spectral convergences of discretization errors are derived in the H1-norm and the L2-norm in each variable. Finally, we present some numerical examples.

Hybrid First-Order System Least Squares Finite Element Methods with Application to Stokes Equations

This paper combines first-order system least squares (FOSLS) with first-order system $LL^*$ (FOSLL$^*$) to create a Hybrid method. The FOSLS approach minimizes the error, ${\bf e}^h = {\bf u}^h - {\bf u}$, over a finite element subspace, ${\cal V}^h$, in the operator norm: $\min_{{\bf u}^h\in{\cal V}^h}\| L ({\bf u}^h-{\bf u})\|$. The FOSLL$^*$ method looks for an approximation in the range of $L^*$, setting ${\bf u}^h = L^*{\bf w}^h$ and choosing ${\bf w}^h \in {\cal W}^h$, a standard finite element space. FOSLL$^*$ minimizes the $L^2$ norm of the error over $L^*({\cal W}^h)$: $\min_{{\bf w}^h\in{\cal W}^h} \| L^*{\bf w}^h - {\bf u}\|$. FOSLS enjoys a locally sharp, globally reliable, and easily computable a posteriori error estimate, while FOSLL$^*$ does not. However, FOSLL$^*$ has the major advantage that it applies to problems that do not exhibit enough smoothness to enable the full advantages that the FOSLS approach otherwise provides. The Hybrid method attempts to retain the best properties of both FOSLS and FOSLL$^*$. This is accomplished by combining the FOSLS functional, the FOSLL$^*$ functional, and an intermediate term that draws them together. The Hybrid method produces an approximation, ${\bf u}^h$, that is nearly the optimal over ${\cal V}^h$ in the graph norm, $\|{\bf e}^h\|_{{\cal G}}^2:= \frac{1}{2}\|{\bf e}^h\|^2 + \|L{\bf e}^h\|^2$. The FOSLS and intermediate terms in the Hybrid functional provide a very effective a posteriori error measure. This paper establishes that the Hybrid functional is coercive and continuous in a graph-like norm with modest coercivity and continuity constants, $c_0 = 1/3$ and $c_1=3$; that both $\|{\bf e}^h \|$ and $\|L {\bf e}^h\|$ converge with rates based on standard interpolation bounds; and that if $LL^*$ has full $H^2$ regularity, the $L^2$ error, $\|{\bf e}^h\|$, converges with a full power of the discretization parameter, $h$, faster than the functional norm. Letting $\tilde{{\bf u}}^h$ denote the optimum over ${\cal V}^h$ in the graph norm, this paper also shows that if superposition is used with nested iteration, then $\| {\bf u}^h-\tilde{{\bf u}}^h\|_{{\cal G}}$ converges two powers of $h$ faster than the functional norm. Numerical tests are provided to confirm the efficiency of the Hybrid method and effectiveness of the a posteriori error measure.

Island Coalescence Using Parallel First-Order System Least Squares on Incompressible Resistive Magnetohydrodynamics

This paper investigates the performance of a parallel Newton, first-order system least-squares (FOSLS) finite-element method with local adaptive refinement and algebraic multigrid (AMG) applied to incompressible, resistive magnetohydrodynamics. In particular, an island coalescence test problem is studied that models magnetic reconnection using a reduced two-dimensional (2D) model of a tokamak fusion reactor. The results show that, using an appropriate temporal and spatial resolution, these methods are capable of resolving the physical instabilities accurately at small computational cost. The time-dependent, nonlinear system of PDEs is solved using work equivalent to about 50--60 simple relaxation sweeps (Gauss--Seidel iterations) per time step. Experiments show that, unless the time step is sufficiently small, nonphysical numerical instabilities may occur. Further, decreasing the time step size does not proportionally increase the cost of the computation, because AMG convergence is improved. In addition, an effective implementation of the methods in parallel keeps load balancing issues to a minimum. Various quantities, such as the reconnection rate and the “sloshing” effect of the plasma instability, are measured to confirm that the correct physics is reproduced.

First-order system least squares and the energetic variational approach for two-phase flow

This paper develops a first-order system least-squares (FOSLS) formulation for equations of two-phase flow. The main goal is to show that this discretization, along with numerical techniques such as nested iteration, algebraic multigrid, and adaptive local refinement, can be used to solve these types of complex fluid flow problems. In addition, from an energetic variational approach, it can be shown that an important quantity to preserve in a given simulation is the energy law. We discuss the energy law and inherent structure for two-phase flow using the Allen–Cahn interface model and indicate how it is related to other complex fluid models, such as magnetohydrodynamics. Finally, we show that, using the FOSLS framework, one can still satisfy the appropriate energy law globally while using well-known numerical techniques.

Efficiency Based Adaptive Local Refinement for First-Order System Least-Squares Formulations

In this paper, we propose new adaptive local refinement (ALR) strategies for first-order system least-squares finite elements in conjunction with algebraic multigrid methods in the context of nested iteration. The goal is to reach a certain error tolerance with the least amount of computational cost and nearly uniform distribution of the error over all elements. To accomplish this, the refinement decision at each refinement level is determined based on optimizing efficiency measures that take into account both error reduction and computational cost. Two efficiency measures are discussed: predicted error reduction and predicted computational cost. These methods are first applied to a two-dimensional (2D) Poisson problem with steep gradients, and the results are compared with the threshold-based methods described in [W. Dörfler, SIAM J. Numer. Anal., 33 (1996), pp. 1106–1124]. Next, these methods are applied to a 2D reduced model of the incompressible, resistive magnetohydrodynamic equations. These equations are used to simulate instabilities in a large aspect-ratio tokamak. We show that, by using the new ALR strategies on this system, we are able to resolve the physics using only 10 percent of the computational cost used to approximate the solutions on a uniformly refined mesh within the same error tolerance.