Decision models enhancing environmental flow sustainability: A strategic approach to water resource management

1.   Introduction and main result

We consider the Stokes system in which the viscosity is a periodically varying function of the space variable with small period $\epsilon > 0$. Many physical phenomena (boiling flows, porous media, oil reservoirs, etc.) lead to mixture of fluids with different viscosities. For incompressible slow or creeping flows, such a situation is modeled by the system (1) for a Stokesian fluid with variable viscosity which is further assumed to be a periodic function. From the point of view of application, it is difficult to realize such a periodic distribution of droplets of one fluid in another without deforming the periodic structure, and (1) may seem as too much of an idealized system. Therefore, we also treat another model, which is a variant of the Stokes system and is physically more relevant. Namely, we consider the so-called incompressible elasticity system (10) which corresponds to a mixture of incompressible elastic phases in a composite material (this situation is quite common for rubber or elastomers).

The first goal of this paper is to study homogenization of the above systems via Bloch Wave Method which is based on the fact that the homogenized operator can be defined using differential properties of the bottom of the so-called Bloch spectrum. The second goal of the paper is to explore this regularity issue which is delicate for the systems under consideration because of the presence of the incompressibility condition. These points are elaborated below.

Through out the paper, we will follow the usual Einstein summation convention with respect to repeated indices. We introduce now our first model. Assuming that the viscosity (denoted by $\mu$ in the sequel) is a periodic function, the goal is to capture the effective viscosity of the mixture. To write down the model we start with a $1$-periodic function $\mu= \mu (y) \in L^\infty (\mathbb{T}^d)$ or equivalently, a $Y$-periodic function where $Y= ]0, 1[^d$ which represents the viscosity of the fluids. Here $\mathbb{T}^d$ is the unit torus in $\mathbb{R}^d.$ We assume $\mu (y) \geq \mu_0 >0$ a.e in $\mathbb{T}^d$. Denote by $\mu^\epsilon= \mu^\epsilon(x) = \mu \left(\frac{x}{\epsilon} \right)$ the corresponding scaled function which is $\epsilon$-periodic. With $f = f (x) \in L^2(\Omega)^d$ representing external force, we consider the Stokes system with no-slip boundary condition in a bounded connected nonempty open set $\Omega \subset \mathbb{R}^d$ having Lipschitz boundary, :

As usual, $u^\epsilon$ and $p^\epsilon$ represent respectively the velocity and pressure fields of the fluid. Well-posedness theory of (1) is classical [15]. We recall some of its elements. To write down the weak formulation, we introduce the spaces

Here $\nu$ denotes unit outward normal to $\partial \Omega.$ Multiplying (1) by $v \in V$ gives the following problem for $u^\epsilon$ which does not involve $p^\epsilon:$ Find $u^\epsilon \in V$ satisfying

The classical Lax-Milgram Lemma (essentially, Riesz Representation Theorem due to the symmetry of our bilinear form) ensures existence and uniqueness of a solution $u^\epsilon\in V$ for (3). To get the pressure field one applies de Rham's Theorem in the following form [15]:

which implies that the pressure $p^{\epsilon}$ in (1) belongs to $L^2(\Omega)$. Since $\Omega$ is a connected set, the pressure is defined up to an additive constant. To guarantee the uniqueness of the pressure, we seek $p$ in the space $L^2_0(\Omega)= \{f\in L^2(\Omega) : \int_\Omega f\ dx = 0\}$ with $L^2$ norm. Moreover, by using Poincaré inequality and inf-sup inequality [15], one shows that the solution $(u^{\epsilon},p^{\epsilon})\in (H^1_{0}(\Omega))^d\times L^2_0(\Omega)$ of (1) are uniformly bounded, namely there exists a constant $C$, independent of $\epsilon$, such that

We are interested here in the homogenization limit of (1), that is the asymptotic limit of the solution $(u^\epsilon, p^\epsilon)$ as $\epsilon \rightarrow 0.$ This problem is very classical and its solution by means of a combination of two-scale asymptotic expansions and the method of oscillating test functions was provided in various references, including [8,16,20]. We recall their main results and follow the notations of [8,chapter Ⅰ,section 10]. The homogenized tensor $(A^{*})_{\alpha \beta}^{kl}$, which represents " effective viscosity", is defined by

in which figure the cell test functions $\{\chi^{k}_{\alpha};\ \alpha, k =1 \ldots d\}\in H^1(\mathbb{T}^d)^d$ and $\{\Pi^{k}_{\alpha};\ \alpha, k =1 \ldots d\}\in L^2(\mathbb{T}^d)$ solutions of the following problem in the unit torus $\mathbb{T}^d$:

We impose $\int_{\mathbb{T}^d} \chi^{k}_\alpha\ dy = \int_{\mathbb{T}^d} \Pi^{k}_\alpha\ dy = 0$ to obtain uniqueness of the solutions. It is easy to see that the above homogenized tensor possesses the following ''simple" symmetry, for any indices $1\leq\alpha,\beta,k,l \leq d$,

which corresponds to the fact that the fourth-order tensor $A^{*}$ is a symmetric linear map from the set of all matrices (or second-order tensors) into itself. Since we follow the notations of [8], the simple symmetry (8) seems a bit awkward since it mixes Latin and Greek indices but it is just the usual symmetry for a pair of indices $(k,\alpha)$ and $(l,\beta)$ in a fourth-order tensor. In other words, (8) holds for a simultaneous permutation of $k,l$ and $\alpha,\beta$. It is straight-forward (see [8]) to check that the tensor $A^{*}$ is positive-definite so that the following system (9) is well-posed.

Theorem 1.1. The homogenized limit of the problem (1) is

More precisely, we have the convergence of solutions:

Note that the simple symmetry (8) does not imply that $A^{*}$ is symmetric in $k,l$ or in $\alpha,\beta$. However, in the homogenized equation (9), since $A^{*}$ is constant, only its symmetric version, obtained by symmetrizing in both $k,l$ and $\alpha,\beta$, plays a role.

Let us next consider the second model of incompressible elasticity :

Here the strain rate tensor is given by

As before, there exists a unique solution $(u_s^{\epsilon},p_s^{\epsilon})$ of the above problem (10) in $(H^1_{0}(\Omega))^d\times L^2_0(\Omega)$ and using Korn's inequality and the inf-sup inequality, the following uniform bound can be proved :

where the constant $C$ does not depend on $\epsilon$. Here the homogenized tensor ${(A^{*}_s)}_{\alpha \beta} ^{kl}$ is given by

where the cell test functions $\widetilde{\chi}^{k}_{\alpha}\in H^1(\mathbb{T}^d),\widetilde{\Pi}^{k}_\alpha\in L^2(\mathbb{T}^d)$ are now solutions in the torus $\mathbb{T}^d$ of

We impose $\int_{\mathbb{T}^d} \widetilde{\chi}^{k}_\alpha = \int_{\mathbb{T}^d} \widetilde{\Pi}^{k}_\alpha =0$. It is known [10] that the above homogenized tensor possesses the following "full" symmetry, for any indices $1\leq\alpha,\beta,k,l \leq d$,

which corresponds to the fact that the fourth-order tensor $A^{*}_s$ is a symmetric linear map from the set of all symmetric matrices into itself (the conditions (14) are the usual symmetry conditions for Hooke's laws in linearized elasticity). The homogenization limit of the problem (10) is again of the form (9) with $A^{*}_s$ replacing $A^{*}$.

The first goal of this paper is to give an alternate proof of Theorem 1.1 using the Bloch Wave Method instead of two-scale asymptotic expansions and the method of oscillating test functions. The notion of Bloch waves is well-known in physics and mathematics [8,11,19,22]. Bloch waves are eigenfunctions of a family of "shifted" spectral problems in the unit cell $Y$ for the corresponding differential operator. Its link with homogenization theory was first explored in [8,13,17]. The key point is that the homogenized operator can be defined in terms of differential properties of the bottom of the Bloch spectrum. The second goal of this paper is to explore this issue which is especially delicate in the case of Stokes equations. Indeed, it was discovered in [7] that the Bloch spectrum for the Stokes equations is not regular enough at the origin because of the incompressibility constraint. Therefore, its differential properties are all the more intricate to establish. Here we complete the task started in [7] and in particular we prove a conjecture of [7] on the homogenization of the Stokes system (1). Since the treatment of the incompressible elasticity system (10) is almost analogous to that of (1), we focus on (1) and we content ourselves in highlighting the main differences for (10) throughout the sequel.

The Bloch wave method for scalar equations and systems without differential constraints (like the incompressibility condition) was studied in [12,13,14,21]. In such cases, this approach gives a spectral representation of the homogenized tensor $A^{*}=(A^{*})^{kl}_{\alpha\beta}$ in terms of the lowest energy Bloch waves and their behaviour for small momenta (what we call the bottom of the spectrum). For instance, the homogenized matrix in the scalar case was found to be equal to one -half of the Hessian of the ground energy (or first eigenvalue) at zero momentum. For a system, several bottom eigenvalues play a role and they are merely directionally differentiable by lack of simplicity. In the present case of the Stokes system, the situation is more complicated. The main characteristic of the Stokes system is the presence of the differential constraint expressing incompressibility of the fluid. One of its effects is that the Bloch energy levels are degenerate and the corresponding eigenfunctions are discontinuous at zero momentum. Even though energy levels are continuous at zero momentum, the second order derivatives are not (cf. Theorem 3.1). Thus, we cannot really make sense of the eigenvalue Hessian at zero momentum. Further, it is not clear if the homogenized tensor can be fully recovered from the Bloch spectral data. In fact, this issue is left open in [7]. In the non-self adjoint case treated in [21], only the symmetric part of the homogenized matrix is determined by Bloch spectral data and this is enough to determine the homogenized operator uniquely. Combining all these difficulties, the homogenization of Stokes system using Bloch waves is an interesting issue which is not a direct extension of previous results. Our work, roughly speaking, shows that Bloch spectral data does not determine the homogenized tensor uniquely, but determines the homogenized operator uniquely. This is in sharp contrast with the linear elasticity system treated in [14] in which the homogenized tensor was uniquely determined from Bloch spectral data. We see thus the effect of differential constraints (the incompressibility condition in the case of Stokes equations) on the homogenization process via Bloch wave method. For further discussion on this point, see Section 4. Bloch wave method of homogenization presented in Section 5 consists of localizing (1), taking its Bloch transform and passing to the limit to get the localized version of homogenized system in the Fourier space. Passage to the limit in the Bloch method is straight forward, though arguments are long. We do not run into the classical difficulty of having a product of two weakly convergent sequences. In fact, we use the Taylor approximation of Bloch spectral elements which gives strongly convergent sequences. This is one of the known features of the method. The required homogenized system is obtained by making a passage to the physical space from the Fourier space. Extracting macro constitutive relation and macro balance equation from the localized homogenized equation in the Fourier space turns out to be not very straight forward because of differential constraints.

Let us end this discussion with two general remarks on Bloch wave method. First one is about the nature of convergence of the homogenization process. It is well-known in the homogenization theory that the convergence in Theorem 1.1 is only weak and not strong. To have strong convergence, we need the so-called correctors [8]. Within Block wave theory, correctors are discussed in [12] for the scalar equation. We do not construct explicitly correctors for Stokes system in this paper, even though all necessary ingredients are presented. Because of the lack of smoothness of the bottom level Bloch spectrum, corrector issue is worth considering separately. The second remark is about non-periodic coefficients. Bloch wave approach to homogenization is well developed only in the case of periodic coefficients. It is known that for some restricted class of locally-periodic/modulated coefficients, new phenomena (like localization) may appear [2,3,4]. We are not aware of a Bloch wave approach for more general coefficients.

The plan of this paper is as follows. In section 2, we recall from [7] the properties of Bloch waves associated with the Stokes operator. It turns out that the Bloch waves and their energies can be chosen to be directionally regular, upon modifying the spectral cell problem at zero momentum. Bloch transform using eigenfunctions lying at the bottom of the spectrum is also introduced in this section. Its asymptotic behaviour for low momenta is also described. Next, Section 3 is devoted to the computation of directional derivatives of Bloch spectral data. Even though these results are essentially borrowed from [7], some new ones are also included because of their need in the sequel. In particular we derive the so-called propagation relation linking the homogenized tensor $A^{*}$ with Bloch spectral data, and the extent to which it determines homogenized tensor is studied in Section 4. Using this information, we prove Theorem 1.1 in Section 5 following the Bloch wave homogenization method.

Note added in proof. At the end of the introduction we claimed that we were not aware of a Bloch wave approach for non-periodic coefficients. We recently learned about a new work in this direction by A. Benoit and A. Gloria, ''Long-time homogenization and asymptotic ballistic transport of classical waves", which will appear in Annales Scientifiques de l'Ecole Normale Supérieure.

2.   Bloch waves

In this section, we introduce Bloch waves associated to the Stokes operator following the lead of [7]. The Bloch waves are defined by considering the shifted (or translated) eigenvalue problem in the torus $\mathbb{T}^d$ parametrized by elements in the dual torus which we take to be $2\pi\mathbb{T}^d$. We denote by $y$ the points of the original torus and by $\eta$ the points of the dual torus. The spectral Bloch problem amounts to find $\lambda = \lambda (\eta) \in \mathbb{R}$, $\phi =\phi (\eta) \in (H^1 (\mathbb{T}^d))^d$, with $\phi \neq 0$ and $\Pi = \Pi (\eta) \in L^2 (\mathbb{T}^d)$, satisfying

The solutions of (15) $\phi,\Pi$ are a priori complex valued, so all functional spaces are complex valued too. Here, we denote

the shifted gradient operator, with $i$ the imaginary root $\sqrt{-1}$. Its action on a vector function $\phi$ yields a matrix: $(D (\eta) \phi )_{kl} = \frac{\partial \phi_l}{\partial y_k} + i \eta_k \phi_l$ for all $k, l=1, \ldots ,d$. The corresponding divergence operation yields a scalar: $D (\eta) \cdot \phi =\frac{\partial \phi_k}{\partial y_k} + i \eta_k \phi_k.$ Analogously, if $\phi$ is a matrix function then its shifted divergence $D (\eta) \cdot \phi$ is a vector function obtained by acting $D (\eta)$ on the column vectors of $\phi.$ At this stage of discussion, spectral problem (15) is stated only formally. Rigourous versions of it with modification will be formulated and used below. For reasons of self-adjointness, its eigenvalues will be real and that is why, we take $\lambda$ to be real, without any loss of generality. Though we may think of eigen-solution $\Pi$ as some sort of complex-valued pressure, there is no guarantee that it represents physical pressure field of a real fluid. Let us remark that the system (15) is only a mathematical model which appears as an useful intermediate step in the study of a physical model.

The main feature of (15) is that the state space keeps varying with $\eta$ due to the differential constraints defined by the incompressibility of the fluid. That is why, the standard spectral theory for elliptic operators does not apply as such; it has to be modified. This is accomplished in [18]. Secondly, it is easily seen that when $\eta=0$, the corresponding eigenvalue $\lambda (0)$ is equal to zero and its multiplicity is $d$. In fact, we can take $e_k, k=1 \ldots d$ as eigenvectors (with corresponding eigen-pressure being zero). Because of this degeneracy, spectral elements of (15) are not guaranteed to be smooth at $\eta=0.$

As discussed in introduction, lack of regularity of the Bloch spectrum at $\eta=0$ is an issue because the representation of the homogenized tensor in terms of Bloch spectral elements is then not clear. This was not the case for the scalar problem, the lack of regularity of the Bloch spectrum at $\eta=0$ does not appear there, in fact the Bloch spectrum is analytic near $\eta=0$, see [13].To overcome the difficulty in the present case, the idea is to consider directional regularity as we approach $\eta=0$ [14]. Accommodating the directional limit at $\eta=0$ requires a modification of the above shifted problem with the addition of a new constraint and corresponding Lagrange multiplier in the equation [7]. Fixing a direction $e \in \mathbb{R}^d, |e| =1$ and taking $\eta=\delta e,$ with $\delta>0$, we consider the modified problem: find $\lambda(\delta) \in \mathbb{R}$, $\phi(.;\delta) \in (H^1 (\mathbb{T}^d))^d$, $q(.;\delta) \in L^2_0(\mathbb{T}^d)$ where $L^2_0(\mathbb{T}^d) = \{ q \in L^2 (\mathbb{T}^d); \ \int_{\mathbb{T}^d}q = 0 \}$ and $q_0(\delta) \in \mathbb{C}$ satisfying

Note that if $\delta \neq 0$ then the relation $e \cdot \int\limits_{\mathbb{T}^d} \phi(.;\delta) =0$ can be easily obtained from $D (\delta e)\cdot\phi(.;\delta) =0$ simply by integration. (However, this is not the case if $\delta =0$.) Hence (16) is the same as (15) provided $\delta \neq 0$ and $\eta = \delta e$. However for $\delta=0,$ (15) is not good because the condition $e \cdot \int\limits_{\mathbb{T}^d} \phi(.;\delta) =0$ is not included. See [7] on the appearance of this new constraint and the corresponding Lagrange multiplier $q_0(\delta)e$.

It is natural to consider the system (16) with $\delta$ small as a perturbation from the following one which corresponds to $\delta =0$. We fix a unit vector $\hat{\eta}\in \mathbb{S}^{d-1}$ and we consider the eigenvalue problem: find $\nu (\hat{\eta})\in\mathbb{R}$, $w(. ,\hat{\eta}) \in (H^1 (\mathbb{T}^d))^d$, $q(.;\hat{\eta}) \in L^2_0(\mathbb{T}^d)$ and $q_0(\hat{\eta}) \in \mathbb{C}$ satisfying

Existence of eigenvalues and eigenvectors for either (16) or (17) is proved in [7]. Let us recall their result, by specializing to the eigenvalue $\nu(\hat{\eta}) =0$ of (17). Note that $\nu(\hat{\eta})=0$ is clearly an eigenvalue of multiplicity $(d-1)$ of (17) with corresponding eigenfunctions being constants, namely $q^0_{m,\hat{\eta}} =0$, $q^0_{0, m,\hat{\eta}} = 0$ and $\phi_{m,\hat{\eta}}^0 (y)$ is a constant unit vector of $\mathbb{R}^d$ orthogonal to $\hat{\eta}$ for $m=1, \ldots , (d-1)$, say $\{\phi^0_{1,\hat{\eta}}, \ldots \phi^0_{d-1,\hat{\eta}}\}$. Doing perturbation analysis of the above situation, the following result was proved in [7].

Theorem 2.1. Fix $\hat{\eta}\in \mathbb{S}^{d-1}$. Consider the first $(d-1)$ eigenvalues of (16). There exists $\delta_0>0$ and exactly $(d-1)$ analytic functions defined in the real interval $|\delta|\leq \delta_0$, $\delta \mapsto ( \lambda_{m,\hat{\eta}}(\delta), \ \phi_{m,\hat{\eta}}(.;\delta), \ q_{m,\hat{\eta}}(.;\delta), \ q_{0, m,\hat{\eta}}(\delta) )$, for $m= 1, \ldots , (d-1)$, with values in $\mathbb{R} \times (H^1 (\mathbb{T}^d))^d \times L^2_0 (\mathbb{T}^d) \times \mathbb{C}$, such that

(ⅰ) $\displaystyle{ \lambda_{m,\hat{\eta}}(\delta) \bigg|_{\delta =0} = 0, \ \phi_{m,\hat{\eta}}(.;\delta) \bigg|_{\delta =0} = \phi^0_{m,\hat{\eta}} , \ q_{m,\hat{\eta}}(.;\delta) \bigg|_{\delta =0} = \ q_{0,m,\hat{\eta}}(\delta) \bigg|_{\delta =0} = 0, }$

(ⅱ) $\displaystyle{ ( \lambda_{m,\hat{\eta}}(\delta), \ \phi_{m,\hat{\eta}}(.;\delta), \ q_{m,\hat{\eta}}(.;\delta) , \ q_{0, m,\hat{\eta} }(\delta) \ ) }$ satisfies (16).

(ⅲ) The set $\displaystyle{ \left\{ \phi_{1,\hat{\eta}}(.;\delta), \ldots . \phi_{(d-1),\hat{\eta}}(.;\delta) \right\} }$ is orthonormal in $(L^2 (\mathbb{T}^d))^d.$

(ⅳ) For each interval $I \subset \mathbb{R}$ with $\overline{I}$ containing exactly the eigenvalue $\nu(\hat{\eta})=0$ of (17) (and no other eigenvalue of (17) then $\left\{\lambda_{1,\hat{\eta}}(\delta) \ldots, \lambda_{(d-1),\hat{\eta}}(\delta)\right\}$ are the only eigenvalues of (16 (counting multiplicities) lying in the interval $I.$

The above theorem says that there are $(d-1)$ smooth curves emanating out of the zero eigenvalue as $\delta$ varies in an interval $(-\delta_0, \delta_0)$. We call them Rellich branches. Using them, for $m= 1, \ldots , (d-1)$, we can define the corresponding $m^{th}$ Bloch transform of $g \in \left(L^2 (\mathbb{R}^d)\right)^d$ via the expression

where $\delta =\delta (\epsilon, \xi) = \epsilon |\xi|$ and $\hat{\eta} = \xi/|\xi|.$ This is well defined provided $\epsilon$ is sufficiently small so that $\epsilon |\xi| \leq \delta_0$. For other $\xi$, we define $B^\epsilon_{m,\hat{\eta}} g(\xi) =0.$

For later purposes we need the Bloch transform for $(H^{-1}(\mathbb{R}^d))^d$ elements also. Let us consider $F \equiv ( g^0 + \sum_{j=1}^d \frac{\partial}{\partial x_j}g^j ) \in (H^{-1}(\mathbb{R}^d))^d$, where $F, g^0, g^1,...,g^d$ are valued in $\mathbb{C}^d$ and $g^j\in ((L^2(\mathbb{R}^d))^d$ for $j=0,1,...,d$. Then we define $B^{\epsilon}_{m,\hat{\eta}} F(\xi)$ in $L^2_{loc}(\mathbb{R}^d_{\xi})$ by

Definition (19) is independent of the representation used for $F\in (H^{-1}(\mathbb{R}^d))^d$ in terms of $\{g^j ,\ j=0,...,d\}$ and is consistent with the previous definition (18) whenever $F\in (L^2(\mathbb{R}^d))^d$.

Remark 1. Due to the property $\nabla\cdot(e^{ix\cdot\xi}\phi^\epsilon_{m,\hat{\eta}})=0$ in $\mathbb{R}^d$, we see from (19) that

In fact, by considering $\nabla \psi = g^0 + \sum_{j=1}^d \frac{\partial}{\partial x_j}g^j\in (H^{-1}(\mathbb{R}^d))^d$ then we can take $g_0 = 0$ and $g_j = \psi e_j$ $j=1,\ldots,d$, in (19) to obtain (20). That is, Bloch transform of gradient field is zero. Therefore the kernel of the Bloch transform $B^\epsilon_{m,\hat{\eta}} : L^2(\mathbb{R}^d)^d \mapsto L^2(\mathbb{R}^d)$ contains the closed subspace $\{\nabla \psi : \psi\in H^1(\mathbb{R}^d)\}$ for each $m=1,\ldots,d-1$. Roughly speaking since Bloch waves satisfy incompressibility condition the Bloch transform on gradient field vanish. Thus we may anticipate that the pressure effects may not be captured in the Bloch method. This impression is not correct. Indeed, as shown Section 5, by means of localization via a cut-off function, we manage to keep the pressure term.

Our next result is concerned with the asymptotic behavior of these Bloch transforms as $\epsilon \rightarrow 0$. In the physical space such convergence often modeled as two scale convergence (or, multi scale convergence) can be found in [5,1]. Here we will be studying such convergence in the Fourier space. Since $\phi_{m,\hat{\eta}} (y; 0)$ is a fixed unit vector $(=\phi^0_{m,\hat{\eta}})$ orthogonal to $\hat{\eta}$ and independent of $y$ (see Theorem 2.1), we have

Theorem 2.2. Let $g^\epsilon$ be a sequence in $(L^2(\mathbb{R}^d))^d$ such that its support is contained in a fixed compact set $K \subset \mathbb{R}^d,$ independent of $\epsilon.$ If $g^\epsilon$ converges weakly to $g$ in $(L^2 (\mathbb{R}^d))^d,$ then we have

where $\hat{g}$ denotes the Fourier transform of $g$ and we recall that $\hat{\eta}=\frac{\xi}{|\xi|}$.

Proof. Let us remark that $B^\epsilon_{m,\hat{\eta}} g^\epsilon (\xi)$ is defined for $\epsilon \leq \frac{\delta_0}{M}$ if $|\xi| \leq M.$ We can write

By using Cauchy-Schwarz, the second term on the above right hand side can be estimated by the quantity

where $C_K$ is a constant depending on $K$ but not on $\epsilon$. Recall that $\delta$ is a function of $(\epsilon, \xi),$ namely $\delta =\epsilon |\xi|$. This quantity is easily seen to converge to zero as $\epsilon \rightarrow 0$ for each fixed $\xi$ because of the directional continuity of $\phi_{m,\hat{\eta}}(.,\delta) \mapsto \phi^0_{m,\hat{\eta}}$ in $(L^2(\mathbb{T}^d))^d$ as $\delta\to0$. We merely use the continuity of the $m^{th}$ Rellich branch at $\delta=0$ with values in $(L^2 (\mathbb{T}^d))^d$. On the other hand, thanks to our normalization, the integral on $K$ is bounded by a constant independent of $(\epsilon, \xi).$ The proof is completed by a simple application of the Dominated Convergence Theorem which guarantees that the second term on the above right hand side converges strongly to $0$ in $L^2_{loc} (\mathbb{R}^d_\xi)$ as $\epsilon\to0$.

Since compactly supported elements are dense in $(L^2(\mathbb{R}^d))^d$, we have the following:

Corollary 1. In the setting of Theorem 2.2, if $g^\epsilon$ be a sequence in $(L^2(\mathbb{R}^d))^d$ such that its support is contained in a fixed compact set $K \subset \mathbb{R}^d,$ independent of $\epsilon$ and $g^\epsilon \rightarrow g$ in $L^2(\mathbb{R}^d)^d$ then we have the following strong convergence

We recall the classical orthogonal decomposition [15] :

Let us denote

By our choice, $\{\phi^0_{1,\hat{\eta}},\ldots,\phi^0_{d-1,\hat{\eta}},\hat{\eta}\}$ forms an orthonormal basis in $\mathbb{R}^d$, and so we can deduce the following :

Proposition 1. If $g\in X^{\perp}$ and $\phi^0_{m,\hat{\eta}}\cdot \widehat{g} = 0$ for all $m=1,\ldots, d-1$, then $g=0$.

Proof. The proof is immediate, as $\{\phi^0_{1,\hat{\eta}},\ldots,\phi^0_{d-1,\hat{\eta}}\}$ forms an orthogonal basis in $\mathbb{R}^{d-1}$ and $\phi^0_{m,\hat{\eta}}\cdot \widehat{g} = 0$ for all $m=1,\ldots, d-1$, so $\widehat{g}(\xi) = c(\xi)\xi$ for some scalar $c\in L^2(\mathbb{R}^d)$. Now if $c\neq 0$, it contradicts the hypothesis $g\in X^{\perp}$. Thus $c=0$. Consequently, $g= 0$.

Corollary 2. In the setting of Theorem 2.2 and Proposition 1 if $g^\epsilon$ be a sequence in $X^{\perp} \subset L^2(\mathbb{R}^d)^d$ such that its support is contained in a fixed compact set $K \subset \mathbb{R}^d,$ independent of $\epsilon$ and $g^\epsilon \rightharpoonup g$ in $L^2(\mathbb{R}^d)^d$ weak, and $\phi^0_{m,\hat{\eta}}\cdot \widehat{g} = 0$ for all $m=1,\ldots, d-1$, then $g=0$.

Proof. The proof simply follows as $X^{\perp}$ is a closed subspace of $L^2(\mathbb{R}^d)^d$, so the limit $g\in X^{\perp}$ and the result follows by applying Proposition 1.

Theorem 2.2 and its corollaries give some sufficient conditions and specify the sense in which the first Bloch transform tends to Fourier transform, which is a sign of homogenization on the Fourier side. We do not believe that these conditions are sharp and exhaust all possibilities. It would be interesting to explore further in this direction.

Remark 2. Bloch waves being incompressible are transversal. Longitudinal direction is missing and it has to be added to get the full basis. Naturally, asymptotics of the Bloch transform contains information of the Fourier transform only in transversal directions. It contains no information in the longitudinal direction. Because of this feature, in the homogenization limit also, there is no information in the longitudinal direction. This is however proved to be enough to complete the homogenization process because the limiting velocity field is incompressible. See Section 5.

3.   Computation of derivatives

In this section, we give the expressions of the derivatives (at $\delta=0$) of the Rellich branches $\{\phi_{m,\hat{\eta}}(y;\delta),q_{m,\hat{\eta}}(y;\delta), q_{0,m,\hat{\eta}}(\delta), \lambda_{m,\hat{\eta}}(\delta)\}$ obtained in Theorem 2.2. These results are essentially borrowed from [7] except for the second order derivative of $q_{0,m,\hat{\eta}}(\delta)$ which is new. We differentiate, with respect to $\delta\in\mathbb{R}$, (16) or equivalently the following system, fixing $m=1,\ldots,d-1$ and $\hat{\eta} =\frac{\xi}{|\xi|}\in \mathbb{S}^{d-1}$,

Zeroth order derivatives : For $m = 1, \ldots ,d-1$ and for a fixed direction $\hat{\eta} \in \mathbb{S}^{d -1}$ we have $\lambda_{m,\hat{\eta}} (0) =0$ and a corresponding eigenfunction is such that $q_{m,\hat{\eta}}(y;0) =0, q_{0, m,\hat{\eta}}(0) =0$ and $\phi_{m,\hat{\eta}}(y;0)$ is a constant unit vector of $\mathbb{R}^d$ orthogonal to $\hat{\eta}$. We give a notation for this constant $\phi_{m,\hat{\eta}}(y;0) = \phi_{m,\hat{\eta}}^0$. We recall that $\{ \phi^0_{1,\hat{\eta}},\ldots, \phi^0_{d-1,\hat{\eta}},\hat{\eta}\}$ is such that they form an orthonormal basis for $\mathbb{R}^d$.

First order derivatives : Let us differentiate (25) once with respect to $\delta$ to obtain (prime denotes derivatives with respect to $\delta$) :

where,

We put $\delta =0$ in (26) and by integrating over $\mathbb{T}^d$, we obtain

Taking scalar product with $\hat{\eta}$, we simply get $\lambda^\prime_{m,\hat{\eta}}(0) =q^\prime_{0,m,\hat{\eta}}(0)=0$ as $\hat{\eta}\perp \phi^0_{m,\hat{\eta}}$.

Using the above information in (26), we find that $(\phi^\prime_{m,\hat{\eta}}(y;0),q^\prime_{m,\hat{\eta}}(y;0))$ is a solution of the following cell problem :

Comparing this with (7), it can be seen that that $\phi^\prime_{m,\hat{\eta}}(y;0)$ is given by (see [7]) :

where $\zeta_{m,\hat{\eta}}\in\mathbb{C}^d$ is a constant vector (independent of $y$), orthogonal to $\hat{\eta}$. In other words, the $y$-dependence of $\phi_{m,\hat{\eta}}^\prime (y;0)$ is completely determined by the cell test function $\chi^r_\alpha(y)$, solution of problem (7).

In a similar manner, the derivative of the eigenpressure $q_{m,\hat{\eta}}(y;0)$ is given by (see [7]):

That is, the $y$-dependence of $q_{m,\hat{\eta}}^\prime (y;0)$ is completely determined by the cell test function $\Pi^r_\alpha(y)$, solution of problem (7).

Second order derivatives : Next we differentiate (26) with respect to $\delta$ to obtain :

where

We consider (30) at $\delta=0$ and by integrating over $\mathbb{T}^d$, we get

or,

where $M (\hat{\eta}, A^\ast)$ is the symmetric matrix whose entries are given by

This is nothing but a contraction of the homogenized tensor $A^{*}$. As a simple consequence of (32), we get

It also follows that $M(\hat{\eta},A^{\ast})\phi_{m,\hat{\eta}}^0 \perp \phi^0_{m^\prime,\hat{\eta}}$ for all $m\neq m^\prime$.

By summarizing the above computations, we have

Theorem 3.1. For $m = 1, \ldots ,d-1$ and for a fixed direction $\hat{\eta} \in \mathbb{S}^{d -1}$ we have

(ⅰ) $\lambda_{m,\hat{\eta}} (0) =0$ and a corresponding eigenfunction is such that $q_{m,\hat{\eta}}(y;0) =0, q_{0, m,\hat{\eta}}(0) =0$ and $\phi_{m,\hat{\eta}}(y;0) = \phi_{m,\hat{\eta}}^0$ a unit vector orthogonal to $\hat{\eta}$.

(ⅱ) $\lambda_{m,\hat{\eta}}^\prime (0) =0$ and $q_{0,m,\hat{\eta}}^\prime(0) = 0$.

(ⅲ) The derivative of the eigenfunction $\phi_{m,\hat{\eta}}(y; {\delta})$ at $\delta =0$ satisfies:

where $\zeta_{m,\hat{\eta}}\in\mathbb{C}^d$ is a constant vector (independent of $y$), orthogonal to $\hat{\eta}$.

(ⅳ) The derivative of the eigenfunction $q_{m,\hat{\eta}}(y;{\delta})$ at $\delta =0$ satisfies:

(v) The second derivative of the eigenvalue $\lambda_{m,\hat{\eta}}({\delta})$ and $q_{0,m,\hat{\eta}}(\delta)$ at $\delta =0$ satisfy the relation

where $M (\hat{\eta}, A^\ast)$ is the symmetric matrix whose entries are given by

Remark 3. The above matrix $M (\hat{\eta}, A^\ast)$ is precisely that which must be positive definite in the Legendre-Hadamard definition of ellipticity. A relation analogous to (34) is called "propagation relation" in [14] in the study of linearized elasticity system and it shows how the homogenized tensor $A^{*}$ enters into the Bloch wave analysis. The above relation (34) generalizes the relation (22) in [7].

Remark 4. In the linearized elasticity system, the propagation relation is an eigenvalue relation. Here, relation (34) can again be seen as an eigenvalue problem, posed in the $(d-1)$-dimensional subspace orthogonal to $\hat{\eta}$. More precisely, $1/2 \lambda_{m,\hat{\eta}}^{\prime \prime} (0)$ is an eigenvalue and $\phi_{m,\hat{\eta}}^0$ (which is orthogonal to $\hat{\eta}$) is an eigenvector of the restriction of the matrix $M(\hat{\eta},A^{\ast})$ to the subspace $\hat{\eta}^\perp$. In (34) $1/2 q^{\prime\prime}_{0,m,\hat{\eta}}(0)$ is the Lagrange multiplier corresponding to the constraint that the eigenvalue problem is posed in the $(d-1)$-dimensional subspace orthogonal to $\hat{\eta}$.

Case of Symmetrized gradient : We recall the incompressible elasticity system (10) with the symmetrized gradient introduced in Section 1

where $E (v) = \frac{1}{2}( \nabla v + \nabla^t v ).$

We introduce Bloch waves associated to the Stokes operator defined in (35).

Find $\lambda_s = \lambda_s (\eta) \in \mathbb{R}, \phi_s =\phi_s (\eta) \in H^1 (\mathbb{T}^d)^d , \ \phi_s \neq 0$ and $\Pi_s = \Pi_s (\eta) \in L^2 (\mathbb{T}^d)$ satisfying

As usual $D (\eta) =\nabla_y+ i\eta$ is the shifted gradient operator and the shifted strain rate tensor is defined by :

As earlier, we modify the spectral problem (36) as follows : Find $\lambda_s(\delta) \in \mathbb{R}, \phi_s(.;\delta) \in H^1 (\mathbb{T})^d, \ q_s(.;\delta) \in L^2_0 (\mathbb{T}^d)$ and $q_{0,s}(\delta) \in \mathbb{C}$ satisfying

As before, we can compute directional derivatives of the solution of (37) and prove a result completely analogous to Theorem 3.1. In particular, we will have the following propagation relation : For $m = 1, \ldots d-1$ and for fixed direction $\hat{\eta}\in \mathbb{S}^{d -1}$ the second derivative of the eigenvalue $\lambda_{s,m,\hat{\eta}}({\delta})$ at $\delta =0$ satisfies the relation

where $M (\hat{\eta}, A_s^\ast)$ is the matrix whose entries are given by

4.   Recovery of homogenized tensor from Bloch waves

In the scalar self-adjoint case, it is known that the homogenized matrix is equal to one-half the Hessian of the first Bloch eigenvalue at zero momentum [13]. In the general (non-symmetric) scalar case, treated in [21], it was shown that only the symmetric part of the homogenized matrix is determined by the Bloch spectrum and it is given again by the same one-half of the Hessian of the first Bloch eigenvalue (which exists by virtue of the Krein-Rutman theorem). The fact that only the symmetric part of the homogenized matrix plays a role is not a big surprise since, the homogenized tensor $A^{*}$ being constant, the differential operator

depends only on the symmetric part of $A^{*}$.

In the case of systems, another phenomenon takes place. For example, the linearized elasticity system (in which there are no differential constraints) was treated in [14] where it was recognized that not only Bloch eigenvalues but also Bloch eigenfunctions at zero momentum are needed to determine the homogenized tensor. More precisely, this connection between Bloch eigenvalues and eigenfunctions, on the one hand, and the homogenized tensor, on the other hand, was expressed via a relation called propagation relation in [14] which uniquely determines the homogenized tensor.

In the case of Stokes system, a new phenomenon arises because of the presence of a differential constraint (the incompressibility condition). Even though there is an analogue of the propagation relation (see (34) above), it does not determine uniquely the homogenized tensor. In fact the propagation relation (34) is unaltered if we add a multiple of $I \otimes I$ (where $I$ is the $d \times d$ identity matrix) to the homogenized tensor. The homogenized Stokes operator clearly remains the same under such an addition since it corresponds to adding a gradient of the velocity divergence which vanishes because of the incompressibility constraint. The authors in [7] conjectured that the homogenized Stokes tensor is uniquely characterized by the propagation relation up to the addition of a term $c (I \otimes I)$ (where $c$ is a constant). We prove this assertion in the case of the Stokes system (10) with a symmetrized gradient. For the other Stokes system (1), the homogenized tensor is not uniquely determined by the propagation relation (34). In this section, we investigate this non-uniqueness. Neverheless, we shall prove that for both Stokes systems the homogenized operators (9), and its equivalent for the symmetric gradient case of (10), are uniquely determined.

Our concern now is the following question: to what extent do the Bloch spectral elements determine the homogenized tensor $A^\ast$ via the propagation relation (34)? Since $\lambda^{\prime\prime}_{m,\hat{\eta}}(0), q^{\prime\prime}_{0,m,\hat{\eta}}(0),\phi^0_{m,\hat{\eta}}$ are known from Bloch spectral data, it follows that $M(\hat{\eta},A^{*})\phi^0_{m,\hat{\eta}}$ is uniquely determined via the relation (34). But it may happen that different tensors $A^\ast$ give rise to the same matrix $M(\hat{\eta},A^{*})$. Three main results are proved in this section and they are stated in the following three propositions.

Proposition 2. Let $A^{*}$ and $B^{*}$ be two fourth order tensors possessing the simple symmetry (8). They satisfy the same propagation relation (34), if and only if

where $I$ is the $d \times d$ identity matrix and $N$ is a fourth order tensor satisfying, on top of the simple symmetry (8), the following anti-symmetry property

Proof. Let us observe that the addition of $c (I \otimes I)$ and $N$, having properties (8) and (40), to $A^\ast$ does not alter the propagation relation (34). Indeed, we have,

Since $\phi_{m,\hat{\eta}}^0$ is orthogonal to $\hat{\eta}$, we deduce

Conversely, let us assume that there are two fourth-order tensors $A^\ast$ and $B^\ast$, possessing the simple symmetry (8) and such that $M(\hat{\eta},A^{*}) \phi_{m,\hat{\eta}}^0 = M(\hat{\eta},B^{*}) \phi_{m,\hat{\eta}}^0$, $m=1,...,d-1$, for all $\hat{\eta}\in\mathbb{S}^{d-1}$. We must then deduce (39). For convenience, the proof is divided into five steps.

Step 1. We begin with showing the matrix $M(\hat{\eta}, A^\ast)$ is symmetric. By interchanging the dummy indices $\alpha$ and $\beta$ and using the simple symmetry (8) of the homogenized coefficients, ${(A^{*})}^{jl}_{\alpha\beta} = {(A^{*})}^{lj}_{\beta\alpha}$, we get

which shows the required symmetry.

Step 2. For $\widetilde{N}= B^\ast -A^\ast$ define $M (\hat{\eta}) = M (\hat{\eta},\widetilde{N}) = M (\hat{\eta} , B^\ast) - M (\hat{\eta} , A^\ast)$. Since $A^\ast$ and $B^\ast$ satisfy (34), it follows that $M(\hat{\eta}) \phi^0_{m,\hat{\eta}} = 0$ for $m=1,...,d-1$. Since the family $\phi^0_{m,\hat{\eta}}$ is a basis of the orthogonal space to $\hat{\eta}$, it implies that $M (\hat{\eta}) = c(\hat{\eta}) \hat{\eta} \otimes \hat{\eta}$ for some scalar $c(\hat{\eta})$. Since $M (\hat{\eta})$ depends quadratically on $\hat{\eta}$, it must be that $c(\hat{\eta})$ is independent of $\hat{\eta}$. Thus, for $c\in\mathbb{R}$, we have $M (\hat{\eta}) = c \, \hat{\eta} \otimes \hat{\eta}$, that is, for any $\hat{\eta}\in\mathbb{S}^{d-1}$,

Step 3. Under condition (42), we verify that

For this purpose, let us take $\hat{\eta} = e_i$ in (42). We obtain $\widetilde{N}_{ii}^{jl}= c \delta_{ij} \delta_{il}$ and so

In particular, (43) is proved. Next, choosing $\hat{\eta} =e_{i} +e_k$ in (42), we get

To check (44), there are several cases to consider.

(ⅰ) $(i\neq j$ and $k\neq j)$. In this case, (44) is a direct consequence of (46) and (47).

(ⅱ) Similarly, for $(k\neq l$ and $i\neq l)$ 44) is a direct consequence of (46) and (47).

(ⅲ) $(i \neq j, \ k=j ).$ In this case,

Now together with $i \neq l$ we have

Then both $j=l$ or $j\neq l$ cases lead to verify (43) and (44) respectively.

(ⅳ) Similarly, for $(k\neq l$ and $i = l)$

Together with $k \neq j$ we have

Then both $i=j$ or $i\neq j$ cases lead to verify (43) and (44) respectively.

Step 4. Now we consider the two remaining cases not covered in (44).

(ⅰ) $(i, k) = (j, l)$. Then from (47) we have

For $i\neq k$ it gives using (46)

(ⅱ) Similarly, for $(k,i) = (j,l)$, together with $i \neq k$ we have

Step 5. Let us set $N = \widetilde{N} - c (I \otimes I).$ Thanks to the properties (43) and (44), we can easily check that $N$ is an anti-symmetric tensor in the sense that it satisfies

From its very definition $N$ also possesses the symmetry $N_{ik}^{jl} = N_{ki}^{lj}.$ Thus $N$ has all the properties listed in (40).

Next we extend Proposition 2 to the Stokes system (10), featuring a symmetric gradient tensor. In this case the propagation relation (34) is replaced by (38) and the homogenized tensor is denoted by $A^{*}_s$.

Proposition 3. The propagation relation (38) characterizes uniquely the tensor $A_s^\ast$, up to the addition of a constant multiple of $I \otimes I$. In other words, $A^{*}_s$ and $B^{*}_s$ satisfy the same propagation relation (38) if and only if, for some $c\in \mathbb{R}$,

Proof. The proof continues from the Step 5 of the previous proof of Proposition 2. We defined $N = \widetilde{N} - c(I\otimes I)$ satisfying (54) i.e.

Now as $\widetilde{N}= B^{*}_s - A^{*}_s$ possess with the symmetry of coefficients of linear elasticity, so we have

This symmetry combined with the anti-symmetry established in the previous step implies that $N =0.$ Note that antisymmetry property holds precisely for the interchange of those pairs of indices for which symmetry property does not hold.

This can be seen as follows: whenever $(i, k) \neq (j, l) \mbox{ and } (k, i ) \neq (j, l)$

Similarly, whenever $(i, k) = (j, l) \mbox{ or } (k, i ) = (j, l)$ together with $i\neq k$; from (52), (53) we have

Then using (56) and (46) we clearly have

Therefore (58), (59) imply that $N =0$ or, $\widetilde{N} = c (I\otimes I)$ and hence $B_s^\ast -A_s^\ast = c (I \otimes I).$

Remark 5. The conclusion of the above proposition was conjectured in [7] and it is proved here to be true whenever we are working with the system (10) with symmetrized gradient. However, it is not true with the full gradient Stokes system (1) as shown by Proposition 2. However, in both of these cases the propagation relation fixes the homogenized operator(9) uniquely, as is stated in the following proposition.

Proposition 4. If (39) is satisfied, then $A^{*}$ and $B^{*}$ give rise to the same homogenized operator (9).

Proof. We have to check that $A^{*}$ and $B^{*}$ define the same Stokes differential operator for divergence-free vector fields. Indeed the Fourier symbol of the operator

is $(A^{*}-B^{*})^{kl}_{\alpha \beta}\xi_\alpha\xi_\beta$ which, by virtue of (42), is equal to $c\xi_k\xi_l$ which is precisely the symbol of the operator $u \to -c \nabla (\nabla\cdot u)$ which vanishes on the space of divergence free functions.

5.   Homogenization result

This section is devoted to a proof of Theorem 1.1, our main homogenization result stated in the first section. It is based on the tools that we have introduced so far. A similar proof is given for the linear elasticity problem in [21]. However, the presence of a pressure and a differential constraint in the Stokes system seriously complexifies the analysis and has a non-trivial effect in the homogenization process. Besides, we also bring some simplifications to the proof given in [21].

We consider a sequence of solutions $(u^{\epsilon},p^{\epsilon})\in (H^1_{0}(\Omega))^d\times L^2_0 (\Omega)$ solving the Stokes system (1). It is classical to derive the following bound [8] :

where $C$ is independent of $\epsilon.$ Then there exist $(u,p)\in (H^1_{0}(\Omega))^d\times L^2_0 (\Omega) \ $ and a subsequence ($u^{\epsilon}$, $p^{\epsilon}$) converging weakly to $(u, p)$ in $(H^1_{0}(\Omega))^d \times L^2_0(\Omega)$. Our aim is to show that $(u,p)$ satisfies the homogenized Stokes system (9). Due to the uniqueness of solutions for the system (9), it follows that the entire sequence $(u^{\epsilon},p^{\epsilon})$ converges to $(u,p)$ weakly in $(H^1_0 (\Omega))^d \times L^2_0 (\Omega).$

There are several steps in the proof. First, we localize the Stokes system (1) by applying a cut-off function technique to the velocity $u$ in order to get the equation (61) in the whole $\mathbb{R}^d$. Next, by taking the Bloch transformation $B_{m,\hat{\eta}}^{\epsilon}$ $(1\leq m \leq d-1 )$ of the equation (61) and passing to the limit, we arrive at the homogenized equation in the Fourier space. Finally, we take the inverse Fourier transform to go back to the physical space which gives our desired result.

Notation. in the sequel L.H.S. stands for left hand side, and R.H.S. for right hand side.

Step 1. Localization of the velocity $u$ : Let $v \in \mathcal{D}(\Omega)$ be arbitrary. Then $vu^{\epsilon}$ and $p^{\epsilon}$ satisfy (for $l=1, \ldots , d)$

where,

Note that, $g_l^{\epsilon}$ and $h_l^{\epsilon}$ correspond to terms containing zero and first order derivatives of $\mu^{\epsilon}$ respectively. In the sequel, we extend $u^\epsilon$ and $p^\epsilon$ by zero outside $\Omega$ and such extensions are denoted by the same letters.

Step 2. Limit of $B_{m,\hat{\eta}}^{\epsilon}$ applied to the L.H.S. of (61) : We consider the following $\epsilon$-scaled spectral problem of (25) as follows : Let $\hat{\eta}= \frac{\xi}{|\xi|} \in \mathbb{S}^{d-1}$, $\delta = \epsilon(\xi \cdot\hat{\eta})$;

They satisfy the following system because of (25) :

Let us first consider the L.H.S. of (61). For $g\in H^1(\mathbb{R}^d)^d$ with compact support in $\Omega$, using the definition Bloch transformation (19) and spectral equation (63), we obtain for $m= 1, \ldots d-1,$

In the previous equation the duality bracket is between $H^1_{comp}(\mathbb{R}^d)^d$ and $H^{-1}_{loc}(\mathbb{R}^d)^d$.

Therefore, $B_{m,\hat{\eta}}^{\epsilon}$ applied to the L.H.S. of (61) ($1\leq m \leq d-1$) is equal to

Below, we treat each term of (64) one by one.

1st term of (64) : By using the Taylor expansion

and then using Theorem 2.2, we get

where we recall that $\phi_{m,\hat{\eta}}^{0}$ is a constant unit vector of $\mathbb{R}^d$ orthogonal to $\hat{\eta}$. Note that $\lambda_{m,\hat{\eta}}^{\prime\prime}(0)$ is linked to $A^{*}$ via the propagation relation (34). Using this relation, the above limit can be written as

2nd term of (64) :

Using the Taylor expansion of $q^\epsilon_{m,\hat{\eta}}(.;\delta)$ :

(prime denotes the derivative with respect to the second variable), with the properties that (cf. Theorem 2.1)

where, $M_{\mathbb{T}^d}(f) = \int_{\mathbb{T}^d} f(y)\ dy$.

Then by using $u^\epsilon\rightarrow u$ strongly in $L^2(\Omega)^d$ from (68) we get

It is also used that, the error term $\mathcal{O}(\epsilon(\xi\cdot\hat{\eta})^2)$ in the above Taylor expansion tends to $0$ in the space $L^2_{loc} (\mathbb{R}^d_{\xi}; L^2_{loc} (\mathbb{R}^d))$. Thus the oscillating eigen-pressure $q^\epsilon_{m,\hat{\eta}}$ does not contribute to the homogenized system.

3rd term of (64) : We use the Taylor expression of $q_{0,m,\hat{\eta}}^\epsilon(\xi)$ with the property $q_{0,m,\hat{\eta}}(0)=q^{\prime}_{0,m,\hat{\eta}}(0)=0$ (cf. Theorem 3.1) to have

So,

4th term of (64) : Finally, we consider the remaining fourth term in (64), and doing integration by parts we get

We use the Taylor expansion

And from (60) as $||p^\epsilon||_{L^2(\Omega)}$ is uniformly bounded, so up to a subsequence we have

Thus by passing to the limit in the R.H.S. of (74), we get

Thus

This property proved for $H^{-1}$ elements is analogous to Theorem 2.1.

Summary so far : Combining the previous results, therefore, by taking the Bloch transformation $B_{m,\hat{\eta}}^{\epsilon}$ of the L.H.S. of (61) ($1\leq m \leq d-1$) and multiplying by $\chi_{\epsilon^{-1}\mathbb{T}^d}$, we see that it converges to

Step 3. Limit of $B_{m,\hat{\eta}}^{\epsilon}$ applied to the R.H.S. of (61) : Applying $B_{m,\hat{\eta}}^{\epsilon}$ to the R.H.S. of (61) ($1\leq m \leq d-1$), we obtain

We treat below each of these terms separately. Passing to the limit in the first term is straightforward (cf. Corollary 1) and we obtain

Limit of $B^{\epsilon}_{m,\hat{\eta}}(g^{\epsilon})$ : We pose $\sigma^{\epsilon}= \mu^{\epsilon}\nabla u^{\epsilon}$ ($\sigma^\epsilon_{l\alpha} = \mu^\epsilon\frac{\partial u^\epsilon_l}{\partial x_\alpha}$) which is a bounded matrix in $(L^2(\Omega))^{d \times d}$ and so there exists a weakly convergent subsequence in\break $(L^2(\Omega))^{d \times d }$. Let $\sigma$ be its limit as well as its extension by zero outside $\Omega$. Then via Theorem 2.2,

Due to the strong convergence of $u^\epsilon$ in $L^2 (\mathbb{R}^d)^d,$ (cf. Corollary 1) we have

Combining the above two convergence results and doing integration by parts, we obtain

Limit of $B^{\epsilon}_{m,\hat{\eta}}(h^{\epsilon})$ : We decompose it into two terms:

We start with the second term. By doing integration by parts, it becomes

Thanks to the strong convergence of $u^\epsilon$ in $L^2 (\mathbb{R}^d)^d$, the above quantity converges in $L^2_{loc}(\mathbb{R}^d_\xi)$ strongly to

Next, we consider the first term of the R.H.S. of (85). After doing integration by parts, one has

In a manner similar to the above arguments, the limit of (87) would be

Now combining (86) and (88) and using the fact

is a constant vector of $\mathbb{C}^d$ independent of $y$, which in turn implies that

we see that $\chi_{{\epsilon}^{-1}\mathbb{T}^{d}}B_{m,\hat{\eta}}^{\epsilon}(h^\epsilon)(\xi)$ converges strongly in $L^2_{loc}(\mathbb{R}^d_\xi)$ to

Step 4. Limit of $B_{m,\hat{\eta}}^{\epsilon}$ applied to (61) : By equating the limiting identities that we have derived in the last two steps, we obtain

The above equation has to be considered as the localized homogenized equation in the Fourier space. The conclusion of Theorem 1.1 will follow as a consequence of this equation.

Step 5. Passage from Fourier space ($\xi)$ to physical space $(x)$ : We note that the L.H.S. and the R.H.S. of (90) can be written as $L(\xi)\cdot\phi_{m,\hat{\eta}}^0$ and $R(\xi)\cdot\phi_{m,\hat{\eta}}^0$, respectively, so that we have

Observe that, the quantity $\left[L(\xi)-R(\xi)\right]$ is independent of $m$. Varying $m=1,..,(d-1)$ and using the fact $\xi \perp \phi_{m,\hat{\eta}}^0, \xi\in\mathbb{R}^d$ and $\{\phi^0_{1, \hat{\eta}} \cdots \phi^0_{d-1, \hat{\eta}}\}$ forms a basis of $\mathbb{R}^{d-1}$, we get

Therefore, for all test functions $w \in (L^2(\mathbb{R}^d))^d$ satisfying $\xi\cdot\hat{w}(\xi)=0$ (i.e. $\text{div}\,w=0$ in $\mathbb{R}^d$) we also have

Now by using the Plancherel's theorem, we have

where $\mathcal{F}^{-1}$ denotes the inverse Fourier transformation.

We easily compute $I (x) = \mathcal{F}^{-1}\left[L(\xi)-R(\xi)\right](x)$ to obtain

which simplifies in

We pose

to write $I_{l}$ in the form

Using (91), it follows from de Rham's theorem that $I$ is a gradient and furthermore this is true whatever be $v\in\mathcal{D}(\Omega)$. This imposes restriction on $F^1,F^2$. In fact, we show using (91) that $F^2_{l\beta}= q\delta_{l\beta}$ and $F^1 = \nabla q$ for some scalar $q\in L^2(\Omega)$ so that $I=v\nabla q + q\nabla v = \nabla (vq)$.

Step 5A. To show $F^2_{l\beta} = q\delta_{l\beta}$ : Let us choose $v = v_0 e^{inx\cdot \omega}$, where $\omega$ is a unit vector in $\mathbb{R}^d$ and $v_0 \in \mathcal{D}(\Omega)$ is fixed. Next, we choose $w = \psi_{\zeta,\omega} \in (L^2(\mathbb{R}^d))^d$ where for any two constant perpendicular vectors $\zeta$ and $\omega$ in $\mathbb{R}^d$, $\psi_{\zeta,\omega}\in (L^2(\mathbb{R}^d))^d$ solves

The existence of such a function $\psi_{\zeta,\omega}$ can be shown as follows. Let $R_0>0$ be such that $\overline{\Omega} \subset B(0,R_0)$ and consider the following boundary value problem

There exists a solution of (94) (see [15,Page No. 24]) since the boundary data satisfies the required compatibility condition (recall that we assume $\zeta\cdot \omega = 0$)

Then extending $\psi_{\zeta,\omega}$ by $0$ outside $B(0,R_0)$ and by $\zeta e^{-inx\cdot\omega}$ in $\Omega$, clearly the extended function $\psi_{\zeta,\omega}$ solves (93).

Now using these $v$ and $w$ in (91), we have

and dividing by $n$ and letting $n\rightarrow \infty$ in the above relation, we get

As $v_0\in \mathcal{D}(\Omega)$ is arbitrary, (95) gives $F^2\ \omega\cdot \zeta = 0$ in $\Omega$. As $F^2$ is symmetric, and further using that $\omega, \zeta$ are arbitrary satisfying $\omega\cdot\zeta=0$, we conclude $F^2_{l\beta} = F^2_{\beta l} = q\delta_{l\beta}$ for some scalar function $q\in L^2(\Omega)$. This means that we have the relation :

Step 5B. To show $F^1 =\nabla q$ : We choose $v\in\mathcal{D}(\Omega)$ and $w= \psi_{e_k,0}$ with $\psi_{e_k,0}$ as in (93) with $\zeta=e_k$ and $\omega=0$. Then using these $v$ and $w$ in (93) and using the conclusion from Step 5A, we have

which implies $(F^1_k -\frac{\partial q}{\partial x_k}) =0$ for $k=1,..,d$ or, $F^1=\nabla q$.

Step 5C Using Step 5A and Step 5B in (92), and considering the relation $F^1 -\nabla F^2 = 0$ in $\Omega$, we get the macro balance equation :

Step 5D In this step, we prove that $q=0$ in $\Omega$ by using the divergence-free condition. Indeed, as $\nabla\cdot u^\epsilon = 0$ in $\Omega$, we have

Passing to the limit $\epsilon \rightarrow 0$, we get

Using this relation in (96) with $\beta=l$, we get

On the other hand, from (6) and (7) we have

Thus for fixed $k,\alpha=1,\ldots,d$ summing over $l$, since $\text{div} \chi^k_\alpha = 0$ in $Y$, we obtain

Using (98) and (99), as $\text{div} u =0$, we deduce

Finally, the macro constitutive law follows as a consequence from (96) :

Step 5E. Since $q=0$, we deduce from Step 5B that $F^1=0$ and from (96) we get the following homogenized Stokes system satisfied by $u,p$ :

This completes the proof of Theorem 1.1.

Acknowledgments

This work has been carried out within a project supported by Indo -French Centre for Applied Maths -UMI, IFCAM. G. A. is a member of the DEFI project at INRIA Saclay Ile-de-France. M.V. is in the current address during the preparation of the revised version.

The authors thank the referees for their useful comments to improve the presentation of this paper.