A framework for approximation of the Stokes equations in an axisymmetric domain

We develop a framework for solving the stationary, incompressible Stokes equations in an axisymmetric domain. By means of Fourier expansion with respect to the angular variable, the three-dimensional Stokes problem is reduced to an equivalent, countable family of decoupled two-dimensional problems. By using decomposition of three-dimensional Sobolev norms we derive natural variational spaces for the two-dimensional problems, and show that the variational formulations are well-posed. We analyze the error due to Fourier truncation and conclude that, for data that are sufficiently regular, it suffices to solve a small number of two-dimensional problems.


Introduction
To determine approximate solutions to fluid flow problems in three-dimensional geometries is a computationally demanding task. In this paper we present a framework for efficiently solving the stationary, incompressible Stokes equations in an axisymmetric domainΩ, which is obtained by rotating its half section Ω around the symmetry axis.
We use Fourier expansions with respect to the angular variable θ, both of the solution and the data, to reduce the three-dimensional Stokes problem to an equivalent, countable family of decoupled two-dimensional problems (set in Ω) for the Fourier coefficients. A natural way to approximate the three-dimensional problem is then to use Fourier truncation and, to obtain a fully discrete scheme, compute approximate solutions to a finite number of the two-dimensional problems.
This is an established technique to approximate boundary value problems that are invariant by rotation. Early error analysis results for second order elliptic problems can be found in [8], and for Poisson's equation in domains with reentrant edges in [6]; both using finite element approximation for the two-dimensional problems. We refer to [5] for additional references to problems described by Laplace or wave equations, the Lamé system, Stokes or Navier-Stokes systems, and Maxwell's equations, and to [6] for early references to algorithms and applications.
We analyze the error due to Fourier truncation and show that, for data that are sufficiently regular with respect to θ, it suffices to solve a small number of two-dimensional problems, which makes the method efficient. Also, the decoupling of the two-dimensional problems makes it suitable for parallel implementation. A further advantage is simplification of the mesh-generation, which is only required for the two-dimensional half section Ω.
An added complexity is that the natural, variational spaces for the Fourier coefficients turn out to be weighted Sobolev spaces, where the weight is either the distance to the symmetry axis, or its inverse. We derive these spaces by decomposing (through a change of variables to cylindrical coordinates) the three-dimensional norms for the relevant spaces L 2 (Ω), L 2 0 (Ω), H 1 (Ω) 3 , H 1 0 (Ω) 3 , and H −1 (Ω) 3 , into sums over all wavenumbers. As a result, we show that the three-dimensional spaces are isometrically isomorphic to a certain subspace of the Cartesian product, over all wavenumbers, of the two-dimensional weighted spaces. The characterizations of the three-dimensional spaces thus obtained are in agreement with the results in [3], where characterizations of H s (Ω) and H s (Ω) 3 by Fourier coefficients for any positive real s are derived; here non-integer order spaces are treated first and the derivation for integer order spaces is then based on Hilbert space interpolation.
As recently shown in [5], the more direct approach for Sobolev spaces of integer order (based on changing to cylindrical coordinates in the three-dimensional norms) results in equivalent norms (compared with the norms in [3]), but where the equivalence constants (unlike in [3]) are independent of the domain. In [5], characterizations of H m (Ω) by Fourier coefficients for any positive integer m are facilitated by the introduction of new differential operators ∂ ζ = 1 √ 2 (∂ x − i∂ y ) and ∂ζ = 1 √ 2 (∂ x + i∂ y ). The results for vector spaces are then derived from a relation between scalar and vector norms for the Fourier coefficients. In this paper, we work with the differential operators ∂ x and ∂ y and, in the vector case, derive the characterization by directly rewriting the H 1 (Ω) 3 -norm. We compare our results with [5] in Appendix A. We also refer to the early results in [8], where this direct approach was used to characterize the scalar spaces H m (Ω), for m = 1, 2, 3. The purpose of the present work is to give a comprehensive presentation directly aimed at the Stokes problem providing, inter alia, detailed derivations of the relevant two-dimensional spaces and norms. Taking as starting-point, in fact, Fourier decompositions of the three-dimensional inner products (·, ·) L 2 (Ω) , (·, ·) (H 1 (Ω)) 3 , and (·, ·) (H 1 0 (Ω)) 3 , additionally enables us to derive a decomposition of the negative norm · (H −1 (Ω)) 3 , and to highlight the relation between the three-dimensional weak formulation of the Stokes problem and the two-dimensional weak formulations for the Fourier coefficients.
Examples of how to build on this framework by discretizing the two-dimensional problems can be found in [3], where spectral methods are used, and [2], where two families of finite elements of order 2 (one with continuous pressure corresponding to the Taylor-Hood element and one with discontinuous pressure) are used. The case with an axisymmetric solution (where only the Fourier coefficient of order 0 is considered, and the angular velocity component is equal to zero), has been treated with finite elements in [1,9,10].
A paper in preparation will be devoted to design and analysis of stabilized finite elements for the two-dimensional problems.
An outline of the paper is as follows: • In Sec. 2, we give some examples of axisymmetric domains and state the stationary, incompressible Stokes equations. • In Sec. 3, we recall some basic formulas, and state the Stokes problem in cylindrical coordinates. • In Sec. 4, we use Fourier expansion with respect to the angular variable to reduce the three-dimensional Stokes problem to a countable family of two-dimensional problems. • In Sec. 5, we derive natural variational spaces for the Fourier coefficients by decomposing the relevant three-dimensional norms into sums over all wavenumbers.
• In Sec. 6, we state variational formulations of the two-dimensional problems and show that these are well-posed. • In Sec. 7, we introduce two families of anisotropic spaces that we need to analyze the error due to Fourier truncation. • In Sec. 8, we prove an error estimate due to Fourier truncation.

Model description
We consider fluid flow in a bounded domainΩ which is invariant by rotation around an axis. We begin by discussing, and give a few examples of, such domains, following the notation in [3]. We then state the stationary, incompressible Stokes equations, which we use to model the flow.
2.1. Axisymmetric domains. An example of an axisymmetric domainΩ is given in Figure 1. The axisymmetric domain is obtained by rotating its half section (meridian domain) Ω around the symmetry axis. We assume that Ω is polygonal. in Figure 2 and Figure 3. Note that the boundary ∂Ω = Γ ∪ Γ 0 of the half section Ω consists of two parts where Γ 0 , the interior of the part of ∂Ω contained in the symmetry axis, is a kind of artificial boundary.
Remark 2.1. In [3], Γ 0 is assumed to be the union of a finite number of segments with positive measure, which means that Ω is not allowed to meet the symmetry axis at isolated points. This assumption, as noted in [5], implies thatΩ (as well as its polygonal half section) is a Lipschitz domain, and it guarantees existence of certain trace operators, needed since [3] uses vanishing traces on Γ 0 in the definition of the Fourier coefficient spaces. The more direct approach for integer order Sobolev spaces used in this paper, and more generally in [5], allows for more general axisymmetric domains whose intersection with the symmetry axis is not necessarily a union of intervals and where the trace operators are not well defined.

2.2.
Stokes problem. We model fluid flow through an axisymmetric domainΩ by the stationary, incompressible Stokes equations where the unknowns are the velocityȗ and the pressurep u =ȗ x e x +ȗ y e y +ȗ z e z ∈ H 1 (Ω) 3 , g =g x e x +g y e y +g z e z ∈ H 1 2 (∂Ω) 3 .
For a vector fieldv =v x e x +v y e y +v z e z defined onΩ, we will letv denote both the vector field itself and its Cartesian component vector (v x ,v y ,v z ) T . Note that the divergence-free property divȗ = 0 implies a necessary compatibility (zero flux) condition ong: wheren denotes the unit outward normal vector toΩ on ∂Ω, and dȂ is (the magnitude of) the area element on ∂Ω.
We recall the standard definitions of the Lebesgue and Sobolev spaces (with all derivatives being taken in the sense of distributions) and the corresponding inner products (·, ·) L 2 (Ω) and (·, ·) H 1 (Ω) .

Remark 2.2.
We consider spaces of complex-valued functions. In Section 4 we will use Fourier expansions of the data and the unknowns, to reduce the threedimensional Stokes problem inΩ to a countable family of two-dimensional problems in the half section Ω. Since the Fourier coefficients (also for real-valued functions) are complex-valued, the two-dimensional variational spaces, that we will define in Section 5, will be spaces of complex-valued functions. Since the interplay between the three-and two-dimensional spaces, which will play an important part in what follows, is developed most naturally when all spaces contain complex-valued functions, we will use that also for the three-dimensional spaces, even though in practice the three-dimensional data and solution will be real-valued.
We also recall the subspaces Remark 2.3. By the dual space H * we mean the linear space of continuous antilinear functionals on H, so that, in particular, On H 1 0 (Ω) we will use the semi-norm which is a norm, equivalent to · H 1 (Ω) , on H 1 0 (Ω), and the corresponding inner product (·, ·) H 1 0 (Ω) . Denoting by γ 0 the linear and continuous trace operator defined on H 1 (Ω), we have

Cylindrical coordinates
Since the domainΩ is axisymmetric, we change from Cartesian coordinates We defineΩ as the product of the half section Ω and (−π, π] and, by analogy, whereΩ andΓ are point sets in cylindrical coordinates corresponding to the domain Ω and its boundary ∂Ω, respectively. 3.1. Basic formulas. We recall the identities (see Figure 4 for an illustration in the xy-plane) x y e x e y r θ e r e θv =v x e x +v y e y = v r e r + v θ e θ P Figure 4. A point P with Cartesian coordinates (x, y) and cylindrical coordinates (r, θ), where x = r cos θ and y = r sin θ. A vectorv, with tail in P , with Cartesian components (3.4) e x = cos θ e r − sin θ e θ , e y = sin θ e r + cos θ e θ , relating the partial derivatives and the orthonormal basis vectors of the Cartesian and cylindrical coordinate systems, and T of a vector fieldv =v x e x +v y e y +v z e z = v r e r + v θ e θ + v z e z expressed in the two coordinate systems. We will write (3.5) From these identities follow the formulas for the gradient and the Laplacian operator acting on a scalar functionv(x, y, z) = v(r, θ, z) and for the divergence and the vector Laplacian operator acting on a vector function v(x, y, z) = R θ v(r, θ, z) where the last two terms in the radial and angular components of ∆v result from the θ-dependence of e r = cos θ e x + sin θ e y and e θ = − sin θ e x + cos θ e y , by the relations ∂ θ e r = e θ and ∂ θ e θ = −e r .

Stokes problem in cylindrical coordinates. Expressing both the data and the unknowns in cylindrical coordinates
u =ȗ x e x +ȗ y e y +ȗ z e z = u r e r + u θ e θ + u z e z , g =g x e x +g y e y +g z e z = g r e r + g θ e θ + g z e z , and g = (g r , g θ , g z ) T onΓ, from (3.6)-(3.9) we can write the Stokes problem (2.1) in cylindrical coordinates where, from (3.7),

Fourier expansion
A natural way to reduce the three-dimensional Stokes problem (3.10) inΩ to a countable family of two-dimensional problems in Ω, is to use Fourier expansion with respect to the angular variable θ, both of the solution where u k = u k r , u k θ , u k z T , and of the data 4.1. Two-dimensional problems. Inserting the Fourier expansions (4.1)-(4.4) into (3.10) results, since the Stokes problem is linear and invariant by rotation (which means that the coefficients of the Stokes operator in cylindrical coordinates do not depend on θ), in uncoupled two-dimensional problems for each Fourier coefficient pair u k , p k , k ∈ Z: in Ω, where ∆ a denotes the axisymmetric part of ∆: Remark 4.1. We will, for all k ∈ Z, show existence and uniqueness of solutions to (4.5) in Section 6.2.
Remark 4.2. By taking the complex conjugate of (4.5), it is easy to see that for real-valued data f and g, in which case (lettingf k , with some ambiguity of notation, denote the complex conjugate of f k ) the pair ū k ,p k solves (4.5) with k replaced by −k. This means that the Fourier coefficients of the solution will also satisfy (corresponding, of course, to a unique, real-valued solution of the three-dimensional Stokes problem for real-valued data) so in the practical case with real-valued data we only need to solve the problems (4.5) for k ≥ 0 and, in addition, the solution for k = 0 will be real-valued. where n = (n r , n z ) T denotes the unit outward normal vector to Ω on Γ, and ds is (the length of ) the line element along Γ.

Variational spaces
In Section 6.1, we will state variational formulations of the two-dimensional problems (4.5). To determine natural variational spaces for the Fourier coefficients defined on the half section Ω, we start by expressing the L 2 (Ω)-and H 1 (Ω) 3 -inner products in cylindrical coordinates, as integrals overΩ, and use Fourier expansions to derive decompositions of the inner products, and associated norms, into sums over all wavenumbers. Based on the structure of the different terms, which are weighted integrals over Ω of the Fourier coefficients and their derivatives, we define weighted Sobolev spaces on Ω.
As a result, we obtain characterizations of the three-dimensional spaces L 2 (Ω) and H 1 (Ω) 3 , in terms of two-dimensional weighted spaces on Ω for all Fourier coefficients. In particular, we show that the three-dimensional spaces are isometrically isomorphic to a subspace of the Cartesian product, over all wavenumbers, of the two-dimensional weighted spaces. As a corollary, we obtain corresponding characterizations of the subspaces L 2 0 (Ω) and H 1 0 (Ω) 3 .
Using the results for H 1 0 (Ω) 3 , we finally derive a characterization (also in terms of spaces for all Fourier coefficients) of its dual H −1 (Ω) 3 .

5.1.
Fourier decomposition of inner products and norms. In cylindrical coordinates, the L 2 (Ω)-inner product of two scalar functions p(x, y, z) = p(r, θ, z),q(x, y, z) = q(r, θ, z), (5.1) is expressed as an integral overΩ, defined by (3.1), with weight r: For two vector functions can, through repeated use of relations (3.3), (3.5) and the Pythagorean trigonometric identity, be expressed in cylindrical coordinates: We now consider Fourier expansions of the family e ikθ +∞ k=−∞ of basis functions, we obtain (lettingq k , as before, denote the complex conjugate of q k ): and, similiarly, in four other representative cases: Inserting (5.7)-(5.8) (and analogous results for the remaining θ-integrals) into (5.2), (5.4) gives (5.9) (p,q) L 2 (Ω) = k∈Z Ω p kqk r dr dz, expressing the L 2 (Ω)-inner product and the H 1 (Ω) 3 -inner product as sums, over all wavenumbers, of weighted integrals over the half section Ω of the Fourier coefficients and their derivatives. The corresponding decompositions of the associated norms are:

5.2.
Weighted Sobolev spaces on Ω. Led by (5.11)-(5.12), where each term is an integral over Ω with weight r or r −1 , we first introduce the spaces Next, we define H 1 1 (Ω) as the space of functions in L 2 1 (Ω) such that their partial derivatives (being taken in the sense of distributions) of order 1 belong to L 2 1 (Ω), equipped with the semi-norm Remark 5.1. The definition can be extended in a natural way to H m 1 (Ω) for an arbitrary integer m ≥ 2 and further, by interpolation, to H s 1 (Ω) for non-integer s > 0.
We will also need the weighted space It can be proved [8,Proposition 4.1] that all functions in V 1 1 (Ω) have a null trace on the part Γ 0 of the boundary contained in the z-axis.
We finally introduce the subspaces consisting of functions in H 1 1 (Ω) and V 1 1 (Ω) that vanish on the part Γ = ∂Ω\Γ 0 of the boundary that is not contained in the z-axis.
All spaces defined above are Hilbert spaces for the inner products associated with the given norms.

Characterization of H 1 (Ω)
3 and H 1 0 (Ω) 3 . Similarly, led by the righthand side terms in the identity (5.12), for all k ∈ Z, we define spaces H 1 (k) (Ω) for the Fourier coefficient triples by the norms To see that (5.17) satisfies the properties of a norm, we consider the function v k ∈ (H 1 (Ω)) 3 corresponding to a single Fourier coefficient and, from (5.12), note that v k From (5.18) and the norm properties of H 1 (Ω) 3 , we immediately get: Also from (5.18), the completeness of H 1 (k) (Ω) is a consequence of the completeness of (H 1 (Ω)) 3 .

Two-dimensional problems
In this section we state variational formulations of the two-dimensional problems for the Fourier coefficients, which are set in the spaces from the previous section. We show existence, uniqueness and stability from inf-sup conditions. We conclude by discussing some features of the special case with axisymmetric data.
6.1. Variational formulation. Introducing, for all k ∈ Z, the sesquilinear forms where div k is defined by (4.6), we consider the following variational formulation of (4.5): Remark 6.1. Note that the term B k is conjugated in the first equation in (6.1), making the sesquilinear form Hermitian.
Remark 6.2. We use the same notation g k for the Fourier coefficients of the Dirichlet boundary datag, and their liftings in H 1 (k) (Ω). Existence of these liftings follows from the three-dimensional trace theorem sinceg ∈ H 1 2 (∂Ω) 3 and thus, again using the same notation, admits a liftingg ∈ H 1 (Ω) 3 .
We have thus proved (see, e.g., [4, Theorem 4.2.3]) the following Theorem 6.5. Forf ∈ H −1 (Ω) 3 andg ∈ H 1 2 (∂Ω) 3 , where g 0 satisfies the compatibility condition (4.7), the variational formulations (6.1) have unique solutions for all k ∈ Z. The solutions are bounded, uniformly in k, by: (k) (Ω) ). (6.5) 6.3. Axisymmetric data. For axisymmetric data, i.e., only f 0 and g 0 are nonzero, it follows immediately from uniqueness that only u 0 and p 0 are non-zero, which means that the solution will also be axisymmetric. Also note that problem (6.1) for k = 0 decouples into two independent problems: one for (u 0 r , u 0 z , p 0 ) and one for u 0 θ . If the data are real-valued, these problems have real-valued solutions.

Anisotropic spaces
We introduce two families H ±1,s (Ω) of anisotropic spaces, where s is a nonnegative real number measuring "extra regularity" in the angular direction, that we (following [3]) will use to derive an estimate of the error due to Fourier truncation in Section 8.
In conclusion, the H 1 (k) (Ω)-norms for the Fourier coefficient spaces that we have derived by directly rewriting the H 1 (Ω) 3 -norm, are the same as those obtained from the relation (A.1). In addition to conveying structural understanding, relation (A.1) also facilitates the general study, by induction, for positive integer order Sobolev spaces carried out in [5]. For our exposition, directly aimed at the Stokes problem, the method we have chosen is straightforward to follow. The same remark can be made as to the difference between introducing the differential operators ∂ ζ = 1 √ 2 (∂ x − i∂ y ) and ∂ζ = 1 √ 2 (∂ x + i∂ y ) in [5], and our choice of working directly with ∂ x and ∂ y .