First-Order System Least Squares for the Stokes Equations, with Application to Linear Elasticity.
Following our earlier work on general second-order scalar equations, here we develop a least-squares functional for the two- and three-dimensional Stokes equations, generalized slightly by allowing a pressure term in the continuity equation. By introducing a velocity flux variable and associated curl and trace equations, we are able to establish ellipticity in an H1 product norm appropriately weighted by the Reynolds number. This immediately yields optimal discretization error estimates for finite element spaces in this norm and optimal algebraic convergence estimates for multiplicative and additive multigrid methods applied to the resulting discrete systems. Both estimates are naturally uniform in the Reynolds number. Moreover, our pressure-perturbed form of the generalized Stokes equations allows us to develop an analogous result for the Dirichlet problem for linear elasticity, where we obtain the more substantive result that the estimates  are uniform in the Poisson ratio.
Introduction
. In earlier work [10], [11], we developed least-squares functionals for
a first-order system formulation of general second-order elliptic scalar partial differential
equations. The functional developed in [11] was shown to be elliptic in the sense that its
homogeneous form applied to the (pressure and velocities) is equivalent to
the norm. This means that the individual variables in the functional are essentially
decoupled (more precisely, their interactions are essentially subdominant). This important
property ensures that standard finite element methods are of H 1 -optimal accuracy in each
variable, and that multiplicative and additive multigrid methods applied to the resulting
discrete equations are optimally convergent.
The purpose of this paper is to extend this methodology to the Stokes equations in two
and three dimensions. To this end, we begin by reformulating the Stokes equations as a
first-order system derived in terms of an additional vector variable, the velocity flux, defined
as the vector of gradients of the Stokes velocities. We first apply a least-squares principle
to this system using L 2 and H \Gamma1 norms weighted appropriately by the Reynolds number,
Re. We then show that the resulting functional is elliptic in a product norm involving Re
and the L 2 and H 1 norms. While of theoretical interest in its own right, we use this result
here primarily as a vehicle for establishing that a modified form of this functional is fully
elliptic in an H 1 product norm scaled by Re.
This appears to be the first general theory of this kind for the Stokes equations in
general dimensions with velocity boundary conditions. Bochev and Gunzburger [6] developed
least-squares functionals for Stokes equations in norms that include stronger Sobolev
terms and mesh weighting, but none are product H 1 elliptic. Chang [13] also used velocity
derivative variables to derive a product H 1 elliptic functional for Stokes equations, but it is
inherently limited to two dimensions. For general dimensions, a vorticity-velocity-pressure
Center for Applied Mathematical Sciences, Department of Mathematics, University of Southern Cali-
fornia, 1042 W. 36th Place, DRB 155, Los Angeles, CA 90089-1113. email: zcai@math.usc.edu
y Program in Applied Mathematics, Campus Box 526, University of Colorado at Boulder, Boulder,
CO 80309-0526. email: tmanteuf@boulder.colorado.edu and stevem@boulder.colorado.edu. This work was
sponsored by the Air Force Office of Scientific Research under grant number AFOSR-91-0156, the National
Science Foundation under grant number DMS-8704169, and the Department of Energy under grant number
DE-FG03-93ER25165.
form (cf.[4] and [20]) proved to be product H 1 elliptic, but only for certain nonstandard
boundary conditions. For the more practical (cf. [17], [22], and [25]) velocity boundary
conditions treated here, the velocity-vorticity-pressure formulation examined by Chang [14]
can be shown by counterexample [3] not to be equivalent to any H 1 product norm, even
with the added boundary condition on the normal component of vorticity. Moreover, this
formulation admits no apparent additional equation, such as the curl and trace constraints
introduced below for our formulation, that would enable such an equivalence. The velocity-
pressure-stress formulation described in [7] has the same shortcomings. (If the vorticity
and deformation stress variables are important, then they can be easily and accurately
reconstructed from the velocity-flux variables introduced in our formulation.)
While our least-squares form requires several new dependent variables, we believe that
the added cost is more than offset by the strengthened accuracy of the discretization and the
speed that the attendant multigrid solution process attains. Moreover, while our modified
functional requires strong regularity conditions, this is to be expected for obtaining full
product H 1 ellipticity in all variables, including velocity fluxes. (We thus obtain optimal
estimates for the derivatives of velocity.) In any case, strengthened regularity is not
necessary for the first functional we introduce.
Our modified Stokes functional is obtained essentially by augmenting the first-order
system with a curl constraint and a scalar (trace) equation involving certain derivatives of
the velocity flux variable, then appealing to a simple L 2 least-squares principle. As in [11] for
the scalar case, the important H 1 ellipticity property that we establish guarantees optimal
finite element accuracy and multigrid convergence rates applied to this Stokes least-squares
functional that are uniform in Re.
One of the more compelling benefits of least squares is the freedom to incorporate
additional equations and impose additional boundary conditions as long as the system is
consistent. In fact, many problems are perhaps best treated with overdetermined (but
first-order systems, as we have here for Stokes. We therefore abandon the so-called
ADN theory (cf. [1] [2]), which is restricted to square systems, in favor of more direct
tools of analysis.
An important aspect of our general formulation is that it applies equally well to the
Dirichlet problem for linear elasticity. This is done by posing the Stokes equations in a
slightly generalized form that includes a pressure term in the continuity equation. Our
development and results then automatically apply to linear elasticity. Most important, our
optimal discretization and solver estimates are uniform in the Lam'e constants.
We emphasize that the discretization and algebraic convergence properties for the generalized
Stokes equations are automatic consequences of the H 1 product norm ellipticity
established here and the finite element and multigrid theories established in Sections 3-5
of [11]. We are therefore content with an abbreviated paper that focuses on establishing
ellipticity, which we do in Section 3. Section 2 introduces the generalized Stokes equations,
the two relevant first-order systems and their functionals, and some preliminary theory.
Concluding remarks are made in Section 4.
2. The Stokes Problem, Its First-Order System Formulation, and Other
Preliminaries.
Let\Omega be a bounded, open, connected domain in ! n
Lipschitz boundary @
\Omega\Gamma The pressure-perturbed form of the generalized stationary Stokes
equations in dimensionless variables may be written as
\Gamma- \Delta u +r
r
where the symbols \Delta, r, and r\Delta stand for the Laplacian, gradient, and divergence operators,
- is the reciprocal of the Reynolds number Re; f is a given vector function; g
is a given scalar function; and ffi is some nonnegative constant
linear elasticity). Without loss of generality, we may assume that
Z
\Omega
Z
\Omega
equation (2.1) can have a solution only when g satisfies (2.2), and we are then
free to ask that p satisfy (2.2). For ffi ? 0, in general we have only that
R\Omega
but this can be reduced to (2.2) simply by replacing p by
and g by 0 in (2.1).)
We consider the (generalized) Stokes equations (2.1) together with the Dirichlet velocity
boundary condition
The slightly generalized Stokes equations in (2.1) allow our results to apply to linear
elasticity. In particular, consider the Dirichlet problem
where u now represents displacements and - and - are the (positive) Lam'e constants. By
here we mean the n-vector of components \Delta u i , that is, \Delta applies to u componentwise.
This is recast in form (2.1)-(2.2) by introducing the pressure variable 1
rescaling f , and by letting
is easy to see that this p
must satisfy (2.2).) An important consequence of the results we develop below is that
standard Rayleigh-Ritz discretization and multigrid solution methods can be applied with
optimal estimates that are uniform in h, -, and -. For example, we obtain optimal uniform
approximation of the gradients of displacements in the H 1 product norm. This in turn
implies analogous H 1 estimates for the stresses, which are easily obtained from the "velocity
fluxes". For related results with a different methodology and weaker norm estimates, see
[16].
Let curl j r\Theta denote the curl operator. (Here and henceforth, we use notation for
the case consider the special case in the natural way by identifying
with the if u is two dimensional, then the curl of u means the
scalar function
r\Theta
1 Perhaps a more physical choice for this artificial pressure would have been
r \Delta u, since it then
becomes the hydrostatic pressure in the incompressible limit. We chose our particular scaling because it
most easily conforms to (2.1). In any case, our results apply to virtually any nonnegative scaling of p, with
no effect on the equivalence constants (provided the norms are correspondingly scaled); see Theorems 3.1
and 3.2.
where u 1 and u 2 are the components of u.) The following identity is immediate:
r\Theta (r\Theta
interpreted as
r ? is the formal adjoint of r\Theta defined by
We will be introducing a new independent variable defined as the n 2 -vector function
of gradients of the u i , It will be convenient to view the original n-vector
functions as column vectors and the new n 2 -vector functions as either block column vectors
or matrices. Thus, given
and denoting u then an operator G defined on scalar functions (e.g.,
r) is extended to n-vectors componentwise:
and
If U i j Gu i is a n-vector function, then we write the matrix
U
U n1 U
We then define the trace operator tr according to
tr
If D is an operator on n-vector functions (e.g., its extension to matrices is
defined by
When each DU i is a scalar function (e.g., then we will want to view the extension
as a mapping to column vectors, so we will use the convention
We also extend the tangential operator n\Theta componentwise:
Finally, inner products and norms on the matrix functions are defined in the natural componentwise
way, e.g.,
Introducing the velocity flux variable
then the Stokes system (2.1) and (2.3) may be recast as the following equivalent first-order
system:
r
Note that the definition of U, the "continuity" condition r \Delta in \Omega\Gamma and the
Dirichlet condition
@\Omega imply the respective properties
r\Theta
Then an equivalent extended system for (2.6) is
r
r\Theta
Let
D(\Omega\Gamma be the linear space of infinitely differentiable functions with compact support
on\Omega and let D 0
(\Omega\Gamma denote the dual space of
D(\Omega\Gamma0 The duality pairing between D
D(\Omega\Gamma is denoted by ! \Delta; \Delta ?. We use the standard notation and definition for the Sobolev
spaces H s
the standard associated inner products are denoted
by (\Delta; \Delta)
s;\Omega and (\Delta; \Delta) s;
@\Omega , and their respective norms by k \Delta k
s;\Omega and k \Delta k s;
@\Omega . (We suppress
the superscript n because dependence of the vector norms on dimension will be clear by
context. We also
omit\Omega from the inner product and norm designation when there is no
risk of confusion.) For
coincides with L
. In this case, the norm and
inner product will be denoted by k \Delta k and (\Delta; \Delta), respectively. As usual, H s
0(\Omega\Gamma is the closure
of
with respect to the norm k \Delta k s and H
\Gammas(\Omega\Gamma is its dual with norm defined by
Define the product spaces H s
\Gammas(\Omega\Gamma with standard
product norms. Let
and
which are Hilbert spaces under the respective norms
and
@\Omega
where denote the respective
unit vectors normal and tangent to the boundary. Finally, define
Z
\Omega
It is well-known that the (weak form of the) boundary value problem (2.1)-(2.2) has a
unique solution (u;
n \Theta L 2
(e.g., see
[21, 22, 17]). Moreover, if the boundary of the
domain\Omega is C or a convex polyhedron,
then the following H 2 -regularity result holds:
(We use C with or without subscripts in this paper to denote a generic positive constant,
possibly different at different occurrences, that is independent of the Reynolds number and
other parameters introduced in this paper, but may depend on the
domain\Omega or the constant
.) Bound (2.9) is established for the case the case for general
and the case ffi ? 0 follows from the well-known linear
elasticity bound kuk 2 is the (unscaled) source term in (3.19) and
oe is the stress tensor. We will need (2.9) to establish full H 1 product ellipticity of one of
our reformulations of (2.1)-(2.2); see Theorem 3.2.
The following lemma is an immediate consequence of a general functional analysis result
due to Ne-cas [24] (see also [17]).
Lemma 2.1. For any p in L 2
, we have
Proof. See [24] for a general proof.
A curl result analogous to Green's theorem for divergence follows from [17] (Theorem
2.11 in Chapter I):
(r\Theta z;
Z
@\Omega
ds
for z 2 H(curl ; \Omega\Gamma and OE 2 H
Finally, we summarize results from [17] that we will need for G 2 in the next section. The
first inequality follows from Theorems 3.7-3.9 in [17], while the second inequality follows
from Lemmas 3.4 and 3.6 in [17].
Theorem 2.1. Assume that the
domain\Omega is a bounded convex polyhedron or has C
boundary. Then for any vector function v in either H 0 (div;
If, in addition, the domain is simply connected, then
kr
3. First-Order System Least Squares. In this section, we consider least-squares
functionals based on system (2.6) and its extension (2.8). Our primary objective here is to
establish ellipticity of these least-squares functionals in the appropriate Sobolev spaces.
Our first least-squares functional is defined in terms of appropriate weights and norms
of the residuals for system (2.6):
Note the use of the H \Gamma1 norm in the first term here. Our second functional is defined as a
weighted sum of the L 2 norms of the residuals for system (2.8):
Let
n \Theta L 2
n \Theta (H
where
@\Omega g:
Note that V 2 ae V 1 . For 2, the first-order system least-squares variational problem
for the Stokes equations is to minimize the quadratic functional G i (U; u;
find (U; u; p) 2 V i such that
Theorem 3.1. There exists a constant C independent of - such that for any (U; u; p) 2
and
Proof. Upper bound (3.5) is straightforward from the triangle and Cauchy-Schwarz
inequalities. We proceed to show the validity of (3.4) for (U; u; p)
g. Then (3.4) would follow for (U; u; p) 2 V 1 by continuity. For
any (U; u; p)
n , we have
(r
Hence, by Lemma 2.1, we have
From (3.6) and the Poincar'e-Friedrichs inequality on u we have
Using the "-inequality, 2ab - 1
for the first two products yields
Again from (3.6) and the Poincar'e-Friedrichs inequality on u we have
Using the "-inequality on the first three products and (3.7), we then have
Again using the "-inequality, we find that
Using (3.8) in (3.6) and (3.7), we now have that
The theorem now follows from these bounds, (3.8), and the Poincar'e-Friedrichs inequality
on u.
The next two lemmas will be useful in the proof of Theorem 3.2.
Lemma 3.1. (Poincar'e-Friedrichs-type inequality) Suppose that the assumptions of
Theorem 2.1 hold. Let p 2 H
where C depends only
on\Omega . Further, let q 2 (H 1
kr
where C depends only
on\Omega .
Proof. Equation
R
at some point in \Omega\Gamma The first result now
follows from the standard Poincar'e-Friedrichs inequality. The second result follows from
the fact that
R
Lemma 3.2. Under the assumptions of Theorem 2.1 with simply
connected\Omega , for any
(\Omega\Gamma we have:
2(\Omega\Gamma and each OE
is such that \DeltaOE i 2 L
@\Omega , then
jr
jr
2(\Omega\Gamma and each
@\Omega ,
then
jr
jr
Proof. The assumptions of Theorem 2.1 are sufficient to guarantee H 2 -regularity
of the Laplace equation on \Omega\Gamma that is, the second inequality in the equation
Note that tr (r ? OE OE. Then, from the above and the triangle inequality, we
have
jr
jr
which is (3.11).
applied to each column of r\Theta\Phi
imply that
since each OE i is divergence free. Eqn. (3.12) now follows from the triangle inequality as for
the case 2.
Theorem 3.2. Assume that the
domain\Omega is a bounded convex polyhedron or has C
boundary and that regularity bound (2.9) holds. Then, there exists a constant C independent
of - such that for any (U; u;
and
Proof. Upper bound (3.14) is straightforward from the triangle and Cauchy-Schwarz
inequalities. To prove (3.13), note that the H \Gamma1 norm of a function is always bounded by
its Hence, by Theorem 3.1, we have
From Theorem 2.1 and (3.9), we have
It thus suffices to show that
We will prove (3.17) only for the case because the proof for First
we assume that the
domain\Omega is simply connected. Since n \Theta U = 0 on @ the following
decomposition is admitted (see Theorems 3.4 in [17]):
and \Phi is columnwise divergence free with n \Theta (r\Theta
@
Here we choose q to satisfy
By taking the curl of both sides of this decomposition, it is easy to see that
(3.
so that k\Delta \Phik is bounded and Lemma 3.2 applies. Hence,
kr\Theta Uk 2
(by equation (3.18))
(by Lemma 3.2)
assumption (2.9) with
(by equation (3.18))
This proves (3.17) and, hence, the theorem for simply
connected\Omega\Gamma The proof for
general\Omega\Gamma that is, when we assume only that
@\Omega is C 1;1 , now follows by
an argument similar to the proof of Theorem 3.7 in [17].
We now show that the last two terms in the definition of G 2 are necessary for the
bound (3.13) to hold, even with the extra boundary condition n \Theta U = 0. We consider the
Stokes equations, so that first that we omit the term kr\ThetaUk 2 but include
the term krtr Uk 2 . We offer a two-dimensional counterexample; a three-dimensional
counterexample can be constructed in a similar manner. Let
Choose any ! 2
D(\Omega\Gamma such that \Deltar! 6= 0 and define
Clearly, n \Theta U = 0. It is easy to show that
r
However,
(r\Theta
by construction. Thus,
which cannot bound kUk 2
1 . That is, since
D(\Omega\Gamma is arbitrary, we may choose it so
oscillatory that kUk 1 =kUk is as large as we like. This prevents the bound (3.13) from
holding.
Next suppose we include the kr\Theta Uk 2 term but omit the krtr Uk 2 term. Now set
choose q i to satisfy
2. Then
We also know that
where C is independent of k. Now set
2. Then n \Theta U
where C is independent of k. On the other hand, we have
which again prevents the bound (3.13) from holding.
4. Concluding Remarks. Full regularity assumption (2.9) is needed in Theorem 3.2
only to obtain full H 1 product ellipticity of augmented functional G 2 in (3.2). This somewhat
restrictive assumption is not necessary for functional G 1 in (3.1), which supports
an efficient practical algorithm (the H \Gamma1 norm in (3.1) can be replaced by a discrete inverse
norm or a simpler mesh weighted norm; see [5] and [8] for analogous inverse norm
algorithms) and which has the weaker norm equivalence assured by Theorem 3.1.
Nevertheless, the principal result of this paper is Theorem 3.2, which establishes full H 1
product ellipticity of least-squares functional G 2 for the generalized Stokes system. Since we
have assumed full H 2 -regularity of the original Stokes (linear elasticity) equations, we may
then use this result to establish optimal finite element approximation estimates and optimal
multiplicative and additive multigrid convergence rates. This can be done in precisely the
same way that these results were established for general second-order elliptic equations (see
[11], Sections 3-5). We therefore omit this development here. However, it is important to
recognize that the ellipticity property is independent of the Reynolds parameter - (Lam'e
constants - and -). This automatically implies that the optimal finite element discretization
error estimates and multigrid convergence factor bounds are uniform in - and -). At
first glance, it might appear that the scaling of some of the H 1 product norm components
might create a scale dependence of our discretization and algebraic convergence estimates.
However, the results in [11] are based only on assumptions posed in an unscaled H 1 product
norm, in which the individual variables are completely decoupled; and since the constant -
appears only as a simple factor in individual terms of the scaled H 1 norm, these assumptions
are equally valid in this case. On the other hand, for problems where the necessary H 1
scaling is not (essentially) constant, extension of the theory of Section 3-5 of [11] is not
straightforward. Such is the case for convection-diffusion equations, which will be treated
in a forthcoming paper.



--R

Estimates near the boundary for solutions of elliptic partial differential equations satisfying general boundary conditions II


Analysis of least-squares finite element methods for the Navier-Stokes equations
Accuracy of least-squares methods for the Navier-Stokes equations
Analysis of least-squares finite element methods for the Stokes equations

A least-squares approach based on a discrete minus one inner product for first order system
On the existence


Schwarz alternating procedure for elliptic problems discretized by least-squares mixed finite elements
A mixed finite element method for the Stokes problem: an acceleration-pressure formulation Appl
An error estimate of the least squares finite element methods for the Stokes problem in three dimensions Math.
The Finite Element Method for Elliptic Problems
analysis of some Galerkin least squares methods for the linear elasticity equationns
Finite Element Methods for
Elliptic Problems in Nonsmooth Domains
Finite Element Methods for Viscous Incompressible Flows
Theoretical study of the incompressible Navier-Stokes equations by the least-squares method
A regularity result for the Stokes problem in a convex polygon
The Mathematical Theory of Viscous Incompressible Flow
Variational multigrid theory


--TR

--CTR
Hongxing Rui , Seokchan Kim , Sang Dong Kim, A remark on least-squares mixed element methods for reaction-diffusion problems, Journal of Computational and Applied Mathematics, v.202 n.2, p.230-236, May, 2007
Suh-Yuh Yang, Analysis of piecewise linear approximations to the generalized stokes problem in the velocity-stress-pressure formulation, Journal of Computational and Applied Mathematics, v.147 n.1, p.53-73, 1 October 2002
M. M. J. Proot , M. I. Gerrtisma, Least-squares spectral elements applied to the Stokes problem, Journal of Computational Physics, v.181 n.2, p.454-477, 20 September 2002
On the velocity-vorticity-pressure least-squares finite element method for the stationary incompressible Oseen problem, Journal of Computational and Applied Mathematics, v.182 n.1, p.211-232, 1 October 2005
Yves Tourigny, The Optimisation of the Mesh in First-Order Systems Least-Squares Methods, Journal of Scientific Computing, v.24 n.2, p.219-245, August    2005
J. P. Pontaza , J. N. Reddy, Space-time coupled spectral/hp least-squares finite element formulation for the incompressible Navier-Stokes equations, Journal of Computational Physics, v.197 n.2, p.418-459, 1 July 2004
J. P. Pontaza , J. N. Reddy, Spectral/hp least-squares finite element formulation for the Navier-Stokes equations, Journal of Computational Physics, v.190 n.2, p.523-549, 20 September
Youngmi Choi , Hyung-Chun Lee , Byeong-Chun Shin, A least-squares/penalty method for distributed optimal control problems for Stokes equations, Computers & Mathematics with Applications, v.53 n.11, p.1672-1685, June, 2007
P. Bolton , R. W. Thatcher, On mass conservation in least-squares methods, Journal of Computational Physics, v.203 n.1,
Sang Dong Kim , Byeong Chun Shin, H-1 least-squares method for the velocity-pressure-stress formulation of Stokes equations, Applied Numerical Mathematics, v.40 n.4, p.451-465, March 2002
J. J. Heys , T. A. Manteuffel , S. F. McCormick , J. W. Ruge, First-order system least squares (FOSLS) for coupled fluid-elastic problems, Journal of Computational Physics, v.195
J. A. Sethian , Jon Wilkening, A numerical model of stress driven grain boundary diffusion, Journal of Computational Physics, v.193 n.1, p.275-305, January 2004
J. J. Heys , T. A. Manteuffel , S. F. McCormick , L. N. Olson, Algebraic multigrid for higher-order finite elements, Journal of Computational Physics, v.204
stokes equations;multigrid;least squares