# Full asymptotic expansion of the permeability matrix of a dilute periodic porous medium

F Feppon

► To cite this version:

F Feppon. Full asymptotic expansion of the permeability matrix of a dilute periodic porous medium. 2024. hal-04588126

HAL Id: hal-04588126

https://hal.science/hal-04588126

Preprint submitted on 25 May 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution - ShareAlike 4.0 International License
# Full asymptotic expansion of the permeability matrix of a dilute periodic porous medium

F Feppon

► To cite this version:

F Feppon. Full asymptotic expansion of the permeability matrix of a dilute periodic porous medium. 2024. hal-04588126v2

HAL Id: hal-04588126

https://hal.science/hal-04588126v2

Preprint submitted on 20 Nov 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
FULL ASYMPTOTIC EXPANSION OF THE PERMEABILITY MATRIX
OF A DILUTE PERIODIC POROUS MEDIUM

F. FEPPON¹, *

¹ NUMA Unit, Department of Computer Science, KU Leuven, Belgium.

**ABSTRACT.** We compute full asymptotic expansions of the drag force generated by a laminar fluid flowing through a periodic array of small solid particles. Our derivation considers obstacles with arbitrary shape in arbitrary space dimension. Equivalently, we compute the full asymptotic expansion of the effective permeability matrix of the medium with respect to the size of the obstacles. In the first step of our analysis, we use hydrodynamics layer potential theory to obtain the asymptotic expansion of the velocity and pressure fields of the flow across the periodic array. The terms of these expansions can be computed through a procedure requiring to solve a cascade of exterior and interior problems. The second step of our derivation is to integrate this expansion over the perforated unit cell. The derivation requires to compute some Hadamard finite part integrals and brings into play specific tensors involving integrals of the fundamental solution and its derivatives over the faces of the unit cell. We verify that our expansions agree to the leading order with the expressions found by Hasimoto [23] in the case of spherical obstacles in two and three dimensions.

**Keywords.** Porous media, Asymptotic analysis, Hydrodynamics layer potentials, Integral equations, Darcy's law.

**AMS Subject classifications.** 76M45, 35B40, 35B10, 45M05

CONTENTS

<table><tr><td>1. Introduction</td><td>1</td></tr><tr><td>2. Preliminaries</td><td>4</td></tr><tr><td>2.1. Hydrodynamics layer potentials</td><td>5</td></tr><tr><td>2.2. Exterior Stokes problems and layer potential representations</td><td>7</td></tr><tr><td>2.3. Stokes systems with periodic boundary conditions</td><td>10</td></tr><tr><td>2.4. Periodic Green matrix</td><td>12</td></tr><tr><td>3. Full asymptotic expansions of the solutions to the perforated periodic Stokes system</td><td>16</td></tr><tr><td>4. Full asymptotic expansions of the permeability matrix</td><td>23</td></tr><tr><td>Appendix A. High-order error estimates without potential theory</td><td>28</td></tr><tr><td>Appendix B. Some box integrals identities</td><td>32</td></tr><tr><td>Appendix C. Moments of the unit hypersphere</td><td>33</td></tr><tr><td>References</td><td>34</td></tr></table>

1. INTRODUCTION

In this work, we propose a complete asymptotic analysis of a Stokes flow past a dilute periodic porous medium and of its induced drag force, in any dimension $d \ge 2$ and for arbitrarily shaped obstacles. More precisely, we consider a steady incompressible laminar fluid flowing into a $d$-dimensional periodic array filled with small isolated obstacles. The periodic lattice and its perforated periodic unit cell $P\setminus(\eta T)$ are represented on Figure 1. $P = (-1/2, 1/2)^d$ is the $d$-dimensional unit hypercube and $T \subset \mathbb{R}^d$ is a smooth bounded set with $C^\infty$ boundary such that $P\setminus(\eta T)$ is connected. The positive number $\eta > 0$ is a small rescaling parameter which can be thought of as the radius of an obstacle $\eta T$. We consider a dilute suspension of particles, namely the low volume fraction regime in which the size of the obstacles is much smaller than the lattice spacing:

$$ \eta \to 0. $$

Due to the periodicity of the lattice, it is natural to assume that the velocity of the fluid and its pressure gradient $\nabla p$ are 1-periodic in the $d$ space directions. Classically, the periodicity of the pressure gradient $\nabla p$ implies that $p$ can be decomposed as the sum of a linear and a periodic term:

$$ p(x) = \overline{\nabla P} \cdot x + p_{\text{per}}(x) \text{ with } p_{\text{per}}(x + e_i) = p_{\text{per}}(x) \text{ for all } 1 \le i \le d, $$

where $\overline{\nabla P} = (\overline{\nabla P}_i)_{1\le i \le d} \in \mathbb{R}^d$ is a constant vector, $p_{\text{per}}$ is a periodic scalar field and $(e_i)_{1\le i \le d}$ stands for the canonical orthonormal basis of $\mathbb{R}^d$ (the axes of the periodic lattice). The motion of the free fluid inside the

<footnote>* Corresponding author. Email: florian.feppon@kuleuven.be.</footnote>
FULL ASYMPTOTIC EXPANSION OF THE PERMEABILITY MATRIX
OF A DILUTE PERIODIC POROUS MEDIUM

F. FEPPON¹, *

¹ NUMA Unit, Department of Computer Science, KU Leuven, Belgium.

**ABSTRACT.** We compute full asymptotic expansions of the permeability matrix of a laminar fluid flowing through a periodic array of small solid particles. Our analysis applies to obstacles of arbitrary shape and arbitrary space dimensions. In the first step, we use hydrodynamics layer potential theory to obtain the asymptotic expansion of the velocity and pressure fields across the periodic array. The terms of these expansions can be computed through a procedure involving a cascade of exterior and interior problems. In the second step, we deduce the full asymptotic expansions of the permeability matrix. The derivation requires evaluating Hadamard finite part integrals and tensors depending on the values of the fundamental solution or its derivatives on the faces of the unit cell. We verify that our expansions agree to the leading order with the expressions found by Hasimoto [24] in the case of spherical obstacles in two and three dimensions.

**Keywords.** Porous media, Asymptotic analysis, Hydrodynamics layer potentials, Integral equations, Darcy's law.

**AMS Subject classifications.** 76M45, 35B40, 35B10, 45M05

CONTENTS

<table><tr><td>1. Introduction</td><td>1</td></tr><tr><td>2. Preliminaries</td><td>4</td></tr><tr><td>2.1. Hydrodynamics layer potentials</td><td>5</td></tr><tr><td>2.2. Exterior Stokes problems and layer potential representations</td><td>7</td></tr><tr><td>2.3. Stokes systems with periodic boundary conditions</td><td>11</td></tr><tr><td>2.4. Periodic Green matrix</td><td>12</td></tr><tr><td>3. Full asymptotic expansions of the solutions to the perforated periodic Stokes system</td><td>16</td></tr><tr><td>4. Full asymptotic expansions of the permeability matrix</td><td>23</td></tr><tr><td>Appendix A. High-order error estimates without potential theory</td><td>28</td></tr><tr><td>Appendix B. Some box integrals identities</td><td>32</td></tr><tr><td>Appendix C. Moments of the unit hypersphere</td><td>33</td></tr><tr><td>References</td><td>34</td></tr></table>

1. INTRODUCTION

In this work, we propose a complete asymptotic analysis of a Stokes flow past a dilute periodic porous medium and of its induced drag force, in any dimension $d \ge 2$ and for arbitrarily shaped obstacles. More precisely, we consider a steady incompressible laminar fluid flowing into a $d$-dimensional periodic array filled with small isolated obstacles. The periodic lattice and its perforated periodic unit cell $P\setminus(\eta T)$ are represented on Figure 1. $P = (-1/2, 1/2)^d$ is the $d$-dimensional unit hypercube and $T \subset \mathbb{R}^d$ is a smooth bounded set with $C^\infty$ boundary such that $P\setminus(\eta T)$ is connected. The positive number $\eta > 0$ is a small rescaling parameter which can be thought of as the radius of an obstacle $\eta T$. We consider a dilute suspension of particles, namely the low volume fraction regime in which the size of the obstacles is much smaller than the lattice spacing:

$$ \eta \to 0. $$

Due to the periodicity of the lattice, it is natural to assume that the velocity of the fluid and its pressure gradient $\nabla p$ are 1-periodic in the $d$ space directions. Classically, the periodicity of the pressure gradient $\nabla p$ implies that $p$ can be decomposed as the sum of a linear and a periodic term:

$$ p(x) = \overline{\nabla P} \cdot \mathbf{x} + p_{\text{per}}(x) \text{ with } p_{\text{per}}(x + e_i) = p_{\text{per}}(x) \text{ for all } 1 \le i \le d, $$

<footnote>* Corresponding author. Email: florian.feppon@kuleuven.be.</footnote>

<footnote>**Data Availability** Data sharing is not applicable to this article as no datasets were generated or analysed during the current study.</footnote>

<footnote>**Conflict of interest** The author declares he has no conflict of interest.</footnote>
<figure><img src="image_2.png" /><figcaption>FIGURE 1. Physical setting: a laminar fluid with average velocity **U** is flowing through an array of small isolated solid obstacles ηT.</figcaption></figure>

array is driven by the opposite of the mean pressure gradient $-\nabla \overline{P}$ which determines the velocity $\mathbf{u}_{\text{per}}$ and the periodic component $p_{\text{per}}$ of the pressure as the solutions to the following Stokes problem in the perforated periodic unit cell $P\setminus(\eta T)$:

$$ 
\left\{ 
\begin{array}{l} 
-\Delta \mathbf{u}_{\text{per}} + \nabla p_{\text{per}} = -\overline{\nabla P} \text{ in } P\setminus(\eta T), \\ 
\text{div}(\mathbf{u}_{\text{per}}) = 0 \text{ in } P\setminus(\eta T), \\ 
\mathbf{u}_{\text{per}} = 0 \text{ on } \partial(\eta T), \\ 
(\mathbf{u}_{\text{per}}, p_{\text{per}}) \text{ is P-periodic.} 
\end{array} 
\right. 
\qquad (1.1) 
$$

By linearity, $(\mathbf{u}_{\text{per}}, p_{\text{per}})$ can be decomposed as

$$ \mathbf{u}_{\text{per}}(x) = - \sum_{j=1}^{d} \overline{\nabla P_j} \boldsymbol{\chi}_j(x), \quad p_{\text{per}} = - \sum_{j=1}^{d} \overline{\nabla P_j} \alpha_j(x), \qquad (1.2) $$

where for all $1 \le j \le d$, $(\boldsymbol{\chi}_j, \alpha_j)$ is the solution to the periodic cell problem

$$ 
\left\{ 
\begin{array}{l} 
-\Delta \boldsymbol{\chi}_j + \nabla \alpha_j = e_j \text{ in } P\setminus(\eta T), \\ 
\text{div}(\boldsymbol{\chi}_j) = 0 \text{ in } P\setminus(\eta T), \\ 
\boldsymbol{\chi}_j = 0 \text{ on } \partial(\eta T), \\ 
(\boldsymbol{\chi}_j, \alpha_j) \text{ is P-periodic.} 
\end{array} 
\right. 
\qquad (1.3) 
$$

By integrating (1.2) over $P\setminus(\eta T)$, it is found that the mean flow velocity $\mathbf{U} := \int_{P\setminus(\eta T)} \mathbf{u}_{\text{per}}(x) dx$ is related to the mean pressure gradient $\overline{\nabla P}$ through the identity

$$ \mathbf{U} = -\boldsymbol{\chi}^* \overline{\nabla P}, \qquad (1.4) $$

where $\boldsymbol{\chi}^*$ is the so-called “permeability matrix”, whose entries are defined by

$$ \boldsymbol{\chi}_{ij}^* := \int_{P\setminus(\eta T)} \boldsymbol{\chi}_j \cdot \mathbf{e}_i \, dx, \quad 1 \le i, j \le d. \qquad (1.5) $$

Equation (1.4) can be interpreted as the Darcy's law governing the apparent porous medium constituted by the periodic lattice. It also shows that the asymptotic expansion of the viscous force field $\mathbf{F} = -\overline{\nabla P}$ for a given prescribed average velocity $\mathbf{U}$ as $\eta \to 0$ can be obtained from the asymptotic expansion of the permeability matrix $\boldsymbol{\chi}^*$, whose complete derivation the objective of this paper.

The asymptotic analysis of the drag force or permeability matrix of Stokes flows past periodic arrays of small solid obstacles has a rather long story over the last 70 years. It is of relevance in industrial applications involving fluids flowing through suspensions of solid particles [35, 19, 10], but also in the quantitative homogenization of the steady (Navier-)Stokes equations with the method of two-scale expansions [33, 2].

Although there had already been attempts by Burgers [8], Brinkman [7] and Debye and Bueche [12] to propose physical models describing the fluid properties in porous media, Hasimoto [23] was the first to derive a mathematically correct asymptotic expansion for the drag force induced by a periodic array of spheres, in two
<figure><img src="image_2.png" /></figure>

<figure><img src="image_3.png" /><figcaption>(B) A fluid flowing through a porous medium filled with a periodic array $\mathbb{Z}^d + \eta T$ of solid particles $\eta T$ of radius $\eta$.</figcaption></figure>

<caption>FIGURE 1. Physical setting: a laminar fluid with average velocity $\mathbf{U}$ is flowing through an array of small isolated solid obstacles $\eta T$.</caption>

where $\overline{\nabla P} = (\overline{\nabla P}_i)_{1\le i \le d} \in \mathbb{R}^d$ is a constant vector, $p_{\text{per}}$ is a periodic scalar field and $(e_i)_{1\le i \le d}$ stands for the canonical orthonormal basis of $\mathbb{R}^d$ (the axes of the periodic lattice). The motion of the free fluid inside the array is driven by the opposite of the mean pressure gradient $-\overline{\nabla P}$, which determines the velocity $\mathbf{u}_{\text{per}}$ and the periodic component $p_{\text{per}}$ of the pressure as the solutions to the following Stokes problem in the perforated periodic unit cell $P\setminus(\eta T)$:

$$ 
\left\{
\begin{array}{l}
-\Delta \mathbf{u}_{\text{per}} + \nabla p_{\text{per}} = -\overline{\nabla P} \text{ in } P\setminus(\eta T), \\
\text{div}(\mathbf{u}_{\text{per}}) = 0 \text{ in } P\setminus(\eta T), \\
\mathbf{u}_{\text{per}} = 0 \text{ on } \partial(\eta T), \\
(\mathbf{u}_{\text{per}}, p_{\text{per}}) \text{ is P-periodic.}
\end{array}
\right.
\qquad (1.1)
$$

By linearity, $(\mathbf{u}_{\text{per}}, p_{\text{per}})$ can be decomposed as

$$ \mathbf{u}_{\text{per}}(x) = - \sum_{j=1}^{d} \overline{\nabla P_j} \mathbf{X}_j(x), \quad p_{\text{per}} = - \sum_{j=1}^{d} \overline{\nabla P_j} \alpha_j(x), \qquad (1.2) $$

where for all $1 \le j \le d$, $(\mathbf{X}_j, \alpha_j)$ is the solution to the periodic cell problem

$$ 
\left\{
\begin{array}{l}
-\Delta \mathbf{X}_j + \nabla \alpha_j = e_j \text{ in } P\setminus(\eta T), \\
\text{div}(\mathbf{X}_j) = 0 \text{ in } P\setminus(\eta T), \\
\mathbf{X}_j = 0 \text{ on } \partial(\eta T), \\
(\mathbf{X}_j, \alpha_j) \text{ is P-periodic.}
\end{array}
\right.
\qquad (1.3)
$$

By integrating (1.2) over $P\setminus(\eta T)$, it is found that the mean flow velocity $\mathbf{U} := \int_{P\setminus(\eta T)} \mathbf{u}_{\text{per}}(x) dx$ is related to the mean pressure gradient $\overline{\nabla P}$ through the identity

$$ \mathbf{U} = -\mathcal{X}^* \overline{\nabla P}, \qquad (1.4) $$

where $\mathcal{X}^*$ is the so-called “permeability matrix”, whose entries are defined by

$$ \mathcal{X}_{ij}^* := \int_{P\setminus(\eta T)} \mathbf{X}_j \cdot \mathbf{e}_i \, dx, \quad 1 \le i, j \le d. \qquad (1.5) $$

Equation (1.4) can be interpreted as the Darcy's law governing the apparent porous medium constituted by the periodic lattice. It also shows that the asymptotic expansion of the viscous force field $\mathbf{F} = -\overline{\nabla P}$ for a given prescribed average velocity $\mathbf{U}$ as $\eta \to 0$ can be obtained from the asymptotic expansion of the permeability matrix $\mathcal{X}^*$, whose complete derivation the objective of this paper.

The asymptotic analysis of the drag force or permeability matrix of Stokes flows past periodic arrays of small solid obstacles has a rather long story over the last 70 years. It is of relevance in industrial applications involving fluids flowing through suspensions of solid particles [39, 20, 10], but also in the quantitative homogenization of the steady (Navier-)Stokes equations with the method of two-scale expansions [37, 2, 27, 31, 28].
and three dimensions. His derivation was based on the computation of an approximation of the fluid velocity by a linear combination of the periodic Stokes fundamental solution and its derivatives. With the present notation conventions, Hasimoto found

$$ \chi^* = \left( \frac{1}{6\pi\eta} - \frac{c}{6\pi} + \frac{2}{9}\eta^2 - \frac{1}{6\pi} \left( \frac{16\pi^2}{45} + 630b^2 \right) \eta^6 + O(\eta^8) \right) \text{I}, \quad \text{if } d=3, \qquad (1.6) $$

where $c \approx 2.8373$ and $b \approx 0.19457$ and I is the identity matrix, and

$$ \chi^* = \left( -\frac{1}{4\pi} \log(\eta) - \frac{\tilde{c}}{4\pi} + \frac{1}{4}\eta^2 + O(\eta^4) \right) \text{I}, \quad \text{if } d=2, \qquad (1.7) $$

with $\tilde{c} \approx 1.3105$. Hasimoto noted that (1.6) can be understood as a higher order correction of the classical formula $F = 6\pi\eta U$ describing the force that acts on a single isolated sphere of radius $\eta$. This reflects the fact that in the dilute regime $\eta \to 0$, the obstacles behave at first order as isolated. Sangani and Acrivos subsequently refined Hasimoto's procedure to obtain higher order terms [34], however their method remained applicable only to spherical obstacles.

Periodic arrays with obstacles of arbitrary shape were first considered in the work of Allaire [1] who identified the first term of the asymptotic expansion of the effective permeability matrix in arbitrary dimension $d \ge 2$. He obtained

$$ \chi^* \sim \begin{cases} \eta^{2-d} M^{-1} & \text{if } d \ge 3, \\ -\frac{1}{4\pi} \log(\eta) \text{I} & \text{if } d = 2, \end{cases} \qquad (1.8) $$

where, in the case $d \ge 3$, $M^{-1}$ is the Brinkmann inverse permeability which is defined from the solutions to exterior Stokes problems (see (2.26) and (2.27) below). In essence, he justified that the drag force acting on each solid obstacle, spherical or not, coincides at the leading order with the one associated to an isolated particle in $\mathbb{R}^d$. However, Allaire's method relies on weak convergence techniques which do not yield quantitative rates of convergence nor provide the next terms of the asymptotic expansion.

Further progress in this direction was achieved by Musolino and Mityushev [30] who studied the asymptotic behavior of the longitudinal permeability of a periodic array of cylinders. This is a simplified setting where the velocity has only a single component in the transverse direction which solves a two-dimensional perforated Poisson equation. The longitudinal permeability corresponds then to a single scalar coefficient

$$ \chi^* := \int_{P\setminus(\eta T)} \chi \, dx, \qquad (1.9) $$

where $\chi$ is the solution to the perforated Dirichlet problem

$$ \left\{ \begin{array}{ll} -\Delta\chi = 1 & \text{in } P\setminus(\eta T), \\ \chi = 0 & \text{on } \partial(\eta T), \\ \chi \text{ is } P\text{-periodic.} \end{array} \right. \qquad (1.10) $$

In the spirit of numerous works applying integral equation techniques to the asymptotic analysis of a various range of physical systems [36, 4, 11, 5, 3], they derived a layer potential representation of the solution $\chi$ to prove the existence of an asymptotic expansion as $\eta \to 0$ for $\chi$ and the coefficient $\chi^*$. They identified that the coefficient $\chi^*$ is the sum of a logarithmic term and a power series in the square of $\eta$; however without providing a very explicit procedure for computing all the terms of this expansion. Independently, Jing [26] obtained with a comparable approach the rate of convergence of the leading order asymptotic expansion of the coefficient $\chi^*$ (the analogue of (1.8) for the Poisson system (1.10)) in dimension $d \ge 3$.

Then, in our recent paper [18], we derived full asymptotic expansions for the solution $\chi$ to the scalar problem (1.10) in both cases $d = 2$ and $d \ge 3$. However, we didn't use these asymptotic series to deduct full asymptotic expansions of the longitudinal permeability $\chi^*$.

Additionally to these works, numerical methods based on layer potential techniques have been developed in order to obtain fast numerical solutions to periodic Stokes problems, see [22] and the references therein. Finally, let us note that a related problem is the computation of the drag force for Stokes flows in random arrays of solid particles, which has been considered in the case of spherical obstacles in [9, 25].

Motivated by the above considerations, the goal of this paper is to derive full asymptotic expansions of the permeability matrix $\chi^*$ defined in (1.5) in arbitrary dimension $d \ge 2$ and for arbitrarily shaped obstacles $T$. The outcome of this computation, which is the main result of this paper, is stated in Proposition 4.2. We note that the full asymptotic expansion of the scalar coefficient $\chi^*$ could be obtained identically with the techniques presented in this work and would result in analogous formulas.

The paper outlines as follows. Section 2 introduces the necessary preliminary background material: we recall basic results on hydrodynamics layer potential theory and exterior Stokes problems. We also propose an existence and uniqueness theory for Stokes problems with some special periodic boundary conditions, and
Although there had already been attempts by Burgers [8], Brinkman [7] and Debye and Bueche [12] to propose physical models describing the fluid properties in porous media, Hasimoto [24] was the first to derive a mathematically correct asymptotic expansion for the drag force induced by a periodic array of spheres, in two and three dimensions. His derivation was based on the computation of an approximation of the fluid velocity by a linear combination of the periodic Stokes fundamental solution and its derivatives. With the present notation conventions, Hasimoto found

$$ \chi^* = \left( \frac{1}{6\pi\eta} - \frac{c}{6\pi} + \frac{2}{9}\eta^2 - \frac{1}{6\pi} \left( \frac{16\pi^2}{45} + 630b^2 \right) \eta^6 + O(\eta^8) \right) \text{I}, \quad \text{if } d=3, \quad (1.6) $$

where $c \approx 2.8373$ and $b \approx 0.19457$ and I is the identity matrix, and

$$ \chi^* = \left( -\frac{1}{4\pi} \log(\eta) - \frac{\tilde{c}}{4\pi} + \frac{1}{4}\eta^2 + O(\eta^4) \right) \text{I}, \quad \text{if } d=2, \quad (1.7) $$

with $\tilde{c} \approx 1.3105$. Hasimoto noted that (1.6) can be understood as a higher order correction of the classical formula $F = 6\pi\eta U$ describing the force that acts on a single isolated sphere of radius $\eta$. This reflects the fact that in the dilute regime $\eta \to 0$, the obstacles behave at first order as isolated. Sangani and Acrivos subsequently refined Hasimoto's procedure to obtain higher order terms [38], however their method remained applicable only to spherical obstacles.

Periodic arrays with obstacles of arbitrary shape were first considered in the work of Allaire [1] who identified the first term of the asymptotic expansion of the effective permeability matrix in arbitrary dimension $d \ge 2$. He obtained

$$ \chi^* \sim \begin{cases} \eta^{2-d} M^{-1} & \text{if } d \ge 3, \\ -\frac{1}{4\pi} \log(\eta) \text{I} & \text{if } d = 2, \end{cases} \quad (1.8) $$

where, in the case $d \ge 3$, $M^{-1}$ is the Brinkmann inverse permeability which is defined from the solutions to exterior Stokes problems (see (2.26) and (2.27) below). In essence, he justified that the drag force acting on each solid obstacle, spherical or not, coincides at the leading order with the one associated to an isolated particle in $\mathbb{R}^d$. However, Allaire's method relies on weak convergence techniques which do not yield quantitative rates of convergence nor provide the next terms of the asymptotic expansion.

Further progress in this direction was achieved by Musolino and Mityushev [34] who studied the asymptotic behavior of the longitudinal permeability of a periodic array of cylinders. This is a simplified setting where the velocity has only a single component in the transverse direction which solves a two-dimensional perforated Poisson equation. The longitudinal permeability corresponds then to a single scalar coefficient

$$ \chi^* := \int_{P \setminus (\eta T)} \chi \, dx, \quad (1.9) $$

where $\chi$ is the solution to the perforated Dirichlet problem

$$ \left\{ \begin{array}{l} -\Delta\chi = 1 \text{ in } P \setminus (\eta T), \\ \chi = 0 \text{ on } \partial(\eta T), \\ \chi \text{ is } P\text{-periodic.} \end{array} \right. \quad (1.10) $$

In the spirit of numerous works applying integral equation techniques to the asymptotic analysis of a various range of physical systems [40, 4, 11, 5, 3], they derived a layer potential representation of the solution $\chi$ to prove the existence of an asymptotic expansion as $\eta \to 0$ for $\chi$ and the coefficient $\chi^*$. They identified that the coefficient $\chi^*$ is the sum of a logarithmic term and a power series in the square of $\eta$. Later, Mityushev and Adler [33] derived a procedure for computing all the terms of this expansion, however restricted to the case of cylindrical obstacles. Independently, Jing [27] obtained with a comparable approach the rate of convergence of the leading order asymptotic expansion of the coefficient $\chi^*$ (the analogue of (1.8) for the Poisson system (1.10)), for arbitrary obstacles and in dimension $d \ge 3$. These results have been recently extended to the Stokes permeability matrix (1.5) in [28] for $d \ge 2$.

Independently, we derived full asymptotic expansions for the solution $\chi$ to the scalar problem (1.10) in both cases $d=2$ and $d \ge 3$ [19]. However, we didn't use these asymptotic series to deduct full asymptotic expansions of the longitudinal permeability $\chi^*$.

Additionally to these works, let us mention that numerical methods based on layer potential techniques have been developed in order to obtain fast numerical solutions to periodic Stokes problems [23]. A related problem is also the computation of the drag force for Stokes flows in random arrays of solid particles, which has been considered in the case of spherical obstacles in [9, 26]. Finally, a lot of attention has been recently given to the rigorous derivation of Einstein's effective viscosity formula for random arrays of moving particles with stochastic frameworks, see [14] and the references therein.

Motivated by the above considerations, the goal of this paper is to derive full asymptotic expansions of the permeability matrix $\chi^*$ defined in (1.5) in arbitrary dimension $d \ge 2$ and for arbitrarily shaped obstacles $T$.
we derive useful identities about some tensors $(\mathfrak{G}^k)_{k \ge 0}$ and $(\mathfrak{R}^k)_{k \ge 0}$ associated to the periodic fundamental solution of the Stokes equations that play an important role in the derivation of our main result.

Section 3 unfolds the first step of our analysis where we derive full asymptotic expansions of the solutions $(\mathcal{X}_j, \alpha_j)$ to the periodic cell problem (1.3) by applying a methodology similar to that of our previous work [18]. Similarly to the results obtained for the scalar problem (1.10), we find that these expansions can be explicitly computed by solving a cascade of interior and exterior problems which are derived in Proposition 3.4. Remarkably, the obtained asymptotic series are summable, which yields error estimates with constants independent of the order chosen for truncating the expansion (Proposition 3.5).

The second step, performed in Section 4, is to integrate the expansion of $\mathcal{X}_j$ over $P \setminus (\eta T)$ to deduce, in Proposition 4.2, full asymptotic asymptotic expansions of the permeability matrix $\mathcal{X}^*$. An important technical aspect of this procedure is the evaluation of some finite part integrals as $\eta \to 0^+$ over the domain with growing boundary $\eta^{-1}P \setminus T$. Using explicit solutions available in the case of spherical obstacles, we then verify that our expansions agree at the leading order with Hasimoto's formula (1.6) and (1.7) (Proposition 4.3).

Finally, Appendix A presents a well-posedness theory for perforated periodic Stokes problems and establishes related a priori estimates. This theory provides a second proof of the error estimates holding for the asymptotic expansion of $(\mathcal{X}_j, \alpha_j)$, which is more elementary in the sense that it does not rely on potential theory, but which yields a slightly weaker bound (featuring a constant that is not uniform with respect to the order of the approximation). Appendices B and C then present a few technical results regarding the evaluation of some moment integrals on the hypercube $P$ and the hypersphere $B(0, 1)$ respectively.

## 2. PRELIMINARIES

Throughout the paper, the bold face notation is used for vectors of $\mathbb{R}^d$ (for instance $(e_i)_{1 \le i \le d}$ are the constant vector fields aligned with the canonical orthonormal basis of $\mathbb{R}^d$) but the standard roman notation for those vectors thought as points, e.g. $x \in P \setminus (\eta T)$. For any velocity and pressure field $v, p$, we denote by $\sigma(v, p)$ the associated Newton stress tensor:

$$ \sigma(v, p) := \nabla v - pI, \qquad (2.1) $$

where $I$ is the $d \times d$ identity matrix. Thus $\sigma(v, p)n$ is the co-normal derivative associated to the differential expression $-\Delta v + \nabla p$. Note that this $\sigma$ must not be confused with the notation $d\sigma$ also used in this paper to denote to the measure arising in surface integrals.

For any $k \in \mathbb{N}$ and matrix, vector, or scalar field $A$, we denote by $\nabla^k A$ the $k$-th order matrix-valued, vector-valued, or scalar valued tensor of partial derivatives:

$$ \nabla^k A = (\partial_{i_1 \dots i_k}^k A)_{1 \le i_1 \dots i_k \le d}. $$

For a $d$-dimensional coordinate vector $y \in \mathbb{R}^d$, we denote by $y^k$ the $k$-th order scalar tensor with components

$$ y_{i_1 \dots i_k}^k = y_{i_1} y_{i_2} \dots y_{i_k}, \quad 1 \le i_1 \dots i_k \le d. $$

The “natural” contraction between two $k$-th order tensors $A^k$ and $B^k$ is denoted by

$$ A^k \cdot B^k = \sum_{1 \le i_1 \dots i_k \le d} A_{i_1 \dots i_k}^k B_{i_1 \dots i_k}^k, $$

where $A_{i_1 \dots i_k}^k B_{i_1 \dots i_k}^k$ is the natural product operation between $A_{i_1 \dots i_k}^k$ and $B_{i_1 \dots i_k}^k$: for instance a matrix-vector product in case of the contraction between a matrix valued tensor and a vector valued tensor, an inner product in the case of two vector valued tensors, a scalar multiplication if one of the two tensors is scalar valued. For example, $\nabla^k G(y) \cdot t^k \phi(t)$ in (2.16) denotes the tensor contraction

$$ \nabla^k G(y) \cdot t^k \phi(t) = \sum_{1 \le i_1 \dots i_k \le d} t_{i_1} \dots t_{i_k} \partial_{i_1 \dots i_k}^k G(y) \phi(t), $$

where a matrix-vector product occurs between the matrix $\partial_{i_1 \dots i_k}^k G(y)$ and the vector $\phi(t)$.

Given a Lipschitz domain $T$, we denote by $n$ or $\mathbf{n}$ the normal pointing outward $T$ and by $[u]$ the jump of a function $u$ across $\partial T$:

$$ [u] := u|_{+} - u|_{-} \text{ with } u|_{+}(y) := \lim_{t \to 0} u(y + tn), \quad u|_{-}(y) := \lim_{t \to 0} u(y - tn), \quad y \in \partial T. $$

Finally, the Kronecker notation $\delta_{d=2}$ is used to mean a constant which is equal to 1 if $d = 2$ and zero otherwise.
The outcome of this computation, which is the main result of this paper, is stated in Proposition 4.2. We note that the full asymptotic expansion of the scalar coefficient $\chi^*$ could be obtained identically with the techniques presented in this work and would result in analogous formulas. This work presents some similarities with [28], where layer potential techniques are also used for the asymptotic analysis of the permeability matrix $\mathcal{X}^*$. The main differences lies in that the authors of [28] rely on a double layer potential representation of the correc- tors $(\mathcal{X}_i)_{1\le i\le d}$, while this work uses a characterization based on the single layer potential (Proposition 3.2). Moreover, the focus of [28] is to derive quantitative periodic homogenization results. Hence, layer potential representations were used to bound the difference between $\mathcal{X}^*$ and its leading order asymptotic (1.8). The focus of this work is to show that it is actually possible to identify all higher order terms through well-defined inductive computational procedures.

The paper outlines as follows. Section 2 introduces the necessary preliminary background material: we recall basic results on hydrodynamics layer potential theory and exterior Stokes problems. We also propose an existence and uniqueness theory for Stokes problems with specific sets of periodic boundary conditions, and we derive useful identities about some tensors $(\mathfrak{B}^k)_{k\ge 0}$ and $(\mathfrak{R}^k)_{k\ge 0}$ associated to the periodic fundamental solution that play an important role in our derivation.

Section 3 unfolds the first step of our analysis, where we derive full asymptotic expansions of the solutions $(\mathcal{X}_j, \alpha_j)$ to the periodic cell problem (1.3) by applying a methodology similar to that of our previous work [19]. Similarly to the results obtained for the scalar problem (1.10), we find that these expansions can be explicitly computed by solving a cascade of interior and exterior problems which are derived in Proposition 3.4. Remark- ably, the obtained asymptotic series are summable, which yields error estimates with constants independent of the order chosen for truncating the expansion (Proposition 3.5).

The second step, performed in Section 4, is to integrate the expansion of $\mathcal{X}_j$ over $P\setminus(\eta T)$ to deduce, in Proposition 4.2, full asymptotic asymptotic expansions of the permeability matrix $\mathcal{X}^*$. An important technical aspect of this procedure is the evaluation of some finite part integrals as $\eta \to 0^+$ over the domain with growing boundary $\eta^{-1}P\setminus T$. Using explicit solutions available in the case of spherical obstacles, we verify that our expansions agree at the leading order with Hasimoto’s formula (1.6) and (1.7) (Proposition 4.3).

This paper contains three additional appendices. Appendix A presents a well-posedness theory for perforated periodic Stokes problems and establishes related a priori estimates. This theory provides a second proof of the error estimates holding for the asymptotic expansion of $(\mathcal{X}_j, \alpha_j)$, which does not rely on potential theory but yields a slightly weaker bound. Appendices B and C present technical results about the evaluation of some moment integrals on the hypercube $P$ and the hypersphere $B(0, 1)$ respectively.

## 2. PRELIMINARIES

Throughout the paper, the bold face notation is used for vectors of $\mathbb{R}^d$ (for instance $(e_i)_{1\le i\le d}$ are the constant vector fields aligned with the canonical orthonormal basis of $\mathbb{R}^d$), but the standard roman notation for those vectors thought as points, e.g. $x \in P\setminus(\eta T)$. For any velocity and pressure field $v$, $p$, we denote by $\sigma(v, p)$ the associated Newton stress tensor:

$$ \sigma(v, p) := \nabla v - pI, \qquad (2.1) $$

where $I$ is the $d \times d$ identity matrix. Thus $\sigma(v, p)n$ is the co-normal derivative associated to the differential expression $-\Delta v + \nabla p$. Note that this $\sigma$ must not be confused with the notation $d\sigma$ also used in this paper to denote to the measure arising in surface integrals.

For any $k \in \mathbb{N}$ and matrix, vector, or scalar field $A$, we denote by $\nabla^k A$ the $k$-th order matrix-valued, vector-valued, or scalar valued tensor of partial derivatives:

$$ \nabla^k A = (\partial_{i_1 \dots i_k}^k A)_{1 \le i_1 \dots i_k \le d}. $$

For a $d$-dimensional coordinate vector $y \in \mathbb{R}^d$, we denote by $y^k$ the $k$-th order scalar tensor with components

$$ y_{i_1 \dots i_k}^k = y_{i_1} y_{i_2} \dots y_{i_k}, \quad 1 \le i_1 \dots i_k \le d. $$

The “natural” contraction between two $k$-th order tensors $A^k$ and $B^k$ is denoted by

$$ A^k \cdot B^k = \sum_{1 \le i_1 \dots i_k \le d} A_{i_1 \dots i_k}^k B_{i_1 \dots i_k}^k, $$

where $A_{i_1 \dots i_k}^k B_{i_1 \dots i_k}^k$ is the natural product operation between $A_{i_1 \dots i_k}^k$ and $B_{i_1 \dots i_k}^k$: for instance a matrix-vector product in case of the contraction between a matrix valued tensor and a vector valued tensor, an inner product in the case of two vector valued tensors, a scalar multiplication if one of the two tensors is scalar valued. For example, $\nabla^k G(y) \cdot t^k \phi(t)$ in (2.16) denotes the tensor contraction

$$ \nabla^k G(y) \cdot t^k \phi(t) = \sum_{1 \le i_1 \dots i_k \le d} t_{i_1} \dots t_{i_k} \partial_{i_1 \dots i_k}^k G(y) \phi(t), $$

where a matrix-vector product occurs between the matrix $\partial_{i_1 \dots i_k}^k G(y)$ and the vector $\phi(t)$.
## 2.1. Hydrodynamics layer potentials

Hydrodynamics layer potentials are defined from the knowledge of the fundamental solution of the Stokes system [13, 32, 28, 20], whose expression can be found from the fundamental solution of the Laplace operator denoted by $\Gamma$:

$$ \Gamma(x) = \begin{cases} \frac{1}{2\pi} \log|x| & \text{if } d=2, \\ -\frac{1}{(d-2)|\partial B(0,1)||x|^{d-2}} & \text{if } d \ge 3, \end{cases} \qquad (2.2) $$

where $|\partial B(0, 1)|$ is the area of the $d$-dimensional unit sphere.

**Proposition 2.1.** For any $1 \le i \le d$, a fundamental solution $(G_i, \psi_i)$ satisfying

$$ \left\{ \begin{array}{l} \Delta G_i - \nabla \psi_i = \delta_0 e_i \text{ in } \mathbb{R}^d, \\ \text{div}(G_i) = 0 \text{ in } \mathbb{R}^d, \end{array} \right. \qquad (2.3) $$

in the sense of distributions is given by

$$ G_i(x) = \frac{1}{2}\Gamma(x)e_i - \frac{1}{2}x_i\nabla\Gamma(x), \quad \psi_i(x) = -\partial_i\Gamma(x), \quad x \in \mathbb{R}^d \setminus \{0\}, \qquad (2.4) $$

that is more explicitly

$$ \left\{ \begin{array}{l} G_i(x) = \frac{1}{4\pi} \left( \log|x|e_i - \frac{x_i x}{|x|^2} \right), \\ \psi_i(x) = -\frac{1}{2\pi} \frac{x_i}{|x|^2}, \end{array} \right. \quad \text{if } d=2, $$

and

$$ \left\{ \begin{array}{l} G_i(x) = -\frac{1}{2(d-2)|\partial B(0,1)|} \left( \frac{e_i}{|x|^{d-2}} + (d-2) \frac{x_i x}{|x|^d} \right), \\ \psi_i(x) = -\frac{1}{|\partial B(0,1)|} \frac{x_i}{|x|^d}, \end{array} \right. \quad \text{if } d \ge 3. $$

*Proof.* Following [3, 13], computing the divergence of the first equation of (2.3) yields $\Delta\psi_i = -\partial_i\delta_0$, which leads to $\psi_i = -\partial_i\Gamma$. Substituting back into (2.3), we obtain

$$ \Delta G_i = \delta_0 e_i - \sum_{j=1}^{d} \partial_{ij} \Gamma e_j. \qquad (2.5) $$

On the other hand, a simple computation yields

$$ \Delta(x_i \partial_j \Gamma) = 2 \partial_{ij} \Gamma + x_i \partial_j \delta_0 = 2 \partial_{ij} \Gamma - \delta_{ij} \delta_0 = 2 \partial_{ij} \Gamma - \delta_{ij} \Delta \Gamma, $$

whence

$$ \Delta \left( \frac{1}{2} x_i \partial_j \Gamma + \frac{1}{2} \delta_{ij} \Gamma \right) = \partial_{ij} \Gamma. $$

Substituting into (2.5), we find that a solution to (2.3) is given by

$$ G_i(x) = \Gamma(x)e_i - \frac{1}{2} \sum_{j=1}^{d} (x_i \partial_j \Gamma(x) + \delta_{ij} \Gamma(x))e_j = \frac{1}{2} \Gamma(x)e_i - \frac{1}{2} x_i \nabla \Gamma(x), $$

which is (2.4).

In all what follows, we denote by $G(x)$ the fundamental matrix of the Stokes system, also called “Oseen matrix”, whose columns are the vector fields $G_i(x)$. Similarly, we denote by $\psi(x)$ the fundamental vector whose components are the scalar pressure fields $\psi_i(x)$:

$$ G(x) = [G_1(x) \quad \dots \quad G_d(x)], \quad \psi(x) = (\psi_i(x))_{1 \le i \le d}. $$

The next lemma states important identities satisfied by the fundamental solution $(G_i, \psi_i)$.

**Lemma 2.1.** (i) For any smooth open bounded domain $T \in \mathbb{R}^d$ with $0 \in T$:

$$ \int_{\partial T} \sigma(G_i, \psi_i) n \, d\sigma = e_i \quad \text{and} \quad \int_{\partial T} \sigma(\nabla^k G_i, \nabla^k \psi_i) n \, d\sigma(x) = 0 \quad \text{for any } k \ge 1. $$

(ii) For any orthogonal symmetry $S$ (satisfying $S = S^T = S^{-1}$), it holds

$$ S G_i \circ S = \sum_{j=1}^{d} S_{ij} G_j, \quad \psi_i \circ S = \sum_{j=1}^{d} S_{ij} \psi_j. \qquad (2.6) $$
Given a Lipschitz domain $T$, we denote by $n$ or $\mathbf{n}$ the normal pointing outward $T$ and by $\llbracket u \rrbracket$ the jump of a function $u$ across $\partial T$:

$$ \llbracket u \rrbracket := u|_{+} - u|_{-} \text{ with } u|_{+}(y) := \lim_{t \to 0} u(y + tn), \quad u|_{-}(y) := \lim_{t \to 0} u(y - tn), \quad y \in \partial T. $$

Finally, the Kronecker notation $\delta_{d=2}$ is used to mean a constant which is equal to 1 if $d=2$ and zero otherwise.

## 2.1. Hydrodynamics layer potentials

Hydrodynamics layer potentials are defined from the knowledge of the fundamental solution of the Stokes system [13, 36, 30, 21], whose expression can be found from the fundamental solution of the Laplace operator denoted by $\Gamma$:

$$ \Gamma(x) = \begin{cases} \frac{1}{2\pi} \log|x| & \text{if } d=2, \\ -\frac{1}{(d-2)|\partial B(0,1)||x|^{d-2}} & \text{if } d \ge 3, \end{cases} \qquad (2.2) $$

where $|\partial B(0,1)|$ is the area of the $d$-dimensional unit sphere.

**Proposition 2.1.** For any $1 \le i \le d$, a fundamental solution $(G_i, \psi_i)$ satisfying

$$ \left\{ \begin{array}{l} \Delta G_i - \nabla \psi_i = \delta_0 e_i \text{ in } \mathbb{R}^d, \\ \text{div}(G_i) = 0 \text{ in } \mathbb{R}^d, \end{array} \right. \qquad (2.3) $$

in the sense of distributions is given by

$$ G_i(x) = \frac{1}{2}\Gamma(x)e_i - \frac{1}{2}x_i\nabla\Gamma(x), \quad \psi_i(x) = -\partial_i\Gamma(x), \quad x \in \mathbb{R}^d\setminus\{0\}, \qquad (2.4) $$

that is more explicitly

$$ \left\{ \begin{array}{l} G_i(x) = \frac{1}{4\pi} \left( \log|x|e_i - \frac{x_i x}{|x|^2} \right), \\ \psi_i(x) = -\frac{1}{2\pi} \frac{x_i}{|x|^2}, \end{array} \right. \quad \text{if } d=2, $$

and

$$ \left\{ \begin{array}{l} G_i(x) = -\frac{1}{2(d-2)|\partial B(0,1)|} \left( \frac{e_i}{|x|^{d-2}} + (d-2) \frac{x_i x}{|x|^d} \right), \\ \psi_i(x) = -\frac{1}{|\partial B(0,1)|} \frac{x_i}{|x|^d}, \end{array} \right. \quad \text{if } d \ge 3. $$

*Proof.* Following [3, 13], computing the divergence of the first equation of (2.3) yields $\Delta\psi_i = -\partial_i\delta_0$, which leads to $\psi_i = -\partial_i\Gamma$. Substituting back into (2.3), we obtain

$$ \Delta G_i = \delta_0 e_i - \sum_{j=1}^{d} \partial_{ij} \Gamma e_j. \qquad (2.5) $$

On the other hand, a simple computation yields

$$ \Delta(x_i \partial_j \Gamma) = 2 \partial_{ij} \Gamma + x_i \partial_j \delta_0 = 2 \partial_{ij} \Gamma - \delta_{ij} \delta_0 = 2 \partial_{ij} \Gamma - \delta_{ij} \Delta \Gamma, $$

whence

$$ \Delta \left( \frac{1}{2} x_i \partial_j \Gamma + \frac{1}{2} \delta_{ij} \Gamma \right) = \partial_{ij} \Gamma. $$

Substituting into (2.5), we find that a solution to (2.3) is given by

$$ G_i(x) = \Gamma(x)e_i - \frac{1}{2} \sum_{j=1}^{d} (x_i \partial_j \Gamma(x) + \delta_{ij} \Gamma(x))e_j = \frac{1}{2} \Gamma(x)e_i - \frac{1}{2} x_i \nabla \Gamma(x), $$

which is (2.4).

In all what follows, we denote by $G(x)$ the fundamental matrix of the Stokes system, also called “Oseen matrix”, whose columns are the vector fields $G_i(x)$. Similarly, we denote by $\psi(x)$ the fundamental vector whose components are the scalar pressure fields $\psi_i(x)$:

$$ G(x) = [G_1(x) \quad \dots \quad G_d(x)], \quad \psi(x) = (\psi_i(x))_{1 \le i \le d}. $$

The next lemma states important identities satisfied by the fundamental solution $(G_i, \psi_i)$.

**Lemma 2.1.** (i) For any smooth open bounded domain $T \in \mathbb{R}^d$ with $0 \in T$:

$$ \int_{\partial T} \sigma(G_i, \psi_i) \mathbf{n} \, d\sigma = e_i \quad \text{and} \quad \int_{\partial T} \sigma(\nabla^k G_i, \nabla^k \psi_i) \mathbf{n} \, d\sigma(x) = 0 \quad \text{for any } k \ge 1. $$
Proof. (i) Let $\epsilon > 0$ be a small positive real number. Integrating the identity $\Delta(\nabla^k G_i) - \nabla(\nabla^k \psi_i) = 0$ inside $T \setminus B(0, \epsilon)$ and using the divergence theorem, we obtain

$$0 = \int_{T \setminus B(0, \epsilon)} (\Delta(\nabla^k G_i) - \nabla(\nabla^k \psi_i)) dx = \int_{\partial T} \sigma(\nabla^k G_i, \nabla^k \psi_i) n \, d\sigma - \int_{\partial B(0, \epsilon)} \sigma(\nabla^k G_i, \nabla^k \psi_i) n \, d\sigma.$$

Performing a change of variable and using the fact that $\sigma(\nabla^k G_i, \nabla^k \psi_i)n = \sum_{j=1}^d \partial_j(\nabla^k G_i)n_j - \nabla^k \psi_i n$ is homogeneous of degree $1-d-k$, we obtain

$$\int_{\partial T} \sigma(\nabla^k G_i, \nabla^k \psi_i) n \, d\sigma = \epsilon^{d-1} \epsilon^{1-d-k} \int_{\partial B(0,1)} \sigma(\nabla^k G_i, \nabla^k \psi_i) n \, d\sigma = \epsilon^{-k} \int_{\partial B(0,1)} \sigma(\nabla^k G_i, \nabla^k \psi_i) n \, d\sigma. \quad (2.7)$$

For $k=0$, an explicit computation using (2.4) yields

$$\int_{\partial T} \sigma(G_i, \psi_i) n \, d\sigma = \int_{\partial B(0,1)} \left( \frac{1}{2} \frac{\partial \Gamma}{\partial n} e_i - \frac{1}{2} (n \cdot e_i) \nabla \Gamma - \frac{1}{2} x_i \nabla^2 \Gamma n + \partial_i \Gamma n \right) \, d\sigma.$$

On the sphere, $n = x$, $\nabla\Gamma(x) = x/|\partial B(0,1)|$ and $\nabla^2\Gamma(x) = (I - dxx^T)/|\partial B(0,1)|$. Using $\int_{\partial B(0,1)} xx^T \, d\sigma(x) = \frac{|\partial B(0,1)|}{d}I$, we obtain

$$\int_{\partial T} \sigma(G_i, \psi_i) n \, d\sigma = \frac{1}{2} e_i - \frac{1}{2d} e_i - \frac{1}{2d} e_i + \frac{1}{2} e_i + \frac{1}{d} e_i = e_i.$$

For $k \ge 1$, we see that the left-hand side of (2.7) is a finite quantity as $\epsilon \to 0$. This is possible only if $\int_{\partial B(0,1)} \sigma(\nabla^k G_i, \nabla^k \psi_i) n \, d\sigma = 0$, which implies the result.

(ii) The spherical symmetry property can be verified directly on the expression (2.4), using $\Gamma(Sx) = \Gamma(x)$. $\square$

For any Lipschitz bounded open set $T \subset \mathbb{R}^d$, we denote by $S_T$ the (Stokes) single layer potential on $\partial T$:

$$\forall x \in \mathbb{R}^d, \forall \phi \in H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d), \quad S_T[\phi](x) := \int_{\partial T} G(x-y)\phi(y)\sigma(y). \quad (2.8)$$

We also define $\mathcal{P}_T[\phi]$ the potential associated to the pressure:

$$\forall x \in \mathbb{R}^d, \forall \phi \in H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d), \quad \mathcal{P}_T[\phi](x) := \text{p.v.} \int_{\partial T} \psi(x-y) \cdot \phi(y) \sigma(y), \quad (2.9)$$

where p.v. indicates a Cauchy principal value integral when $x \in \partial T$. The main properties of these potentials are recalled in the next proposition [16, 32, 37].

**Proposition 2.2.** (i) For any $\phi \in H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d)$, $(S_T[\phi], \mathcal{P}_T[\phi])$ solves the Stokes system

$$\left\{ \begin{array}{l} \Delta S_T[\phi] - \nabla \mathcal{P}_T[\phi] = 0 \quad \text{in } \mathbb{R}^d \setminus \partial T, \\ \text{div}(S_T[\phi]) = 0 \quad \text{in } \mathbb{R}^d \setminus \partial T. \end{array} \right. \quad (2.10)$$

$$(ii) S_T[n] = 0 \text{ and } \mathcal{P}_T[n] = \begin{cases} -1 & \text{in } T, \\ -\frac{1}{2} & \text{on } \partial T, \\ 0 & \text{in } \mathbb{R}^d \setminus \bar{T}. \end{cases}$$

(iii) We have the following jump identities:

$$[S_T[\phi]] = 0, \quad (2.11)$$

$$[\nabla S_T[\phi] n] = \phi - (\phi \cdot n)n, \quad (2.12)$$

$$[\mathcal{P}_T[\phi]] = -\phi \cdot n, \quad (2.13)$$

In particular,

$$\phi = [\sigma(S_T[\phi], \mathcal{P}_T[\phi])n]. \quad (2.14)$$

(iv) The following continuity bound hold true:

$$\|\nabla S_T[\phi]\|_{L^2(\mathbb{R}^d, \mathbb{R}^{d \times d})} + \|\mathcal{P}_T[\phi]\|_{L^2(\mathbb{R}^d \setminus T)} + \|\mathcal{P}_T[\phi]\|_{L^2(T)/\mathbb{R}} \le C \|\phi\|_{H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d)}. \quad (2.15)$$

(v) $S_T[\phi]$ and $\mathcal{P}_T[\phi]$ have the following asymptotic behavior at infinity:

$$S_T[\phi](y) = \sum_{k=0}^{+\infty} \frac{(-1)^k}{k!} \nabla^k G(y) \cdot \int_{\partial T} t^k \phi(t) \, d\sigma(t) \quad \text{as } |y| \to +\infty, \quad (2.16)$$

$$\mathcal{P}_T[\phi](y) = \sum_{k=0}^{+\infty} \frac{(-1)^k}{k!} \nabla^k \psi(y) \cdot \int_{\partial T} t^k \phi(t) \, d\sigma(t) \quad \text{as } |y| \to +\infty, \quad (2.17)$$

where these series converge for $|y|$ sufficiently large.
(ii) For any orthogonal symmetry S (satisfying $S = S^T = S^{-1}$), it holds

$$ S\mathbf{G}_i \circ S = \sum_{j=1}^{d} S_{ij}\mathbf{G}_j, \quad \psi_i \circ S = \sum_{j=1}^{d} S_{ij}\psi_j. \qquad (2.6) $$

Proof. (i) Let $\epsilon > 0$ be a small positive real number. Integrating the identity $\Delta(\nabla^k \mathbf{G}_i) - \nabla(\nabla^k \psi_i) = 0$ inside $T \setminus B(0, \epsilon)$ and using the divergence theorem, we obtain

$$ 0 = \int_{T \setminus B(0, \epsilon)} (\Delta(\nabla^k \mathbf{G}_i) - \nabla(\nabla^k \psi_i)) \, dx = \int_{\partial T} \sigma(\nabla^k \mathbf{G}_i, \nabla^k \psi_i) \mathbf{n} \, d\sigma - \int_{\partial B(0, \epsilon)} \sigma(\nabla^k \mathbf{G}_i, \nabla^k \psi_i) \mathbf{n} \, d\sigma. $$

Performing a change of variable and using the fact that $\sigma(\nabla^k \mathbf{G}_i, \nabla^k \psi_i)\mathbf{n} = \sum_{j=1}^d \partial_j(\nabla^k \mathbf{G}_i)n_j - \nabla^k \psi_i \mathbf{n}$ is homogeneous of degree $1-d-k$, we obtain

$$ \int_{\partial T} \sigma(\nabla^k \mathbf{G}_i, \nabla^k \psi_i) \mathbf{n} \, d\sigma = \epsilon^{d-1} \epsilon^{1-d-k} \int_{\partial B(0,1)} \sigma(\nabla^k \mathbf{G}_i, \nabla^k \psi_i) \mathbf{n} \, d\sigma = \epsilon^{-k} \int_{\partial B(0,1)} \sigma(\nabla^k \mathbf{G}_i, \nabla^k \psi_i) \mathbf{n} \, d\sigma. \quad (2.7) $$

For $k=0$, an explicit computation using (2.4) yields

$$ \int_{\partial T} \sigma(\mathbf{G}_i, \psi_i) \mathbf{n} \, d\sigma = \int_{\partial B(0,1)} \left( \frac{1}{2} \frac{\partial \Gamma}{\partial n} e_i - \frac{1}{2} (\mathbf{n} \cdot \mathbf{e}_i) \nabla \Gamma - \frac{1}{2} x_i \nabla^2 \Gamma \mathbf{n} + \partial_i \Gamma \mathbf{n} \right) \, d\sigma. $$

On the sphere, $\mathbf{n} = \mathbf{x}$, $\nabla\Gamma(\mathbf{x}) = \mathbf{x}/|\partial B(0,1)|$ and $\nabla^2\Gamma(\mathbf{x}) = (I - dxx^T)/|\partial B(0,1)|$. Using $\int_{\partial B(0,1)} xx^T \, d\sigma(\mathbf{x}) = \frac{|\partial B(0,1)|}{d}I$, we obtain

$$ \int_{\partial T} \sigma(\mathbf{G}_i, \psi_i) \mathbf{n} \, d\sigma = \frac{1}{2} e_i - \frac{1}{2d} e_i - \frac{1}{2d} e_i + \frac{1}{2} e_i + \frac{1}{d} e_i = e_i. $$

For $k \ge 1$, we see that the left-hand side of (2.7) is a finite quantity as $\epsilon \to 0$. This is possible only if $\int_{\partial B(0,1)} \sigma(\nabla^k \mathbf{G}_i, \nabla^k \psi_i) \mathbf{n} \, d\sigma = 0$, which implies the result.

(ii) The spherical symmetry property can be verified directly on the expression (2.4), using $\Gamma(Sx) = \Gamma(x)$. $\square$

For any Lipschitz bounded open set $T \subset \mathbb{R}^d$, we denote by $S_T$ the (Stokes) single layer potential on $\partial T$:

$$ \forall x \in \mathbb{R}^d, \forall \phi \in H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d), \quad S_T[\phi](x) := \int_{\partial T} G(x-y)\phi(y)\sigma(y). \qquad (2.8) $$

We also define $\mathcal{P}_T[\phi]$ the potential associated to the pressure:

$$ \forall x \in \mathbb{R}^d, \forall \phi \in H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d), \quad \mathcal{P}_T[\phi](x) := \text{p.v.} \int_{\partial T} \psi(x-y) \cdot \phi(y) \sigma(y), \qquad (2.9) $$

where p.v. indicates a Cauchy principal value integral when $x \in \partial T$. The main properties of these potentials are recalled in the next proposition [17, 36, 41].

**Proposition 2.2.** (i) For any $\phi \in H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d)$, $(S_T[\phi], \mathcal{P}_T[\phi])$ solves the Stokes system

$$ \left\{ \begin{array}{l} \Delta S_T[\phi] - \nabla \mathcal{P}_T[\phi] = 0 \quad \text{in } \mathbb{R}^d \setminus \partial T, \\ \text{div}(S_T[\phi]) = 0 \quad \text{in } \mathbb{R}^d \setminus \partial T. \end{array} \right. \qquad (2.10) $$

$$ \text{(ii) } S_T[\mathbf{n}] = 0 \text{ and } \mathcal{P}_T[\mathbf{n}] = \begin{cases} -1 & \text{in } T, \\ -\frac{1}{2} & \text{on } \partial T, \\ 0 & \text{in } \mathbb{R}^d \setminus \bar{T}. \end{cases} $$

(iii) We have the following jump identities:

$$ \llbracket S_T[\phi] \rrbracket = 0, \qquad (2.11) $$

$$ \llbracket \nabla S_T[\phi] \mathbf{n} \rrbracket = \phi - (\phi \cdot \mathbf{n})\mathbf{n}, \qquad (2.12) $$

$$ \llbracket \mathcal{P}_T[\phi] \rrbracket = -\phi \cdot \mathbf{n}, \qquad (2.13) $$

In particular,

$$ \phi = \llbracket \sigma(S_T[\phi], \mathcal{P}_T[\phi]) \mathbf{n} \rrbracket. \qquad (2.14) $$

(iv) The following continuity bound hold true:

$$ \|\nabla S_T[\phi]\|_{L^2(\mathbb{R}^d, \mathbb{R}^{d \times d})} + \|\mathcal{P}_T[\phi]\|_{L^2(\mathbb{R}^d \setminus T)} + \|\mathcal{P}_T[\phi]\|_{L^2(T)/\mathbb{R}} \le C \|\phi\|_{H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d)}. \qquad (2.15) $$
Proof. Properties (i) and (iii) are standard and can be found in [16]. For (ii), $S_T[n](x) = 0$ is a straightforward consequence of $\text{div}(G) = 0$ in $\mathbb{R}^d$. Remembering that $\psi = -\nabla\Gamma$, $\mathcal{P}_T[n]$ can be rewritten as

$$ \mathcal{P}_T[n](y) = \int_{\partial T} \psi(x-y) \cdot n(y) \, d\sigma(y) = -\int_{\partial T} \nabla\Gamma(x-y) \cdot n(y) \, d\sigma(y) = -\mathcal{D}_T[1](y), $$

where $\mathcal{D}_T$ denotes the standard Laplace double layer potential. It is then a classical fact [29] that $\mathcal{D}_T[1]$ is equal to 1 inside $T$, 0 outside and 1/2 on $\partial T$. The point (iv) is a consequence of (2.10) and (2.14). The asymptotic expansions of (2.16) and (2.17) are consequences of the fact that $G$ (or $\nabla G$ if $d=2$) and $\psi$ are homogeneous kernels. $\square$

## 2.2. Exterior Stokes problems and layer potential representations

In what follows, we denote by $\mathcal{D}^{1,2}(\mathbb{R}^d, \mathbb{R}^d)$ the so-called Deny-Lions [1] or “homogeneous” space [28, 31]:

$$ \mathcal{D}^{1,2}(\mathbb{R}^d, \mathbb{R}^d) := \{v \text{ measurable } | \nabla v \in L^2(\mathbb{R}^{d \times d})\}. $$

We recall in this section basic existence results for exterior Stokes problems [28, 1, 2, 20] as well as the far field asymptotic expansions of their solutions. The cases $d \ge 3$ and $d=2$ need to be treated separately.

**Proposition 2.3.** Let $d \ge 3$ and $g \in H^{\frac{1}{2}}(\partial D)$ satisfying

$$ \int_{\partial T} g \cdot n \, d\sigma = 0. \qquad (2.18) $$

(i) There exists a unique solution $(v, p) \in \mathcal{D}^{1,2}(\mathbb{R}^d \times \mathbb{R}^d) \times L^2(\mathbb{R}^d)$ to the problem

$$ \left\{ \begin{array}{l} -\Delta v + \nabla p = 0 \quad \text{in } \mathbb{R}^d \setminus \partial T, \\ \text{div}(v) = 0 \quad \text{on } \partial T, \\ v = g \quad \text{on } \partial T, \\ \int_T p(y) \, dy = 0. \end{array} \right. \qquad (2.19) $$

(ii) The solution $(v, p)$ can be represented as a Stokes potential

$$ v = S_T[[\sigma(v, p)n]], \quad p = \mathcal{P}_T[[\sigma(v, p)n]] \quad \text{in } \mathbb{R}^d \setminus \partial T. \qquad (2.20) $$

(iii) $(v, p)$ have the following asymptotic expansion at infinity:

$$ v(y) = \sum_{k=0}^{+\infty} \frac{(-1)^k}{k!} \nabla^k G(y) \cdot \int_{\partial T} [[\sigma(v, p)n]] t^k \, d\sigma(t), $$

$$ p(y) = \sum_{k=0}^{+\infty} \frac{(-1)^k}{k!} \nabla^k \psi(y) \cdot \int_{\partial T} [[\sigma(v, p)n]] t^k \, d\sigma(t). $$

Proof. These results are rather standard [28, 16, 2], but it is not easy to find a reference proposing a unified statement, so we provide a short proof for the convenience of the reader.

(i) The existence result for the interior problem is classical, an existence proof for the exterior problem can be found in [28, 16] with layer potential methods, and a variational proof is available in [2].

(ii) Let us denote by $(\hat{v}, \hat{p})$ the functions defined by

$$ \hat{v} := S_T[[\sigma(v, p)n]], \quad \hat{p} = \mathcal{P}_T[[\sigma(v, p)n]]. \qquad (2.21) $$

Note that $v - \hat{v} \in \mathcal{D}^{1,2}(\mathbb{R}^d \setminus T)$. and $\Delta(v - \hat{v}) - \nabla(p - \hat{p}) = 0$ in a distributional sense. Multiplying this identity by $v - \hat{v}$, integrating by parts and using $[[\sigma(v, p)n]] = [[\sigma(\hat{v}, \hat{p})n]]$ (which follows from the jump identities (2.11) to (2.12)), we infer

$$ \int_{\mathbb{R}^d} |\nabla(v - \hat{v})|^2 dy = 0. $$

We deduce that $v - \hat{v}$ is constant. Due to the decay at infinity, this constant is necessary zero, thus $v = \hat{v}$. Then, $\nabla(p - \hat{p}) = 0$ and we also deduce that $p - \hat{p}$ differ by a constant. Since $p \in L^2(\mathbb{R}^d)$ and $\hat{p}$ converges to zero at infinity, this constant must also be zero. Hence $p = \hat{p}$ and the representation formula (2.20) is proved. (iii) follows directly from (2.20). $\square$

**Remark 2.1.** Without the condition $\int_T p dy$ in (2.19), the pressure is defined only up to the addition of a constant function $c1_T$ supported in the obstacle. As it is visible from the identity (ii) of Proposition 2.2, the layer potential representation (2.20) remains consistent when changing $p$ into $p + c1_T$.

In two dimensions, the existence and uniqueness result for exterior Stokes problems reads as follows.

**Proposition 2.4.** Assume $d=2$ and consider $g \in H^{\frac{1}{2}}(\partial D)$ satisfying $\int_{\partial T} g \cdot n \, d\sigma = 0$.
(v) $S_T[\phi]$ and $P_T[\phi]$ have the following asymptotic behavior at infinity:

$$S_T[\phi](y) = \sum_{k=0}^{+\infty} \frac{(-1)^k}{k!} \nabla^k G(y) \cdot \int_{\partial T} t^k \phi(t) \, d\sigma(t) \text{ as } |y| \to +\infty, \qquad (2.16)$$

$$P_T[\phi](y) = \sum_{k=0}^{+\infty} \frac{(-1)^k}{k!} \nabla^k \psi(y) \cdot \int_{\partial T} t^k \phi(t) \, d\sigma(t) \text{ as } |y| \to +\infty, \qquad (2.17)$$

where these series converge for $|y|$ sufficiently large.

*Proof.* Properties (i) and (iii) are standard and can be found in [17]. For (ii), $S_T[n](x) = 0$ is a straightforward consequence of $\text{div}(G) = 0$ in $\mathbb{R}^d$. Remembering that $\psi = -\nabla\Gamma$, $P_T[n]$ can be rewritten as

$$P_T[n](y) = \int_{\partial T} \psi(x-y) \cdot n(y) \, d\sigma(y) = - \int_{\partial T} \nabla\Gamma(x-y) \cdot n(y) \, d\sigma(y) = -D_T[1](y),$$

where $D_T$ denotes the standard Laplace double layer potential. It is then a classical fact [32] that $D_T[1]$ is equal to 1 inside $T$, 0 outside and $1/2$ on $\partial T$. The point (iv) is a consequence of (2.10) and (2.14). The asymptotic expansions of (2.16) and (2.17) are consequences of the fact that $G$ (or $\nabla G$ if $d=2$) and $\psi$ are homogeneous kernels. $\square$

## 2.2. Exterior Stokes problems and layer potential representations

In what follows, we denote by $\mathcal{D}^{1,2}(\mathbb{R}^d, \mathbb{R}^d)$ the so-called Deny-Lions [1] or “homogeneous” space [30, 35], which is defined as the completion of the space of smooth compactly supported vector fields with respect to the $L^2$-norm of the gradient:

$$\mathcal{D}^{1,2}(\mathbb{R}^d, \mathbb{R}^d) := \overline{\mathcal{C}_c^\infty(\overline{\mathbb{R}^d, \mathbb{R}^d})}^{\|\nabla\cdot\|_{L^2(\mathbb{R}^d, \mathbb{R}^d \times d)}}.$$

We recall in this section basic existence results for exterior Stokes problems [30, 1, 2, 21] as well as the far field asymptotic expansions of their solutions. The cases $d \ge 3$ and $d=2$ need to be treated separately.

**Proposition 2.3.** Let $d \ge 3$ and $g \in H^{\frac{1}{2}}(\partial D)$ satisfying

$$\int_{\partial T} g \cdot n \, d\sigma = 0. \qquad (2.18)$$

(i) There exists a unique solution $(v, p) \in \mathcal{D}^{1,2}(\mathbb{R}^d \times \mathbb{R}^d) \times L^2(\mathbb{R}^d)$ to the problem

$$\left\{
\begin{array}{l}
-\Delta v + \nabla p = 0 \quad \text{in } \mathbb{R}^d \setminus \partial T, \\
\text{div}(v) = 0 \quad \text{on } \partial T, \\
v = g \quad \text{on } \partial T, \\
\int_T p(y) \, dy = 0.
\end{array}
\right.
\qquad (2.19)$$

(ii) The solution $(v, p)$ can be represented as a Stokes potential

$$v = S_T[[\sigma(v, p)n]], \quad p = P_T[[\sigma(v, p)n]] \quad \text{in } \mathbb{R}^d \setminus \partial T. \qquad (2.20)$$

(iii) $(v, p)$ have the following asymptotic expansion at infinity:

$$v(y) = \sum_{k=0}^{+\infty} \frac{(-1)^k}{k!} \nabla^k G(y) \cdot \int_{\partial T} [[\sigma(v, p)n]] t^k \, d\sigma(t),$$

$$p(y) = \sum_{k=0}^{+\infty} \frac{(-1)^k}{k!} \nabla^k \psi(y) \cdot \int_{\partial T} [[\sigma(v, p)n]] t^k \, d\sigma(t).$$

*Proof.* These results are rather standard [30, 17, 2], but it is not easy to find a reference proposing a unified statement, so we provide a short proof for the convenience of the reader.

(i) The existence result for the interior problem is classical, an existence proof for the exterior problem can be found in [30, 17] with layer potential methods, and a variational proof is available in [2].

(ii) Let us denote by $(\hat{v}, \hat{p})$ the functions defined by

$$\hat{v} := S_T[[\sigma(v, p)n]], \quad \hat{p} = P_T[[\sigma(v, p)n]]. \qquad (2.21)$$

Note that $v - \hat{v} \in \mathcal{D}^{1,2}(\mathbb{R}^d \setminus T)$. and $\Delta(v - \hat{v}) - \nabla(p - \hat{p}) = 0$ in a distributional sense. Multiplying this identity by $v - \hat{v}$, integrating by parts and using $[[\sigma(v, p)n]] = [[\sigma(v, \hat{p})n]]$ (which follows from the jump identities (2.11) to (2.12)), we infer

$$\int_{\mathbb{R}^d} |\nabla(v - \hat{v})|^2 dy = 0.$$
(i) There exists a unique solution $(v, p) \in D^{1,2}(\mathbb{R}^2 \times \mathbb{R}^2) \times L^2(\mathbb{R}^2)$ to the problem

$$ egin{cases} -\Delta v + \nabla p = 0 \text{ in } \mathbb{R}^2 \setminus \partial T, \\ \operatorname{div}(v) = 0 \text{ on } \partial T, \\ v = g \text{ on } \partial T, \\ \int_T p(y) \, dy = 0. \end{cases} \tag{2.22} $$

(ii) The “drag force” generated by the fluid surrounding the obstacle (as well as inside the obstacle) vanishes:

$$ \int_{\partial T} \sigma(v, p) n \Big|_{+} \, d\sigma(y) = 0. \tag{2.23} $$

(iii) There exists a constant vector $\mathcal{F}^{\infty}v$ such that the solution $(v, p)$ can be represented as:

$$ v = \mathcal{F}^{\infty}v + S_T [[\sigma(v, p)n]], \quad p = P_T [[\sigma(v, p)n]] \text{ in } \mathbb{R}^2 \setminus \partial T. \tag{2.24} $$

(iv) $(v, p)$ has the following asymptotic expansion at infinity:

$$ v(y) = \mathcal{F}^{\infty}v + \sum_{k=1}^{+\infty} \frac{(-1)^k}{k!} \nabla^k G(y) \cdot \int_{\partial T} [[\sigma(v, p)n]] t^k d\sigma(t), $$

$$ p(y) = \sum_{k=1}^{+\infty} \frac{(-1)^k}{k!} \nabla^k \psi(y) \cdot \int_{\partial T} [[\sigma(v, p)n]] t^k d\sigma(t). $$

*Proof.* (i) The proof of the existence result relying on a variational argument as in [2] still applies as it does not depend on the dimension.

(ii) In two dimensions, constant vector fields belong to the space $D^{1,2}(\mathbb{R}^2, \mathbb{R}^2)$. Thus, multiplying (2.22) by $e_i$ and integrating by parts over $\mathbb{R}^d \setminus T$ yields

$$ 0 = \int_{\mathbb{R}^d \setminus T} (-\Delta v + \nabla p) \cdot e_i \, dy = - \int_{\partial T} e_i \cdot \sigma(v, p) n \, d\sigma(y), \quad \forall 1 \le i \le d, $$

which implies (2.23).

(iii) As in the proof of Proposition 2.3, we denote by $(\hat{v}, \hat{p})$ the functions defined by

$$ \hat{v} := S_T [[\sigma(v, p)n]], \quad \hat{p} = P_T [[\sigma(v, p)n]]. \tag{2.25} $$

By using the property (2.23), we obtain that these functions satisfy the far field estimates

$$ \hat{v}(y) = O(|y|^{-1}) \text{ and } p(y) = O(|y|^{-2}) \text{ as } |y| \to +\infty. $$

This shows, in particular, that $\hat{v} \in D^{1,2}(\mathbb{R}^2, \mathbb{R}^2)$. Reasoning then as in the proof of Proposition 2.3, we can again obtain that $v - \hat{v}$ is a constant, that we denote by $\mathcal{F}^{\infty}v$. Similarly, we obtain that $p$ and $\hat{p}$ differ by a constant which must be zero by invoking $p - \hat{p} \in L^2(\mathbb{R}^d)$.

(iv) follows from (2.24). □

We now introduce a set of exterior solutions $(\mathbf{\Phi}_i, \sigma_i)_{1\le i\le d}$ which play important roles in the remainder of our analysis. In dimension $d \ge 3$, we first introduce the solutions $(\mathbf{\Psi}_j, \tau_j)_{1\le j\le d}$ to the exterior problem (2.19) with $\mathbf{g} = \mathbf{e}_i$:

$$ egin{cases} -\Delta \mathbf{\Psi}_j + \nabla \tau_j = 0 \text{ in } \mathbb{R}^d \setminus \partial T, \\ \operatorname{div}(\mathbf{\Psi}_j) = 0 \text{ in } \mathbb{R}^d \setminus \partial T, \\ \mathbf{\Psi}_j = \mathbf{e}_j \text{ on } \partial T, \\ \mathbf{\Psi}_j(y) + |y|\tau_j(y) = O(|y|^{2-d}) \text{ as } |y| \to +\infty. \end{cases} \tag{2.26} $$

The symmetric positive definite matrix $M$ defined by

$$ M_{ij} := \int_{\mathbb{R}^d \setminus T} \nabla \mathbf{\Psi}_i : \nabla \mathbf{\Psi}_j \, dy = - \int_{\partial T} e_j \cdot \sigma(\mathbf{\Psi}_i, \tau_i) n \, d\sigma(y), \quad 1 \le i, j \le d, \tag{2.27} $$

is called the Brinkmann inverse permeability matrix. It is the analogue of the capacity concept for the Stokes operator [3]. We then introduce the vector and pressure fields $(\mathbf{\Phi}_i)_{1\le i\le d}$ and $(\sigma_i)_{1\le i\le d}$ defined by

$$ \mathbf{\Phi}_i := \sum_{j=1}^{d} M_{ij}^{-1}(e_j - \mathbf{\Psi}_j), \quad \sigma_i := - \sum_{j=1}^{d} M_{ij}^{-1} \tau_j. \tag{2.28} $$

In dimension $d = 2$, the auxiliary functions $(\mathbf{\Phi}_i, \sigma_i)_{1\le i\le d}$ are defined as the solutions to exterior Stokes problems with a logarithmic growth at infinity.

**Proposition 2.5.** Assume $d = 2$.
We deduce that $v - \hat{v}$ is constant. Due to the decay at infinity, this constant is necessary zero, thus $v = \hat{v}$. Then, $\nabla(p - \hat{p}) = 0$ and we also deduce that $p - \hat{p}$ differ by a constant. Since $p \in L^2(\mathbb{R}^d)$ and $\hat{p}$ converges to zero at infinity, this constant must also be zero. Hence $p = \hat{p}$ and the representation formula (2.20) is proved. (iii) follows directly from (2.20). $\square$

**Remark 2.1.** Without the condition $\int_T p dy$ in (2.19), the pressure is defined only up to the addition of a constant function $c1_T$ supported in the obstacle. As it is visible from the identity (ii) of Proposition 2.2, the layer potential representation (2.20) remains consistent when changing $p$ into $p + c1_T$.

In two dimensions, the existence and uniqueness result for exterior Stokes problems reads as follows.

**Proposition 2.4.** Assume $d = 2$ and consider $g \in H^{\frac{1}{2}}(\partial D)$ satisfying $\int_{\partial T} g \cdot n \, d\sigma = 0$.

(i) There exists a unique solution $(v, p) \in D^{1,2}(\mathbb{R}^2 \times \mathbb{R}^2) \times L^2(\mathbb{R}^2)$ to the problem

$$ egin{cases} 
-\Delta v + \nabla p = 0 \quad \text{in } \mathbb{R}^2 \setminus \partial T, \\ 
\operatorname{div}(v) = 0 \quad \text{on } \partial T, \\ 
v = g \quad \text{on } \partial T, \\ 
\int_T p(y) \, dy = 0. 
\end{cases} 
\quad (2.22) $$

(ii) The “drag force” generated by the fluid surrounding the obstacle (as well as inside the obstacle) vanishes:

$$ \int_{\partial T} \sigma(v, p)n \Big|_{+} \, d\sigma(y) = 0. \quad (2.23) $$

(iii) There exists a constant vector $\mathcal{F}^{\infty}v$ such that the solution $(v, p)$ can be represented as:

$$ v = \mathcal{F}^{\infty}v + S_T[[\sigma(v, p)n]], \quad p = \mathcal{P}_T[[\sigma(v, p)n]] \quad \text{in } \mathbb{R}^2 \setminus \partial T. \quad (2.24) $$

(iv) $(v, p)$ has the following asymptotic expansion at infinity:

$$ \begin{aligned} v(y) &= \mathcal{F}^{\infty}v + \sum_{k=1}^{+\infty} \frac{(-1)^k}{k!} \nabla^k G(y) \cdot \int_{\partial T} [[\sigma(v, p)n]] t^k \, d\sigma(t), \\ p(y) &= \sum_{k=1}^{+\infty} \frac{(-1)^k}{k!} \nabla^k \psi(y) \cdot \int_{\partial T} [[\sigma(v, p)n]] t^k \, d\sigma(t). \end{aligned} $$

**Proof.** (i) The proof of the existence result relying on a variational argument as in [2] still applies as it does not depend on the dimension.

(ii) In two dimensions, constant vector fields belong to the space $D^{1,2}(\mathbb{R}^2, \mathbb{R}^2)$. Thus, multiplying (2.22) by $e_i$ and integrating by parts over $\mathbb{R}^d \setminus T$ yields

$$ 0 = \int_{\mathbb{R}^d \setminus T} (-\Delta v + \nabla p) \cdot e_i \, dy = - \int_{\partial T} e_i \cdot \sigma(v, p) n \, d\sigma(y), \quad \forall 1 \le i \le d, $$

which implies (2.23).

(iii) As in the proof of Proposition 2.3, we denote by $(\hat{v}, \hat{p})$ the functions defined by

$$ \hat{v} := S_T[[\sigma(v, p)n]], \quad \hat{p} = \mathcal{P}_T[[\sigma(v, p)n]]. \quad (2.25) $$

By using the property (2.23), we obtain that these functions satisfy the far field estimates

$$ \hat{v}(y) = O(|y|^{-1}) \text{ and } p(y) = O(|y|^{-2}) \text{ as } |y| \to +\infty. $$

This shows, in particular, that $\hat{v} \in D^{1,2}(\mathbb{R}^2, \mathbb{R}^2)$. Reasoning then as in the proof of Proposition 2.3, we can again obtain that $v - \hat{v}$ is a constant, that we denote by $\mathcal{F}^{\infty}v$. Similarly, we obtain that $p$ and $\hat{p}$ differ by a constant which must be zero by invoking $p - \hat{p} \in L^2(\mathbb{R}^d)$.

(iv) follows from (2.24). $\square$

We now introduce a set of exterior solutions $(\mathbf{\Phi}_i, \sigma_i)_{1 \le i \le d}$ which play important roles in the remainder of our analysis. In dimension $d \ge 3$, we first introduce the solutions $(\mathbf{\Psi}_j, \tau_j)_{1 \le j \le d}$ to the exterior problem (2.19) with $g = e_i$:

$$ egin{cases} 
-\Delta \mathbf{\Psi}_j + \nabla \tau_j = 0 \quad \text{in } \mathbb{R}^d \setminus \partial T, \\ 
\operatorname{div}(\mathbf{\Psi}_j) = 0 \quad \text{in } \mathbb{R}^d \setminus \partial T, \\ 
\mathbf{\Psi}_j = e_j \quad \text{on } \partial T, \\ 
\mathbf{\Psi}_j(y) + |y|\tau_j(y) = O(|y|^{2-d}) \quad \text{as } |y| \to +\infty. 
\end{cases} 
\quad (2.26) $$

The symmetric positive definite matrix $M$ defined by

$$ M_{ij} := \int_{\mathbb{R}^d \setminus T} \nabla \mathbf{\Psi}_i : \nabla \mathbf{\Psi}_j \, dy = - \int_{\partial T} e_j \cdot \sigma(\mathbf{\Psi}_i, \tau_i) n \, d\sigma(y), \quad 1 \le i, j \le d, \quad (2.27) $$
(i) There exist unique solutions $(\mathbf{\Phi}_i, \sigma_i)_{1\le i\le 2}$ solutions to the exterior Stokes problem

$$
\left\{
\begin{array}{l}
-\Delta \mathbf{\Phi}_i + \nabla \sigma_i = 0 \quad \text{in } \mathbb{R}^2 \setminus \partial T, \\
\operatorname{div}(\mathbf{\Phi}_i) = 0 \quad \text{in } \mathbb{R}^2 \setminus \partial T, \\
\mathbf{\Phi}_i = 0 \quad \text{on } \partial T, \\
\displaystyle\int_{T} \sigma_i = 0, \\
\mathbf{\Phi}_i(y) \sim \frac{1}{4\pi} \log |y| e_i \quad \text{as } |y| \to +\infty, \\
\sigma_i(y) \sim -\frac{1}{2\pi} \frac{y_i}{|y|^2} \quad \text{as } |x| \to +\infty,
\end{array}
\right.
\qquad (2.29)
$$

with the conditions $\mathbf{\Phi}_i - \mathbf{G}_i \in \mathcal{D}_{1,2}(\mathbb{R}^2 \setminus \bar{T}, \mathbb{R}^2)$ and $\sigma_i - \psi_i \in L^2(\mathbb{R}^2 \setminus \bar{T})$.

(ii) Independently on the obstacle T, the drag force generated by the flow $(\mathbf{\Phi}_i, \sigma_i)$ is equal to $e_i$:

$$
\int_{\partial T} \sigma(\mathbf{\Phi}_i, \sigma_i) \mathbf{n} \, d\sigma = e_i.
$$

(iii) There exists a constant $\mathcal{F}^\infty \mathbf{\Phi}_i$ such that $(\mathbf{\Phi}_i, \sigma_i)$ admits the following layer potential representation:

$$
\mathbf{\Phi}_i = S_T [\sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n}|_+] + \mathcal{F}^\infty \mathbf{\Phi}_i, \quad \sigma_i = \mathcal{P}_T [\sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n}|_+] \quad \text{in } \mathbb{R}^2 \setminus \partial T. \quad (2.30)
$$

(iv) $(\mathbf{\Phi}_i, \sigma_i)$ admit the following far field asymptotic expansions

$$
\mathbf{\Phi}_i(y) = \mathbf{G}_i(y) + \mathcal{F}^\infty \mathbf{\Phi}_i + \sum_{k=1}^{+\infty} \nabla^k G(y) \cdot \int_{\partial T} \sigma(\mathbf{\Phi}_i, \sigma_i) \mathbf{n}|_+ t^k \, d\sigma(t) \text{ as } |y| \to +\infty,
$$

$$
\sigma_i(y) = \psi_i(y) + \sum_{k=1}^{+\infty} \nabla^k \psi(y) \cdot \int_{\partial T} \sigma(\mathbf{\Phi}_i, \sigma_i) \mathbf{n}|_+ t^k \, d\sigma(t) \text{ as } |y| \to +\infty.
$$

Proof. (i) $\mathbf{\Phi}_i$ and $\sigma_i$ can be constructed by setting $\mathbf{\Phi}_i := \mathbf{G}_i + \mathbf{v}_i$ and $\sigma_i := \psi_i + p_i$ where $(\mathbf{v}_i, p_i)$ is the solution to the exterior problem (2.22) with $g := -\mathbf{G}_i$ on $\partial T$.

(ii) Due to (2.23), it holds $\int_{\partial T} \sigma(\mathbf{\Phi}_i, \sigma_i) \mathbf{n}|_+ \, d\sigma = \int_{\partial T} \sigma(\mathbf{G}_i, \psi_i) \mathbf{n} \, d\sigma$. This quantity is equal to $e_i$ in virtue of the point (i) of Lemma 2.1.

(iii) Let us denote by $(\tilde{\mathbf{\Phi}}_i, \hat{\sigma}_i)$ the functions defined by $\tilde{\mathbf{\Phi}}_i := S_T[\sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n}|_+]$ and $\hat{\sigma}_i := \mathcal{P}_T[\sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n}|_+]$. From the layer potential representation, it is straightforward to see that $\nabla \tilde{\mathbf{\Phi}}_i(y) = \nabla \mathbf{G}_i(y) + O(|y|^{-2})$ and $\hat{\sigma}_i = \psi_i(y) + O(|y|^{-2})$, thus $\mathbf{\Phi}_i - \tilde{\mathbf{\Phi}}_i \in \mathcal{D}^{1,2}(\mathbb{R}^2 \setminus T, \mathbb{R}^2)$ and $\sigma_i - \hat{\sigma}_i \in L^2(\mathbb{R}^2 \setminus T)$. Multiplying again the identity $\Delta(v - \hat{v}) - \nabla(\sigma_i - \hat{\sigma}_i) = 0$ by $\mathbf{\Phi}_i - \tilde{\mathbf{\Phi}}_i$ and integrating by parts, we obtain as in the point (ii) of Proposition 2.3 that $\int_{\mathbb{R}^2 \setminus T} |\nabla(\mathbf{\Phi}_i - \tilde{\mathbf{\Phi}}_i)|^2 dy = 0$, i.e. $\mathbf{\Phi}_i - \tilde{\mathbf{\Phi}}_i$ is a constant, that we denote by $\mathcal{F}^\infty \mathbf{\Phi}_i$. Similarly $\sigma_i - \hat{\sigma}_i$ is a constant $c$. Since $\sigma_i - \hat{\sigma}_i \in L^2(\mathbb{R}^2 \setminus T, \mathbb{R}^2)$, we obtain $c = 0$.

(iv) These asymptotic expansions follow from the layer potential representation (2.30).

In all what follows, for any function $v$ admitting an asymptotic expansion at infinity of the form

$$
v(y) = v^{\infty} + \sum_{k=0}^{+\infty} v^k \cdot \nabla^k G(y) \text{ as } |y| \to +\infty, \qquad (2.31)
$$

for a constant $v^\infty \in \mathbb{R}^d$ and $k$-th order vector valued tensors $(v^k)_{k \in \mathbb{R}}$, we will denote by $\mathcal{F}^\infty v$ and $\mathcal{F}^k v$ the coefficients of this expansion:

$$
\mathcal{F}^{\infty} \mathbf{v} := \mathbf{v}^{\infty} \text{ and } \mathcal{F}^{k} \mathbf{v} := \mathbf{v}^{k}. \qquad (2.32)
$$

In the context where $(v, p)$ is the solution to an exterior Stokes problem, we have thus

$$
\mathcal{F}^k v = \frac{(-1)^k}{k!} \int_{\partial T} \llbracket \sigma(v, p) n \rrbracket t^k d\sigma(t) \text{ for all } k \in \mathbb{N}. \quad (2.33)
$$

Moreover, it can be obtained through the Green identity that

$$
\mathcal{F}^{\infty} v = \sum_{i=1}^{d} \left( \int_{\partial T} v \cdot \sigma(\mathbf{\Phi}_i, \sigma_i) n|_+ \, d\sigma \right) e_i, \quad (2.34)
$$

where these formulas are valid independently of $d = 2$ or $d = 3$ thanks to the definition of the functions $(\mathbf{\Phi}_i)_{1\le i\le d}$. We have also in both cases $d = 2$ and $d \ge 3$:

$$
\mathcal{F}^0 \mathbf{\Phi}_i = e_i, \qquad (2.35)
$$

and $\mathcal{F}^\infty \mathbf{\Phi}_i = M^{-1} e_i$ if $d \ge 3$ due to (2.28). In dimension $d = 2$, the constants $\mathcal{F}^\infty \mathbf{\Phi}_i$ are the analogue of the logarithmic capacity concept for the Laplace operator [3].
is called the Brinkmann inverse permeability matrix. It is the analogue of the capacity concept for the Stokes operator [3]. We then introduce the vector and pressure fields $(\mathbf{\Phi}_i)_{1\le i\le d}$ and $(\sigma_i)_{1\le i\le d}$ defined by

$$ \mathbf{\Phi}_i := \sum_{j=1}^{d} M_{ij}^{-1}(e_j - \mathbf{\Psi}_j), \quad \sigma_i := -\sum_{j=1}^{d} M_{ij}^{-1}\tau_j. \qquad (2.28) $$

In dimension $d = 2$, the auxiliary functions $(\mathbf{\Phi}_i, \sigma_i)_{1\le i\le d}$ are defined as the solutions to exterior Stokes problems with a logarithmic growth at infinity.

**Proposition 2.5.** Assume $d = 2$.

(i) There exist unique solutions $(\mathbf{\Phi}_i, \sigma_i)_{1\le i\le 2}$ solutions to the exterior Stokes problem

$$ \left\{ \begin{array}{l} -\Delta \mathbf{\Phi}_i + \nabla \sigma_i = 0 \quad \text{in } \mathbb{R}^2 \setminus \partial T, \\ \operatorname{div}(\mathbf{\Phi}_i) = 0 \quad \text{in } \mathbb{R}^2 \setminus \partial T, \\ \mathbf{\Phi}_i = 0 \quad \text{on } \partial T, \\ \\ \int_{T} \sigma_i = 0, \\ \\ \mathbf{\Phi}_i(y) \sim \frac{1}{4\pi} \log |y| e_i \quad \text{as } |y| \to +\infty, \\ \\ \sigma_i(y) \sim -\frac{y_i}{2\pi |y|^2} \quad \text{as } |x| \to +\infty, \end{array} \right. \qquad (2.29) $$

with the conditions $\mathbf{\Phi}_i - \mathbf{G}_i \in \mathcal{D}_{1,2}(\mathbb{R}^2 \setminus \bar{T}, \mathbb{R}^2)$ and $\sigma_i - \psi_i \in L^2(\mathbb{R}^2 \setminus \bar{T})$.

(ii) Independently on the obstacle T, the drag force generated by the flow $(\mathbf{\Phi}_i, \sigma_i)$ is equal to $e_i$:

$$ \int_{\partial T} \sigma(\mathbf{\Phi}_i, \sigma_i) \mathbf{n} \, d\sigma = e_i. $$

(iii) There exists a constant $\mathcal{F}^\infty \mathbf{\Phi}_i$ such that $(\mathbf{\Phi}_i, \sigma_i)$ admits the following layer potential representation:

$$ \mathbf{\Phi}_i = S_T[\sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n}|_+] + \mathcal{F}^\infty \mathbf{\Phi}_i, \quad \sigma_i = P_T[\sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n}|_+] \quad \text{in } \mathbb{R}^2 \setminus \partial T. \qquad (2.30) $$

(iv) $(\mathbf{\Phi}_i, \sigma_i)$ admit the following far field asymptotic expansions

$$ \mathbf{\Phi}_i(y) = \mathbf{G}_i(y) + \mathcal{F}^\infty \mathbf{\Phi}_i + \sum_{k=1}^{+\infty} \nabla^k G(y) \cdot \int_{\partial T} \sigma(\mathbf{\Phi}_i, \sigma_i) \mathbf{n}|_+ t^k \, d\sigma(t) \quad \text{as } |y| \to +\infty, $$

$$ \sigma_i(y) = \psi_i(y) + \sum_{k=1}^{+\infty} \nabla^k \psi(y) \cdot \int_{\partial T} \sigma(\mathbf{\Phi}_i, \sigma_i) \mathbf{n}|_+ t^k \, d\sigma(t) \quad \text{as } |y| \to +\infty. $$

*Proof.* (i) $\mathbf{\Phi}_i$ and $\sigma_i$ can be constructed by setting $\mathbf{\Phi}_i := \mathbf{G}_i + \mathbf{v}_i$ and $\sigma_i := \psi_i + p_i$ where $(\mathbf{v}_i, p_i)$ is the solution to the exterior problem (2.22) with $g := -\mathbf{G}_i$ on $\partial T$.

(ii) Due to (2.23), it holds $\int_{\partial T} \sigma(\mathbf{\Phi}_i, \sigma_i) \mathbf{n}|_+ \, d\sigma = \int_{\partial T} \sigma(\mathbf{G}_i, \psi_i) \mathbf{n} \, d\sigma$. This quantity is equal to $e_i$ in virtue of the point (i) of Lemma 2.1.

(iii) Let us denote by $(\hat{\mathbf{\Phi}}_i, \hat{\sigma}_i)$ the functions defined by $\hat{\mathbf{\Phi}}_i := S_T[\sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n}|_+]$ and $\hat{\sigma}_i := P_T[\sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n}|_+]$. From the layer potential representation, it is straightforward to see that $\nabla \hat{\mathbf{\Phi}}_i(y) = \nabla \mathbf{G}_i(y) + O(|y|^{-2})$ and $\hat{\sigma}_i = \psi_i(y) + O(|y|^{-2})$, thus $\mathbf{\Phi}_i - \hat{\mathbf{\Phi}}_i \in \mathcal{D}^{1,2}(\mathbb{R}^2 \setminus T, \mathbb{R}^2)$ and $\sigma_i - \hat{\sigma}_i \in L^2(\mathbb{R}^2 \setminus T)$. Multiplying again the identity $\Delta(v - \hat{v}) - \nabla(\sigma_i - \hat{\sigma}_i) = 0$ by $\mathbf{\Phi}_i - \hat{\mathbf{\Phi}}_i$ and integrating by parts, we obtain as in the point (ii) of Proposition 2.3 that $\int_{\mathbb{R}^2 \setminus T} |\nabla(\mathbf{\Phi}_i - \hat{\mathbf{\Phi}}_i)|^2 dy = 0$, i.e. $\mathbf{\Phi}_i - \hat{\mathbf{\Phi}}_i$ is a constant, that we denote by $\mathcal{F}^\infty \mathbf{\Phi}_i$. Similarly $\sigma_i - \hat{\sigma}_i$ is a constant $c$. Since $\sigma_i - \hat{\sigma}_i \in L^2(\mathbb{R}^2 \setminus T, \mathbb{R}^2)$, we obtain $c = 0$.

(iv) These asymptotic expansions follow from the layer potential representation (2.30). $\square$

In all what follows, for any function $\mathbf{v}$ admitting an asymptotic expansion at infinity of the form

$$ \mathbf{v}(y) = \mathbf{v}^\infty + \sum_{k=0}^{+\infty} \mathbf{v}^k \cdot \nabla^k G(y) \quad \text{as } |y| \to +\infty, \qquad (2.31) $$

for a constant $\mathbf{v}^\infty \in \mathbb{R}^d$ and k-th order vector valued tensors $(\mathbf{v}^k)_{k \in \mathbb{R}}$, we will denote by $\mathcal{F}^\infty \mathbf{v}$ and $\mathcal{F}^k \mathbf{v}$ the coefficients of this expansion:

$$ \mathcal{F}^\infty \mathbf{v} := \mathbf{v}^\infty \quad \text{and} \quad \mathcal{F}^k \mathbf{v} := \mathbf{v}^k. \qquad (2.32) $$

In the context where $(\mathbf{v}, p)$ is the solution to an exterior Stokes problem, we have thus

$$ \mathcal{F}^k \mathbf{v} = \frac{(-1)^k}{k!} \int_{\partial T} [\sigma(\mathbf{v}, p) \mathbf{n}] t^k \, d\sigma(t) \quad \text{for all } k \in \mathbb{N}. \qquad (2.33) $$
We conclude this part with some useful invertibility theorems for the single layer potential $S_T$. Let us remark that $S_T$ is not invertible as an operator from $H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d)$ to $H^{\frac{1}{2}}(\partial T, \mathbb{R}^d)$ due to $S_T[n] = 0$. In what follows, we denote for any $s \in \mathbb{R}$ by $V^s$ the space

$$V^s := \{v \in H^s(\partial T, \mathbb{R}^d) \mid \int_{\partial T} v \cdot n \, d\sigma = 0\}.$$

**Proposition 2.6.** (i) If $d \ge 3$, Ker($S_T$) = span($n$) and $S_T$ is an isomorphism from $V^{-\frac{1}{2}}$ to $V^{\frac{1}{2}}$.

(ii) In both cases $d = 2$ or $d \ge 3$, the operator $\mathcal{A}: V^{-\frac{1}{2}} \times \mathbb{R}^d \to V^{\frac{1}{2}} \times \mathbb{R}^d$ defined by

$$\mathcal{A}(\phi, a) = \left( S_T[\phi] + a, \int_{\partial D} \phi \, d\sigma \right) \qquad (2.36)$$

has a bounded inverse.

*Proof.* (i) Assuming $d \ge 3$, we observe that if $S_T[\phi] = 0$ on $\partial T$, then due to (2.10), $S_T[\phi] = 0$ in $\mathbb{R}^d$ and $P_T[\phi] = C1_T$ for some constant $C$. Applying the jump identity (2.14), we obtain that $\phi = Cn$.

To prove that $S_T: V^{-\frac{1}{2}} \to V^{\frac{1}{2}}$ is invertible, it is sufficient to show that it is surjective. Let $g \in V^{\frac{1}{2}}$ and denote by $v$ and $p$ the solution to the Stokes problem (2.19). From the representation formula (2.20), it holds $v = S_T[[\sigma(v, p)n]]$. Thus $S_T[\phi] = g$ with $\phi := [[\sigma(v, p)n]] - \frac{1}{|\partial T|} (\int_{\partial T} [[\sigma(v, p)n]] \cdot n \, d\sigma) n \in V^{-\frac{1}{2}}$.

(iii) Case $d=2$. Let us first show the injectivity. If $(\phi, a) \in \text{Ker}(\mathcal{A})$ and $a$ are such that $S_T[\phi] + a = 0$ and $\int_{\partial D} \phi \, d\sigma = 0$, then the function $v := S_T[\phi]$ and $p := P_T[\phi]$ solve the exterior problem (2.22) with $g = -a$. Furthermore, $v \in D^{1,2}(\mathbb{R}^d, \mathbb{R}^d)$ and $p \in L^2(\mathbb{R}^d)$ due to $\int_{\partial D} \phi \, d\sigma = 0$. It follows that $v = a$ and $p = 0$, but then $\phi = 0$ due to the jump relation (2.14). Consequently, $a = 0$. Thus $\mathcal{A}$ is injective.

To prove that $\mathcal{A}$ is surjective, let us consider arbitrary $g \in V^{\frac{1}{2}}$ and $b \in \mathbb{R}^d$ and consider $(v, p)$ the solution to the exterior Stokes problem (2.22). We pose

$$\phi = \sigma(v, p)n + \sum_{i=1}^{d} b_i \sigma(\Phi_i, \sigma_i)n + cn,$$

where $c$ is a constant adjusted so that $\int_{\partial T} \phi \cdot n \, d\sigma = 0$. We have then on $\partial T$:

$$\begin{aligned} g &= v = S_T[\sigma(v, p)n] + \mathcal{F}^{\infty}v = S_T[\phi] - \sum_{i=1}^{d} b_i S_T[\sigma(\Phi_i, \sigma_i)n] + \mathcal{F}^{\infty}v \\ &= S_T[\phi] - \sum_{i=1}^{d} b_i \mathcal{F}^{\infty}\Phi_i + \mathcal{F}^{\infty}v. \end{aligned}$$

Moreover,

$$\int_{\partial T} \phi \, d\sigma = \sum_{i=1}^{d} b_i \int_{\partial T} \sigma(\Phi_i, \sigma_i) n \, d\sigma = \sum_{i=1}^{d} b_i e_i = b,$$

which shows that $\mathcal{A}(\phi, a) = (g, b)$ with $a = \mathcal{F}^{\infty}v - \sum_{i=1}^{d} b_i \mathcal{F}^{\infty}\Phi_i$, proving the surjectivity.

Case $d \ge 3$. To prove that $\mathcal{A}$ is invertible, we show the existence and uniqueness of a solution to the equation $S_T[\phi] + a = g$ with $g \in V^{\frac{1}{2}}$ and $\int_{\partial T} \phi \, d\sigma = b$. Due to the point (i), we have

$$\phi = (S_T)^{-1}g - \sum_{i=1}^{d} a_i \sigma(\Phi_i, \sigma_i)n + Cn, \qquad (2.37)$$

for a suitable constant $C$ adjusted such that $\int_{\partial T} \phi \cdot n \, d\sigma = 0$. Thus

$$\int_{\partial T} \phi \, d\sigma = \int_{\partial T} (S_T)^{-1}g \, d\sigma + \sum_{i,j=1}^{d} a_i M_{ij} e_j = Ma + \int_{\partial T} (S_T)^{-1}g \, d\sigma.$$

Therefore, $(a, \phi)$ with $a = M^{-1}(b - \int_{\partial T} (S_T)^{-1}g \, d\sigma)$ and $\phi$ determined by (2.37) provides the unique solution to $\mathcal{A}(\phi, a) = (g, b)$. $\square$

## 2.3. Stokes systems with periodic boundary conditions

The asymptotic expansions of the velocity fields $(\mathcal{X}_i)_{1 \le i \le d}$ solving (1.3) feature terms that are solutions to interior Stokes problems posed in the unit cell with some specific periodicity conditions on $\partial P$. In this section, we give an existence and uniqueness theory for these cell problems which are well-posed under some compatibility conditions. Analogue results for Stokes problems posed on the perforated cell $P \setminus (\eta T)$ with periodicity conditions on the boundary of $P$ and a Dirichlet boundary condition on $\eta \partial T$ are stated in Appendix A.
Moreover, it can be obtained through the Green identity that

$$ \mathcal{F}^{\infty} \mathbf{v} = \sum_{i=1}^{d} \left( \int_{\partial T} \mathbf{v} \cdot \sigma(\mathbf{\Phi}_i, \sigma_i) \mathbf{n}|_{+} \mathrm{d}\sigma \right) e_i, \qquad (2.34) $$

where these formulas are valid independently of $d = 2$ or $d = 3$ thanks to the definition of the functions $(\mathbf{\Phi}_i)_{1\le i\le d}$. We have also in both cases $d = 2$ and $d \ge 3$:

$$ \mathcal{F}^0 \mathbf{\Phi}_i = e_i, \qquad (2.35) $$

and $\mathcal{F}^{\infty} \mathbf{\Phi}_i = M^{-1}e_i$ if $d \ge 3$ due to (2.28). In dimension $d = 2$, the constants $\mathcal{F}^{\infty} \mathbf{\Phi}_i$ are the analogue of the logarithmic capacity concept for the Laplace operator [3].

We conclude this part with some useful invertibility theorems for the single layer potential $S_T$. Let us remark that $S_T$ is not invertible as an operator from $H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d)$ to $H^{\frac{1}{2}}(\partial T, \mathbb{R}^d)$ due to $S_T[\mathbf{n}] = 0$. In what follows, we denote for any $s \in \mathbb{R}$ by $V^s$ the space

$$ V^s := \left\{ \mathbf{v} \in H^s(\partial T, \mathbb{R}^d) \mid \int_{\partial T} \mathbf{v} \cdot \mathbf{n} \mathrm{d}\sigma = 0 \right\}. $$

**Proposition 2.6.** (i) If $d \ge 3$, $\mathrm{Ker}(S_T) = \mathrm{span}(\mathbf{n})$ and $S_T$ is an isomorphism from $V^{-\frac{1}{2}}$ to $V^{\frac{1}{2}}$.

(ii) In both cases $d = 2$ or $d \ge 3$, the operator $\mathcal{A}: V^{-\frac{1}{2}} \times \mathbb{R}^d \to V^{\frac{1}{2}} \times \mathbb{R}^d$ defined by

$$ \mathcal{A}(\phi, \mathbf{a}) = \left( S_T[\phi] + \mathbf{a}, \int_{\partial D} \phi \mathrm{d}\sigma \right) \qquad (2.36) $$

has a bounded inverse.

*Proof.* (i) Assuming $d \ge 3$, we observe that if $S_T[\phi] = 0$ on $\partial T$, then due to (2.10), $S_T[\phi] = 0$ in $\mathbb{R}^d$ and $\mathcal{P}_T[\phi] = C1_T$ for some constant $C$. Applying the jump identity (2.14), we obtain that $\phi = Cn$.

To prove that $S_T: V^{-\frac{1}{2}} \to V^{\frac{1}{2}}$ is invertible, it is sufficient to show that it is surjective. Let $\mathbf{g} \in V^{\frac{1}{2}}$ and denote by $\mathbf{v}$ and $p$ the solution to the Stokes problem (2.19). From the representation formula (2.20), it holds $\mathbf{v} = S_T[[\sigma(\mathbf{v}, p)\mathbf{n}]].$ Thus $S_T[\phi] = \mathbf{g}$ with $\phi := [[\sigma(\mathbf{v}, p)\mathbf{n}] - \frac{1}{|\partial T|} (\int_{\partial T} [[\sigma(\mathbf{v}, p)\mathbf{n}] \cdot \mathbf{n} \mathrm{d}\sigma) \mathbf{n}] \in V^{-\frac{1}{2}}$.

(iii) *Case d = 2.* Let us first show the injectivity. If $(\phi, \mathbf{a}) \in \mathrm{Ker}(\mathcal{A})$ and $\mathbf{a}$ are such that $S_T[\phi] + \mathbf{a} = 0$ and $\int_{\partial D} \phi \mathrm{d}\sigma = 0$, then the function $\mathbf{v} := S_T[\phi]$ and $p := \mathcal{P}_T[\phi]$ solve the exterior problem (2.22) with $\mathbf{g} = -\mathbf{a}$. Furthermore, $\mathbf{v} \in \mathcal{D}^{1,2}(\mathbb{R}^d, \mathbb{R}^d)$ and $p \in L^2(\mathbb{R}^d)$ due to $\int_{\partial D} \phi \mathrm{d}\sigma = 0$. It follows that $\mathbf{v} = \mathbf{a}$ and $p = 0$, but then $\phi = 0$ due to the jump relation (2.14). Consequently, $\mathbf{a} = 0$. Thus $\mathcal{A}$ is injective.

To prove that $\mathcal{A}$ is surjective, let us consider arbitrary $\mathbf{g} \in V^{\frac{1}{2}}$ and $\mathbf{b} \in \mathbb{R}^d$ and consider $(\mathbf{v}, p)$ the solution to the exterior Stokes problem (2.22). We pose

$$ \phi = \sigma(\mathbf{v}, p)\mathbf{n} + \sum_{i=1}^{d} b_i \sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n} + c\mathbf{n}, $$

where $c$ is a constant adjusted so that $\int_{\partial T} \phi \cdot \mathbf{n} \mathrm{d}\sigma = 0$. We have then on $\partial T$:

$$ \begin{aligned} \mathbf{g} &= \mathbf{v} = S_T[\sigma(\mathbf{v}, p)\mathbf{n}] + \mathcal{F}^{\infty}\mathbf{v} = S_T[\phi] - \sum_{i=1}^{d} b_i S_T[\sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n}] + \mathcal{F}^{\infty}\mathbf{v} \\ &= S_T[\phi] - \sum_{i=1}^{d} b_i \mathcal{F}^{\infty}\mathbf{\Phi}_i + \mathcal{F}^{\infty}\mathbf{v}. \end{aligned} $$

Moreover,

$$ \int_{\partial T} \phi \mathrm{d}\sigma = \sum_{i=1}^{d} b_i \int_{\partial T} \sigma(\mathbf{\Phi}_i, \sigma_i) \mathbf{n} \mathrm{d}\sigma = \sum_{i=1}^{d} b_i \mathbf{e}_i = \mathbf{b}, $$

which shows that $\mathcal{A}(\phi, \mathbf{a}) = (\mathbf{g}, \mathbf{b})$ with $\mathbf{a} = \mathcal{F}^{\infty}\mathbf{v} - \sum_{i=1}^{d} b_i \mathcal{F}^{\infty}\mathbf{\Phi}_i$, proving the surjectivity.

*Case d $\ge$ 3.* To prove that $\mathcal{A}$ is invertible, we show the existence and uniqueness of a solution to the equation $S_T[\phi] + \mathbf{a} = \mathbf{g}$ with $\mathbf{g} \in V^{\frac{1}{2}}$ and $\int_{\partial T} \phi \mathrm{d}\sigma = \mathbf{b}$. Due to the point (i), we have

$$ \phi = (S_T)^{-1}\mathbf{g} - \sum_{i=1}^{d} a_i \sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n} + C\mathbf{n}, \qquad (2.37) $$

for a suitable constant $C$ adjusted such that $\int_{\partial T} \phi \cdot \mathbf{n} \mathrm{d}\sigma = 0$. Thus

$$ \int_{\partial T} \phi \mathrm{d}\sigma = \int_{\partial T} (S_T)^{-1} \mathbf{g} \mathrm{d}\sigma + \sum_{i,j=1}^{d} a_i M_{ij} \mathbf{e}_i = M \mathbf{a} + \int_{\partial T} (S_T)^{-1} \mathbf{g} \mathrm{d}\sigma. $$
In what follows, we say that a measurable function $v \in L^1(\partial P, \mathbb{R}^d)$ is “P-periodic” if its values match on opposite faces of the hypercube $P$. This property can be written in the following weak form:

$$ \forall \phi \in C_{\text{per}}^{\infty}(P), \quad \int_{\partial P} (\phi \cdot v) n \, d\sigma = 0, $$

which allows to extend the notion of P-periodicity to any function $v \in H^s(\partial P, \mathbb{R}^d)$ for any $s \in \mathbb{R}$. In the next proposition and hereafter, $L^2(P\setminus(\eta T))/\mathbb{R}$ is the quotient space of square integrable functions on $P\setminus(\eta T)$ defined up to an additive constant. We have the following existence and uniqueness result for periodic Stokes problems (in which we impose periodicity conditions for the trace of the solution and its co-normal derivative on opposite matching faces).

**Proposition 2.7.** For any $f \in H_{\text{per}}^{-1}(P, \mathbb{R}^d)$, $g \in L^2(P)$, $h_0 \in H^{\frac{1}{2}}(\partial P, \mathbb{R}^d)$ and $h_1 \in H^{-\frac{1}{2}}(\partial P, \mathbb{R}^d)$, the problem

$$ \left\{ \begin{array}{ll} -\Delta \mathbf{u} + \nabla p = \mathbf{f} & \text{in } P, \\ \operatorname{div}(\mathbf{u}) = \mathbf{g} & \text{in } P, \\ (\sigma(\mathbf{u}, p)\mathbf{n} - \mathbf{h}_0)\mathbf{n}^T & \text{is P-periodic,} \end{array} \right. \qquad (2.38) $$

admits a unique solution $(\mathbf{u}, p) \in H^1(P, \mathbb{R}^d) \times L^2(P)/\mathbb{R}$ with $\mathbf{u}$ and $p$ defined up to constants if and only if the following compatibility conditions are satisfied:

$$ \int_{\partial P} \mathbf{h}_0 \cdot \mathbf{n} \, d\sigma - \int_{P} \mathbf{g} \, dx = 0 \quad \text{and} \quad \int_{P} \mathbf{f} \, dx + \int_{\partial P} \mathbf{h}_1 \, d\sigma = 0. \qquad (2.39) $$

Moreover, the following properties hold:

(i) $\mathbf{u} - \tilde{\mathbf{h}}_0 \in H_{\text{per}}^1(P, \mathbb{R}^d)$ for any extension $\tilde{\mathbf{h}}_0 \in H^1(P, \mathbb{R}^d)$ of $\mathbf{h}_0$ satisfying $\mathbf{h} = \mathbf{h}_0$ on $\partial P$,

(ii) there exists a constant $C > 0$ independent of $F, \mathbf{h}_0$ and $\mathbf{h}_1$ such that

$$ \|\nabla \mathbf{u}\|_{L^2(P, \mathbb{R}^{d \times d})} + \|p\|_{L^2(P)/\mathbb{R}} \le C (\|\mathbf{F}\|_{H^{-1}(P, \mathbb{R}^d)} + \|\mathbf{g}\|_{L^2(P)} + \|\mathbf{h}_0\|_{H^{\frac{1}{2}}(\partial P, \mathbb{R}^d)} + \|\mathbf{h}_1\|_{H^{-\frac{1}{2}}(\partial P, \mathbb{R}^d)}). \quad (2.40) $$

*Proof.* We only sketch the proof which is obtained by adapting the arguments of [18, Proposition 5.3] for periodic Poisson problems. We first assume that a solution $(\mathbf{u}, p)$ to (2.38) is given. Using Stokes' theorem, we have

$$ \int_P \mathbf{g} \, dx = \int_P \operatorname{div}(\mathbf{u}) \, dx = \int_{\partial P} \mathbf{u} \cdot \mathbf{n} \, d\sigma = \int_{\partial P} (\mathbf{u} - \mathbf{h}_0) \cdot \mathbf{n} \, d\sigma + \int_{\partial P} \mathbf{h}_0 \cdot \mathbf{n} \, d\sigma. $$

Since $\int_{\partial P}(\mathbf{u} - \mathbf{h}_0) \cdot \mathbf{n} \, d\sigma = 0$ due to the periodicity condition, we infer that the condition $\int_{\partial P} \mathbf{h}_0 \cdot \mathbf{n} \, d\sigma = \int_P \mathbf{g} \, dx$ is necessary. Multiplying now (2.38) by a test function $v \in H_{\text{per}}^1(P, \mathbb{R}^d)$ and integrating by parts, we obtain

$$ \int_P \mathbf{f} \cdot \mathbf{v} \, dx = \int_P (-\Delta \mathbf{u} + \nabla p) \cdot \mathbf{v} \, dx = \int_P (\nabla \mathbf{u} : \nabla \mathbf{v} - p \operatorname{div}(\mathbf{v})) \, dx - \int_{\partial P} \mathbf{v} \cdot \sigma(u, p) \mathbf{n} \, d\sigma. $$

Moreover, the second P-periodicity condition implies $\int_{\partial P} \sigma(\mathbf{u}, p)\mathbf{n} \cdot \mathbf{v} \, d\sigma = \int_{\partial P} \mathbf{h}_1 \cdot \mathbf{v} \, d\sigma$ and we deduce that $\mathbf{u}$ satisfies the variational formulation.

$$ \int_P \nabla \mathbf{u} : \nabla \mathbf{v} \, dx = \int_P \mathbf{f} \cdot \mathbf{v} \, dx + \int_{\partial P} \mathbf{h}_1 \cdot \mathbf{v} \, d\sigma \quad \text{for any } \mathbf{v} \in H_{\text{per}}^1(P, \mathbb{R}^d) \text{ such that } \operatorname{div}(\mathbf{v}) = 0. \qquad (2.41) $$

Setting $\mathbf{v} = e_i$ for $1 \le i \le d$ implies in particular that the second condition of (2.39) is also necessary.

Let now $\mathbf{f}, \mathbf{g}, \mathbf{h}_0$ and $\mathbf{h}_1$ be given such that the conditions of (2.39) are fulfilled. The equality $\int_{\partial P} \mathbf{h}_0 \cdot \mathbf{n} \, d\sigma = \int_P \mathbf{g} \, dx$ enables to construct a lift $\tilde{\mathbf{h}}_0 \in H^1(P, \mathbb{R}^d)$ satisfying $\tilde{\mathbf{h}}_0 = \mathbf{h}_0$ on $\partial P$ and $\operatorname{div}(\tilde{\mathbf{h}}_0) = \mathbf{g}$ in $P$. Then, the existence of the solution $\mathbf{u}$ follows by applying the Lax-Milgram theorem on a variational formulation similar to (2.41) satisfied by $\mathbf{u} - \tilde{\mathbf{h}}_0$, using Poincaré-Wireinger inequality. The pressure field $p$ can then be reconstructed by standard arguments involving the existence of a continuous right-inverse for the divergence operator mapping $H_{\text{per}}^1(P, \mathbb{R}^d)$ to $L_0^2(P)$ [21]. The bound (2.40) follows from standard energy estimates. $\square$

We then state the following elliptic regularity result.

**Proposition 2.8.** Let $\mathbf{f} \in C_{\text{per}}^{\infty}(P, \mathbb{R}^d)$ be a smooth periodic vector field. Assume that there exists an open set $T \subset P$ and $\mathbf{h} \in C^{\infty}(P\setminus T, \mathbb{R}^d)$ and $q \in C^{\infty}(P\setminus T)$ such that:

(i) $\int_{\partial P} \mathbf{h} \cdot \mathbf{n} \, d\sigma = 0,$

(ii) $-\Delta \mathbf{h} + \nabla q$ is a smooth P-periodic function:

$$ -\Delta \mathbf{h} + \nabla q \in C_{\text{per}}^{\infty}(P\setminus T, \mathbb{R}^d), $$

(iii) $\mathbf{h}$ satisfies the compatibility condition

$$ \int_P \mathbf{f} \, dx + \int_{\partial P} \sigma(\mathbf{h}, q) \mathbf{n} \, d\sigma = 0. $$
Therefore, $(\mathbf{a}, \phi)$ with $\mathbf{a} = M^{-1}(\mathbf{b} - \int_{\partial T}(S_T)^{-1}\mathbf{g} d\sigma)$ and $\phi$ determined by (2.37) provides the unique solution to $\mathcal{A}(\phi, \mathbf{a}) = (\mathbf{g}, \mathbf{b})$. $\square$

## 2.3. Stokes systems with periodic boundary conditions

The asymptotic expansions of the velocity fields $(\mathcal{X}_i)_{1\le i\le d}$ solving (1.3) feature terms that are solutions to interior Stokes problems posed in the unit cell with some specific periodicity conditions on $\partial P$. In this section, we give an existence and uniqueness theory for these cell problems which are well-posed under some compatibility conditions. Analogue results for Stokes problems posed on the perforated cell $P\setminus(\eta T)$ with periodicity conditions on the boundary of $P$ and a Dirichlet boundary condition on $\eta\partial T$ are stated in Appendix A.

In what follows, we say that a measurable function $v \in L^1(\partial P, \mathbb{R}^d)$ is “P-periodic” if its values match on opposite faces of the hypercube $P$. This property can be written in the following weak form:

$$ \forall \phi \in C_{\text{per}}^{\infty}(P), \quad \int_{\partial P} (\phi \cdot \mathbf{v}) \mathbf{n} \, d\sigma = 0, $$

which allows to extend the notion of P-periodicity to any function $\mathbf{v} \in H^s(\partial P, \mathbb{R}^d)$ for any $s \in \mathbb{R}$. In the next proposition and hereafter, $L^2(P\setminus(\eta T))/\mathbb{R}$ is the quotient space of square integrable functions on $P\setminus(\eta T)$ defined up to an additive constant. We have the following existence and uniqueness result for periodic Stokes problems (in which we impose periodicity conditions for the trace of the solution and its co-normal derivative on opposite matching faces).

**Proposition 2.7.** For any $\mathbf{f} \in H_{\text{per}}^{-1}(P, \mathbb{R}^d)$, $\mathbf{g} \in L^2(P)$, $\mathbf{h}_0 \in H_{\text{per}}^{1/2}(\partial P, \mathbb{R}^d)$ and $\mathbf{h}_1 \in H_{\text{per}}^{-1/2}(\partial P, \mathbb{R}^d)$, the problem

$$ \left\{ \begin{array}{ll} -\Delta \mathbf{u} + \nabla p = \mathbf{f} & \text{in } P, \\ \text{div}(\mathbf{u}) = \mathbf{g} & \text{in } P, \\ \mathbf{u} - \mathbf{h}_0 & \text{is P-periodic,} \\ (\sigma(\mathbf{u}, p)\mathbf{n} - \mathbf{h}_1)\mathbf{n}^T & \text{is P-periodic,} \end{array} \right. \qquad (2.38) $$

admits a unique solution $(\mathbf{u}, p) \in H^1(P, \mathbb{R}^d) \times L^2(P)/\mathbb{R}$ with $\mathbf{u}$ and $p$ defined up to constants if and only if the following compatibility conditions are satisfied:

$$ \int_{\partial P} \mathbf{h}_0 \cdot \mathbf{n} \, d\sigma - \int_P \mathbf{g} \, dx = 0 \quad \text{and} \quad \int_P \mathbf{f} \, dx + \int_{\partial P} \mathbf{h}_1 \, d\sigma = 0. \qquad (2.39) $$

Moreover, the following properties hold:

(i) $\mathbf{u} - \tilde{\mathbf{h}}_0 \in H_{\text{per}}^1(P, \mathbb{R}^d)$ for any extension $\tilde{\mathbf{h}}_0 \in H^1(P, \mathbb{R}^d)$ of $\mathbf{h}_0$ satisfying $\mathbf{h} = \mathbf{h}_0$ on $\partial P$,

(ii) there exists a constant $C > 0$ independent of $\mathbf{F}, \mathbf{h}_0$ and $\mathbf{h}_1$ such that

$$ \|\nabla \mathbf{u}\|_{L^2(P, \mathbb{R}^{d \times d})} + \|p\|_{L^2(P)/\mathbb{R}} \le C (\|\mathbf{F}\|_{H^{-1}(P, \mathbb{R}^d)} + \|\mathbf{g}\|_{L^2(P)} + \|\mathbf{h}_0\|_{H^{1/2}(\partial P, \mathbb{R}^d)} + \|\mathbf{h}_1\|_{H^{-1/2}(\partial P, \mathbb{R}^d)}). \quad (2.40) $$

**Proof.** We only sketch the proof which is obtained by adapting the arguments of [19, Proposition 5.3] for periodic Poisson problems. We first assume that a solution $(\mathbf{u}, p)$ to (2.38) is given. Using Stokes' theorem, we have

$$ \int_P \mathbf{g} \, dx = \int_P \text{div}(\mathbf{u}) \, dx = \int_{\partial P} \mathbf{u} \cdot \mathbf{n} \, d\sigma = \int_{\partial P} (\mathbf{u} - \mathbf{h}_0) \cdot \mathbf{n} \, d\sigma + \int_{\partial P} \mathbf{h}_0 \cdot \mathbf{n} \, d\sigma. $$

Since $\int_{\partial P}(\mathbf{u} - \mathbf{h}_0) \cdot \mathbf{n} \, d\sigma = 0$ due to the periodicity condition, we infer that the condition $\int_{\partial P} \mathbf{h}_0 \cdot \mathbf{n} \, d\sigma = \int_P \mathbf{g} \, dx$ is necessary. Multiplying now (2.38) by a test function $\mathbf{v} \in H_{\text{per}}^1(P, \mathbb{R}^d)$ and integrating by parts, we obtain

$$ \int_P \mathbf{f} \cdot \mathbf{v} \, dx = \int_P (-\Delta \mathbf{u} + \nabla p) \cdot \mathbf{v} \, dx = \int_P (\nabla \mathbf{u} : \nabla \mathbf{v} - p \text{div}(\mathbf{v})) \, dx - \int_{\partial P} \mathbf{v} \cdot \sigma(u, p) \mathbf{n} \, d\sigma. $$

Moreover, the second P-periodicity condition implies $\int_{\partial P} \sigma(\mathbf{u}, p)\mathbf{n} \cdot \mathbf{v} \, d\sigma = \int_{\partial P} \mathbf{h}_1 \cdot \mathbf{v} \, d\sigma$ and we deduce that $\mathbf{u}$ satisfies the variational formulation.

$$ \int_P \nabla \mathbf{u} : \nabla \mathbf{v} \, dx = \int_P \mathbf{f} \cdot \mathbf{v} \, dx + \int_{\partial P} \mathbf{h}_1 \cdot \mathbf{v} \, d\sigma \quad \text{for any } \mathbf{v} \in H_{\text{per}}^1(P, \mathbb{R}^d) \text{ such that } \text{div}(\mathbf{v}) = 0. \quad (2.41) $$

Setting $\mathbf{v} = e_i$ for $1 \le i \le d$ implies in particular that the second condition of (2.39) is also necessary.

Let now $\mathbf{f}, \mathbf{g}, \mathbf{h}_0$ and $\mathbf{h}_1$ be given such that the conditions of (2.39) are fulfilled. The equality $\int_{\partial P} \mathbf{h}_0 \cdot \mathbf{n} \, d\sigma = \int_P \mathbf{g} \, dx$ enables to construct a lift $\tilde{\mathbf{h}}_0 \in H^1(P, \mathbb{R}^d)$ satisfying $\tilde{\mathbf{h}}_0 = \mathbf{h}_0$ on $\partial P$ and $\text{div}(\tilde{\mathbf{h}}_0) = \mathbf{g}$ in $P$. Then, the existence of the solution $\mathbf{u}$ follows by applying the Lax-Milgram theorem on a variational formulation similar to (2.41) satisfied by $\mathbf{u} - \tilde{\mathbf{h}}_0$, using Poincaré-Wirtinger inequality. The pressure field $p$ can then be reconstructed by standard arguments involving the existence of a continuous right-inverse for the divergence operator mapping $H_{\text{per}}^1(P, \mathbb{R}^d)$ to $L_0^2(P)$ [22]. The bound (2.40) follows from standard energy estimates. $\square$

We then state the following elliptic regularity result.
Then the unique solution $(\mathbf{u}, p)$ defined up to a constant to

$$ egin{cases} -\Delta \mathbf{u} + \nabla p = \mathbf{f} \quad \text{in } P, \\ \text{div}(\mathbf{u}) = 0 \quad \text{in } P, \\ \mathbf{u} - \mathbf{h} \quad \text{is P-periodic,} \\ \sigma(\mathbf{u}, p)\mathbf{n}\mathbf{n}^T - \sigma(\mathbf{h}, q)\mathbf{n}\mathbf{n}^T \quad \text{is P-periodic,} \end{cases} \tag{2.42} $$

satisfies the following properties:

a) $\mathbf{u} - \mathbf{h} \in C_{\text{per}}^{\infty}(P \setminus T)$ is a smooth P-periodic function on $P \setminus T$,

b) for any $k > d/2 - 1$, there exists a constant $C_k > 0$ which depends only on $k$ such that

$$ \left\| \mathbf{u} - \int_P \mathbf{u} \, d\mathbf{x} \right\|_{L^\infty(P, \mathbb{R}^d)} + \|\nabla \mathbf{u}\|_{L^\infty(P, \mathbb{R}^d \times d)} + \|p\|_{L^\infty(P)/\mathbb{R}} \le C_k (\|\mathbf{f}\|_{H^k(P, \mathbb{R}^d)} + \|\mathbf{h}\|_{H^{k+2}(P \setminus T, \mathbb{R}^d)} + \|q\|_{H^{k+1}(P \setminus T)}). \tag{2.43} $$

Proof. Let $\theta \in C^{\infty}(P)$ be a smooth cut-off function satisfying $\theta = 1$ on an open neighborhood of $\partial P$ and compactly supported in $P \setminus T$. Then $(\tilde{\mathbf{u}}, \tilde{p}) := (\mathbf{u} - \theta\mathbf{h}, p - \theta q)$ solves the periodic Stokes problem

$$ egin{cases} -\Delta \tilde{\mathbf{u}} + \nabla \tilde{p} = \mathbf{f} + \Delta(\theta \mathbf{h}) - \nabla(\theta q), \quad \text{in } P, \\ \text{div}(\tilde{\mathbf{u}}) = \text{div}(\theta \mathbf{h}) \quad \text{in } P, \\ \tilde{\mathbf{u}} \quad \text{is P-periodic,} \\ \sigma(\mathbf{u}, p)\mathbf{n}\mathbf{n}^T \quad \text{is P-periodic.} \end{cases} \tag{2.44} $$

The compatibility conditions (2.39) for this problem are satisfied due to $\int_P \text{div}(\theta \mathbf{h}) \, d\mathbf{x} = \int_{\partial P} \mathbf{h} \cdot \mathbf{n} \, d\sigma = 0$ and assumption (iii). This problem has therefore a unique solution defined up to a constant in $\tilde{\mathbf{u}}$ which can be solved explicitly using Fourier series. Due to the smoothness of $\mathbf{f}$, $\mathbf{h}$ and $\theta$, it is readily obtained from the properties of the Fourier expansion that $(\mathbf{u}, p)$ is smooth in $P$ and the energy estimate (2.43) holds. $\square$

## 2.4. Periodic Green matrix

The asymptotic analysis of the permeability matrix $\mathcal{X}^*$ of (1.5) brings into play some universal constants $c$ and $k$-th order tensors $\mathfrak{E}^k$, $\mathfrak{R}^k$ which are related to the periodic Green function. We define the periodic fundamental solutions of the Stokes system as the vector and scalar fields $(G_{\text{per},i}, \psi_{\text{per},i}(x))_{1 \le i \le d}$ satisfying

$$ egin{cases} \Delta G_{\text{per},i} - \nabla \psi_{\text{per},i} = \left( \sum_{n \in \mathbb{Z}^d} \delta(\cdot - n) - 1 \right) e_i, \\ \text{div}(G_{\text{per},i}) = 0, \quad 1 \le i \le d. \end{cases} \tag{2.45} $$

These functions admit the following trigonometric series representations [23]:

$$ egin{cases} G_{\text{per},i}(x) = -\frac{1}{4\pi^2} \sum_{k \in \mathbb{Z}^d \setminus \{0\}} \left( \frac{e_i}{|k|^2} - \frac{k_i k}{|k|^4} \right) e^{2i\pi k \cdot x}, \\ \psi_{\text{per},i}(x) = -\frac{1}{2i\pi} \sum_{k \in \mathbb{Z}^d \setminus \{0\}} \frac{k_i}{|k|^2} e^{2i\pi k \cdot x}. \end{cases} \tag{2.46} $$

Classically, the periodic Green functions $(G_{\text{per}}(x), \psi_{\text{per}}(x))$ are smooth perturbations of the hydrodynamics Green functions $(G, \psi)$.

**Proposition 2.9.** There exists a constant $c \in \mathbb{R}$ depending only on the dimension $d$ such that

$$ G_{\text{per},i}(x) = G_i(x) + R_i(x) + c e_i, \quad \psi_{\text{per},i}(x) = \psi_i(x) + s_i(x), \quad \text{for all } x \in P \setminus \{0\}, \tag{2.47} $$

for any $1 \le i \le d$, where $(R_i, s_i) \in H^1(P, \mathbb{R}^d) \times L^2(P)$ is the solution to the periodic difference problem

$$ egin{cases} -\Delta R_i + \nabla s_i = e_i \quad \text{in } P, \\ \text{div}(R_i) = 0 \quad \text{in } P, \\ R_i + G_i \quad \text{is P-periodic,} \\ \sigma(R_i, s_i)\mathbf{n}\mathbf{n}^T + \sigma(G_i, \psi_i)\mathbf{n}\mathbf{n}^T \quad \text{is P-periodic,} \\ R_i(0) = 0, \quad s_i(0) = 0. \end{cases} \tag{2.48} $$

The functions $R_i$ and $s_i$ are respectively even and odd.
**Proposition 2.8.** Let $\mathbf{f} \in C_{\text{per}}^{\infty}(P, \mathbb{R}^d)$ be a smooth periodic vector field. Assume that there exists an open set $T \subset P$ and $\mathbf{h} \in C^{\infty}(P\setminus T, \mathbb{R}^d)$ and $q \in C^{\infty}(P\setminus T)$ such that:

(i) $\int_{\partial P} \mathbf{h} \cdot \mathbf{n} \, d\sigma = 0,$

(ii) $-\Delta \mathbf{h} + \nabla q$ is a smooth $P$-periodic function:

$$-\Delta \mathbf{h} + \nabla q \in C_{\text{per}}^{\infty}(P\setminus T, \mathbb{R}^d),$$

(iii) $\mathbf{h}$ satisfies the compatibility condition

$$\int_P \mathbf{f} \, dx + \int_{\partial P} \sigma(\mathbf{h}, q) \mathbf{n} \, d\sigma = 0.$$

Then the unique solution $(\mathbf{u}, p)$ defined up to a constant to

$$\left\{
\begin{array}{ll}
-\Delta \mathbf{u} + \nabla p = \mathbf{f} & \text{in } P, \\
\operatorname{div}(\mathbf{u}) = 0 & \text{in } P, \\
\mathbf{u} - \mathbf{h} & \text{is } P\text{-periodic,}
\end{array}
\right.
\qquad (2.42) \\
\sigma(\mathbf{u}, p)\mathbf{n}\mathbf{n}^T - \sigma(\mathbf{h}, q)\mathbf{n}\mathbf{n}^T \text{ is } P\text{-periodic,}$$

satisfies the following properties:

a) $\mathbf{u} - \mathbf{h} \in C_{\text{per}}^{\infty}(P\setminus T)$ is a smooth $P$-periodic function on $P\setminus T$,

b) for any $k > d/2 - 1$, there exists a constant $C_k > 0$ which depends only on $k$ such that

$$\left\| \mathbf{u} - \int_P \mathbf{u} \, dx \right\|_{L^\infty(P, \mathbb{R}^d)} + \|\nabla \mathbf{u}\|_{L^\infty(P, \mathbb{R}^d \times d)} + \|p\|_{L^\infty(P)/\mathbb{R}} \le C_k (\|\mathbf{f}\|_{H^k(P, \mathbb{R}^d)} + \|\mathbf{h}\|_{H^{k+2}(P\setminus T, \mathbb{R}^d)} + \|q\|_{H^{k+1}(P\setminus T)}). \quad (2.43)$$

*Proof.* Let $\theta \in C^{\infty}(P)$ be a smooth cut-off function satisfying $\theta = 1$ on an open neighborhood of $\partial P$ and compactly supported in $P\setminus T$. Then $(\tilde{\mathbf{u}}, \tilde{p}) := (\mathbf{u} - \theta\mathbf{h}, p - \theta q)$ solves the periodic Stokes problem

$$\left\{
\begin{array}{ll}
-\Delta \tilde{\mathbf{u}} + \nabla \tilde{p} = \mathbf{f} + \Delta(\theta \mathbf{h}) - \nabla(\theta q), & \text{in } P, \\
\operatorname{div}(\tilde{\mathbf{u}}) = \operatorname{div}(\theta \mathbf{h}) & \text{in } P, \\
\tilde{\mathbf{u}} & \text{is } P\text{-periodic,} \\
\sigma(\mathbf{u}, p)\mathbf{n}\mathbf{n}^T & \text{is } P\text{-periodic.}
\end{array}
\right.
\qquad (2.44)$$

The compatibility conditions (2.39) for this problem are satisfied due to $\int_P \operatorname{div}(\theta \mathbf{h}) \, dx = \int_{\partial P} \mathbf{h} \cdot \mathbf{n} \, d\sigma = 0$ and assumption (iii). This problem has therefore a unique solution defined up to a constant in $\tilde{\mathbf{u}}$ which can be solved explicitly using Fourier series. Due to the smoothness of $\mathbf{f}$, $\mathbf{h}$ and $\theta$, it is readily obtained from the properties of the Fourier expansion that $(\mathbf{u}, p)$ is smooth in $P$ and the energy estimate (2.43) holds. $\square$

## 2.4. Periodic Green matrix

The asymptotic analysis of the permeability matrix $\mathcal{X}^*$ of (1.5) brings into play some universal constants $c$ and $k$-th order tensors $\mathfrak{G}^k$, $\mathfrak{R}^k$ which are related to the periodic Green function. We define the periodic fundamental solutions of the Stokes system as the vector and scalar fields $(\mathbf{G}_{\text{per},i}, \psi_{\text{per},i}(x))_{1\le i\le d}$ satisfying

$$\left\{
\begin{array}{ll}
\Delta \mathbf{G}_{\text{per},i} - \nabla \psi_{\text{per},i} = \left( \sum_{n \in \mathbb{Z}^d} \delta(\cdot - n) - 1 \right) e_i, & \\
\operatorname{div}(\mathbf{G}_{\text{per},i}) = 0, & 1 \le i \le d. \qquad (2.45) \\
\mathbf{G}_{\text{per},i}, \psi_{\text{per},i} \text{ are } P\text{-periodic,} \\
\int_P \mathbf{G}_{\text{per},i} \, dx = 0, \quad \int_P \psi_{\text{per},i} \, dx = 0,
\end{array}
\right.
$$

These functions admit the following trigonometric series representations [24]:

$$\left\{
\begin{array}{l}
\mathbf{G}_{\text{per},i}(x) = -\frac{1}{4\pi^2} \sum_{k \in \mathbb{Z}^d \setminus \{0\}} \left( \frac{e_i}{|k|^2} - \frac{k_i \mathbf{k}}{|k|^4} \right) e^{2i\pi k \cdot x}, \\
\psi_{\text{per},i}(x) = -\frac{1}{2i\pi} \sum_{k \in \mathbb{Z}^d \setminus \{0\}} \frac{k_i}{|k|^2} e^{2i\pi k \cdot x}.
\end{array}
\right.
\qquad (2.46)$$

Classically, the periodic Green functions $(G_{\text{per}}(x), \psi_{\text{per}}(x))$ are smooth perturbations of the hydrodynamics Green functions $(G, \psi)$.
Proof. The compatibility condition $\int_P G_i \cdot n d\sigma = 0$ is satisfied due to $\text{div}(G_i) = 0$. Moreover,

$$ \int_P e_i dx + \int_{\partial P} \sigma(G_i, \psi_i)n d\sigma = 0 $$

holds in virtue of Lemma 2.1. Therefore, Proposition 2.8 implies that (2.48) has a unique solution. Moreover, $R_i + G_i$ and $\psi_i + s_i$ can be extended by periodicity into a smooth function on $\mathbb{R}^d \setminus \mathbb{Z}^d$ and satisfy (2.45) by construction. By uniqueness, we deduce that there exists constant vectors $c_i \in \mathbb{R}^d$ and constants $\kappa_i \in \mathbb{R}$ such that

$$ G_{\text{per},i}(x) = G_i(x) + R_i(x) + c_i, \quad \psi_{\text{per},i}(x) = \psi_i(x) + s_i(x) + \kappa_i. \quad (2.49) $$

Let us then consider the cubic cell symmetries $S^i := \text{I} - 2e_i e_i^T$ which changes the coordinate $x_i$ into $-x_i$, and the symmetry with respect to the diagonal axes $S^{ij} := \text{I} - e_i e_i^T - e_j e_j^T + e_i e_j^T + e_j e_i^T$ for $i \neq j$, which exchange $x_i$ and $x_j$ for $1 \le i \neq j \le d$. These symmetries preserve the hypercube $P$. By using the spherical symmetry properties of the Stokes fundamental solution $(G_i, \psi_i)$, it can be verified that

$$ S R_i \circ S = \sum_{l=1}^{d} S_{il} R_l, \quad s_i \circ S = \sum_{l=1}^{d} S_{il} s_l, $$

see e.g. for a method of proof [17, Section 3.3]. In particular, with $S = S^i$, these properties imply

$$ \int_P R_i \cdot e_j dx = 0 \text{ for all } 1 \le j \neq i \le d \text{ and } \int_P s_i dx = 0 \text{ for all } 1 \le i \le d. $$

Integrating (2.49) over the cell $P$, we deduce

$$ c_i \cdot e_j = 0 \text{ for all } 1 \le j \neq i \le d \text{ and } \kappa_i = 0 \text{ for all } 1 \le i \le d. $$

This shows that $c_i = c_i e_i$ for some constant $c_i$. Considering now the symmetries $S^{ij}$, we find

$$ \int_P R_i \cdot e_i dx = \int_P R_j \cdot e_j dx \text{ for all } 1 \le i, j \le d, $$

which imply that the constants $c_i$ are identically equal to a single constant $c$. The fact that $\kappa_i = 0$ for $1 \le i \le d$ can be obtained by observing that $\psi_{\text{per},i}$ and $\psi_i$ are odd functions and evaluating (2.49) at zero. $\square$

In the remainder of the paper, we will denote by $G_{\text{per}}(x)$, $\psi_{\text{per}}(x)$ the periodic Green matrix and pressure vector gathering the periodic fundamental solutions $(G_{\text{per},i}, \psi_{\text{per},i})_{1 \le i \le d}$:

$$ G_{\text{per}}(x) = [G_{\text{per},1}(x) \dots G_{\text{per},d}], \quad \psi_{\text{per}}(x) = (\psi_{\text{per},i}(x))_{1 \le i \le d}. \quad (2.50) $$

Similarly, we denote by $R(x)$ and $s(x)$ the matrix and the vector obtained by gathering the vector and pressure fields $(R_i(x))_{1 \le i \le d}$ and $(s_i(x))_{1 \le i \le d}$:

$$ R(x) = [R_1(x) \dots R_d], \quad s(x) = (s_i(x))_{1 \le i \le d}. \quad (2.51) $$

The constant $c$ can be characterized as the constant arising in the asymptotic expansion of $G_{\text{per}}(x)$ near zero:

$$ G_{\text{per}}(x) = G(x) + c\text{I} + O(|x|^2). \quad (2.52) $$

It is related to a similar constant involved in the asymptotic expansion near zero of the periodic Green function for the Laplace operator: there exists a universal constant $c$ depending only on the dimension such that

$$ \Gamma_{\text{per}}(x) = \Gamma(x) + c - \frac{1}{2d}|x|^2 + O(|x|^4) \text{ as } |x| \to 0, \text{ where } \Gamma_{\text{per}}(x) := -\frac{1}{4\pi^2} \sum_{k \in \mathbb{Z}^d \setminus \{0\}} \frac{e^{ik \cdot x}}{|k|^2}. \quad (2.53) $$

The constant $c$ can be computed accurately thanks to Evald's theta method described in [15, 6, 23]. Its numerical values are given in two and three dimensions by

$$ c \approx \begin{cases} 0.208578 & \text{if } d = 2, \\ 0.225791 & \text{if } d = 3. \end{cases} \quad (2.54) $$

The next proposition shows that the constant $c$ arising in (2.52) can be easily computed from the constant $c$. The proof yields also the value of $\Delta R(0)$, which plays a role in Proposition 4.3 below.

**Proposition 2.10.** The constants $c$ and $c$ arising respectively in (2.52) and (2.53) are related through the formula

$$ c = \left(1 - \frac{1}{d}\right) c + \frac{1}{8\pi} \delta_{d=2}. \quad (2.55) $$

Moreover,

$$ \Delta R(0) = - \left(1 - \frac{1}{d}\right) \text{I}. \quad (2.56) $$
**Proposition 2.9.** There exists a constant $c \in \mathbb{R}$ depending only on the dimension $d$ such that

$$G_{\text{per},i}(x) = G_i(x) + \mathbf{R}_i(x) + c\mathbf{e}_i, \quad \psi_{\text{per},i}(x) = \psi_i(x) + s_i(x), \quad \text{for all } x \in P \setminus \{0\}, \qquad (2.47)$$

for any $1 \le i \le d$, where $(\mathbf{R}_i, s_i) \in H^1(P, \mathbb{R}^d) \times L^2(P)$ is the solution to the periodic difference problem

$$\left\{
\begin{array}{l}
-\Delta \mathbf{R}_i + \nabla s_i = \mathbf{e}_i \text{ in } P, \\
\operatorname{div}(\mathbf{R}_i) = 0 \text{ in } P, \\
\mathbf{R}_i + \mathbf{G}_i \text{ is } P\text{-periodic}, \\
\sigma(\mathbf{R}_i, s_i) \mathbf{n}\mathbf{n}^T + \sigma(\mathbf{G}_i, \psi_i) \mathbf{n}\mathbf{n}^T \text{ is } P\text{-periodic}, \\
\mathbf{R}_i(0) = 0, \quad s_i(0) = 0.
\end{array}
\right.
\qquad (2.48)$$

The functions $\mathbf{R}_i$ and $s_i$ are respectively even and odd.

*Proof.* The compatibility condition $\int_P \mathbf{G}_i \cdot \mathbf{n} \, d\sigma = 0$ is satisfied due to $\operatorname{div}(\mathbf{G}_i) = 0$. Moreover,

$$\int_P \mathbf{e}_i \, dx + \int_{\partial P} \sigma(\mathbf{G}_i, \psi_i) \mathbf{n} \, d\sigma = 0$$

holds in virtue of Lemma 2.1. Therefore, Proposition 2.8 implies that (2.48) has a unique solution. Moreover, $\mathbf{R}_i + \mathbf{G}_i$ and $\psi_i + s_i$ can be extended by periodicity into a smooth function on $\mathbb{R}^d \setminus \mathbb{Z}^d$ and satisfy (2.45) by construction. By uniqueness, we deduce that there exists constant vectors $\mathbf{c}_i \in \mathbb{R}^d$ and constants $\kappa_i \in \mathbb{R}$ such that

$$\mathbf{G}_{\text{per},i}(x) = \mathbf{G}_i(x) + \mathbf{R}_i(x) + \mathbf{c}_i, \quad \psi_{\text{per},i}(x) = \psi_i(x) + s_i(x) + \kappa_i. \qquad (2.49)$$

Let us then consider the cubic cell symmetries $S^i := \text{I} - 2\mathbf{e}_i\mathbf{e}_i^T$ which changes the coordinate $x_i$ into $-x_i$, and the symmetry with respect to the diagonal axes $S^{ij} := \text{I} - \mathbf{e}_i\mathbf{e}_i^T - \mathbf{e}_j\mathbf{e}_j^T + \mathbf{e}_i\mathbf{e}_j^T + \mathbf{e}_j\mathbf{e}_i^T$ for $i \neq j$, which exchange $x_i$ and $x_j$ for $1 \le i \neq j \le d$. These symmetries preserve the hypercube $P$. By using the spherical symmetry properties of the Stokes fundamental solution $(\mathbf{G}_i, \psi_i)$, it is can be verified that

$$S\mathbf{R}_i \circ S = \sum_{l=1}^{d} S_{il}\mathbf{R}_l, \quad s_i \circ S = \sum_{l=1}^{d} S_{il}s_l,$$

see e.g. for a method of proof [18, Section 3.3]. In particular, with $S = S^i$, these properties imply

$$\int_P \mathbf{R}_i \cdot \mathbf{e}_j \, dx = 0 \text{ for all } 1 \le j \neq i \le d \text{ and } \int_P s_i \, dx = 0 \text{ for all } 1 \le i \le d.$$

Integrating (2.49) over the cell $P$, we deduce

$$\mathbf{c}_i \cdot \mathbf{e}_j = 0 \text{ for all } 1 \le j \neq i \le d \text{ and } \kappa_i = 0 \text{ for all } 1 \le i \le d.$$

This shows that $\mathbf{c}_i = c_i\mathbf{e}_i$ for some constant $c_i$. Considering now the symmetries $S^{ij}$, we find

$$\int_P \mathbf{R}_i \cdot \mathbf{e}_i \, dx = \int_P \mathbf{R}_j \cdot \mathbf{e}_j \, dx \text{ for all } 1 \le i, j \le d,$$

which imply that the constants $c_i$ are identically equal to a single constant $c$. The fact that $\kappa_i = 0$ for $1 \le i \le d$ can be obtained by observing that $\psi_{\text{per},i}$ and $\psi_i$ are odd functions and evaluating (2.49) at zero. $\square$

In the remainder of the paper, we will denote by $G_{\text{per}}(x)$, $\psi_{\text{per}}(x)$ the periodic Green matrix and pressure vector gathering the periodic fundamental solutions $(\mathbf{G}_{\text{per},i}, \psi_{\text{per},i})_{1 \le i \le d}$:

$$G_{\text{per}}(x) = [\mathbf{G}_{\text{per},1}(x) \quad \dots \quad \mathbf{G}_{\text{per},d}], \quad \psi_{\text{per}}(x) = (\psi_{\text{per},i}(x))_{1 \le i \le d}. \qquad (2.50)$$

Similarly, we denote by $R(x)$ and $s(x)$ the matrix and the vector obtained by gathering the vector and pressure fields $(\mathbf{R}_i(x))_{1 \le i \le d}$ and $(s_i(x))_{1 \le i \le d}$:

$$R(x) = [\mathbf{R}_1(x) \quad \dots \quad \mathbf{R}_d], \quad s(x) = (s_i(x))_{1 \le i \le d}. \qquad (2.51)$$

The constant $c$ can be characterized as the constant arising in the asymptotic expansion of $G_{\text{per}}(x)$ near zero:

$$G_{\text{per}}(x) = G(x) + c\text{I} + O(|x|^2). \qquad (2.52)$$

It is related to a similar constant involved in the asymptotic expansion near zero of the periodic Green function for the Laplace operator: there exists a universal constant $c$ depending only on the dimension such that

$$\Gamma_{\text{per}}(x) = \Gamma(x) + c - \frac{1}{2d}|x|^2 + O(|x|^4) \text{ as } |x| \to 0, \text{ where } \Gamma_{\text{per}}(x) := -\frac{1}{4\pi^2} \sum_{k \in \mathbb{Z}^d \setminus \{0\}} \frac{e^{ik \cdot x}}{|k|^2}. \qquad (2.53)$$
Proof. From (2.47), we have the following asymptotic expansion for the trace of $G_{\text{per}}$:

$$ \mathrm{Tr}(G_{\mathrm{per}}) = \mathrm{Tr}(G) + d\mathfrak{c} + \frac{1}{2}\nabla^2\mathrm{Tr}(R)(0) \cdot x^2 + O(|x|^4) \text{ as } |x| \to 0. \quad (2.57) $$

By using Euler's identity $x_i \partial_i \Gamma = (2-d)\Gamma + \frac{1}{2\pi}\delta_{d=2}$ and (2.4), we find on the one hand

$$ \mathrm{Tr}(G(x)) = \sum_{i=1}^{d} e_i^T G_i(x) = \frac{d}{2}\Gamma(x) - \frac{1}{2}\sum_{i=1}^{d} x_i \partial_i \Gamma(x) = \frac{d}{2}\Gamma(x) - \frac{2-d}{2}\Gamma(x) - \frac{1}{4\pi}\delta_{d=2} = (d-1)\Gamma(x) - \frac{1}{4\pi}\delta_{d=2}. $$

On the other hand, using (2.53),

$$ \begin{aligned} \mathrm{Tr}(G_{\mathrm{per}}(x)) &= -\frac{1}{4\pi^2} \sum_{k \in \mathbb{Z}^d \setminus \{0\}} \sum_{i=1}^{d} \left( \frac{1}{|k|^2} - \frac{k_i^2}{|k|^4} \right) e^{2i\pi k \cdot x} \\ &= -\frac{(d-1)}{4\pi^2} \sum_{k \in \mathbb{Z}^d \setminus \{0\}} \frac{e^{2i\pi k \cdot x}}{|k|^2} = (d-1) \left( \Gamma(x) + c - \frac{1}{2d} |x|^2 + O(|x|^4) \right) \quad \text{as } |x| \to 0. \end{aligned} \quad (2.58) $$

Equating (2.55) and (2.57) and identifying the constant terms, we obtain $(d-1)c = -\frac{1}{4\pi}\delta_{d=2} + d\mathfrak{c}$ which yields (2.55). Let us then observe that due to the cubic symmetry, only the coefficients $\partial_{ii}R_{jj}(0)$, $\partial_{ij}R_{ij}(0)$ and $\partial_{ii}R_{ii}(0)$ with $1 \le i \ne j \le d$ are not zero, and these do not depend on the chosen pair $\{i, j\}$. Consequently, for $1 \le i_1, i_2 \le d$,

$$ \frac{1}{2} \nabla^2 \mathrm{Tr}(R)(0) \cdot x^2 = \frac{1}{2} \sum_{1 \le i_1, i_2 \le d} \sum_{k=1}^{d} \partial_{i_1 i_2} R_{kk}(0) x_{i_1} x_{i_2} = \frac{1}{2} (\partial_{ii} R_{ii}(0) + (d-1) \partial_{jj} R_{jj}(0)) |x|^2. $$

By identifying the quadratic term in (2.57) and (2.58), we thus obtain:

$$ \frac{1}{2}(\partial_{ii}R_{ii}(0) + (d-1)\partial_{jj}R_{jj}(0)) = -\frac{d-1}{2d}. $$

Finally, using the cubic symmetries as above, we find that $\Delta R(0)$ is a diagonal matrix, whose components are

$$ \Delta R_{lm}(0) = \sum_{i_1=1}^{d} \partial_{i_1 i_1} R_{lm}(0) = \delta_{lm} (\partial_{ii} R_{ii}(0) + (d-1) \partial_{ii} R_{jj}(0)), \quad 1 \le l, m \le d. $$

The two above equalities imply (2.56). $\square$

From (2.54), we deduce the following numerical values for the universal constant $c$:

$$ c \approx \begin{cases} 0.1505273 & \text{if } d=3, \\ 0.1440777 & \text{if } d=2. \end{cases} $$

In the remainder of this part, we introduce two family of matrix valued tensors $(\mathfrak{G}^k)_{k \ge 0}$ and $\mathfrak{R}^k_{k \ge 0}$ that play a role in the derivation of the asymptotic expansions of the permeability matrix $\mathcal{X}^*$. For any $k \in \mathbb{N}$, we define the $k$-th order matrix valued tensor $\mathfrak{G}^k$ by the formula

$$ \mathfrak{G}^k := \int_{\partial P} \nabla^k G(x) \boldsymbol{x} \cdot \boldsymbol{n} \, d\sigma(x) = \frac{1}{2} \int_{\partial P} \nabla^k G(x) \, d\sigma(x). \quad (2.59) $$

The tensors $\mathfrak{R}^k$ are defined as the cell average of the derivatives of the smooth matrix $R$ of (2.51):

$$ \mathfrak{R}^k := \int_P \nabla^k R \, dx. \quad (2.60) $$

Due to the radial symmetry of $G$ and the cubic symmetry of $P$, many components of the tensor $\mathfrak{G}^k$ and $\mathfrak{R}^k$ turn out to be zero.

**Lemma 2.2.** The following properties hold for the matrix valued tensors $\mathfrak{G}^k$ and $\mathfrak{R}^k$ defined in (2.59) and (2.60):

(i) $\mathfrak{G}^k_{i_1 \dots i_k, lm} = 0$ and $\mathfrak{R}^k_{i_1 \dots i_k, lm} = 0$ whenever there exists an integer $j \in \{1, 2, \dots, k\}$ occurring an odd number of times in the indices $i_1, \dots, i_k, l, m$. In particular, the tensors $\mathfrak{G}^k$ and $\mathfrak{R}^k$ of odd order $k$ are zero: for any $k \in \mathbb{N}$, $\mathfrak{G}^{2k+1} = 0$.

(ii) For any permutation $\zeta \in \mathfrak{S}_d$ of $\{1, \dots, d\}$, it holds

$$ \mathfrak{G}_{\zeta(i_1)\dots\zeta(i_k), \zeta(l)\zeta(m)}^k = \mathfrak{G}_{i_1\dots i_k, lm}^k, \quad \mathfrak{R}_{\zeta(i_1)\dots\zeta(i_k), \zeta(l)\zeta(m)}^k = \mathfrak{R}_{i_1\dots i_k, lm}^k. $$

Proof. (i) and (ii). For any orthogonal symmetry $S$ satisfying $S^T = S = S^{-1}$, it holds $G \circ S = SGS$ (Lemma 2.1). By differentiating $k$ times and applying the chain rule, we find

$$ \sum_{1 \le j_1 \dots j_k \le d} S_{i_1 j_1} \dots S_{i_k j_k} \partial_{j_1 \dots j_k}^k G \circ S = S \partial_{i_1 \dots i_k}^k GS. $$
The constant c can be computed accurately thanks to Evald's theta method described in [16, 6, 24]. Its numerical values are given in two and three dimensions by

$$c \approx \begin{cases} 0.208578 & \text{if } d = 2, \\ 0.225791 & \text{if } d = 3. \end{cases} \qquad (2.54)$$

The next proposition shows that the constant c arising in (2.52) can be easily computed from the constant c. The proof yields also the value of $\Delta R(0)$, which plays a role in Proposition 4.3 below.

**Proposition 2.10.** The constants c and c arising respectively in (2.52) and (2.53) are related through the formula

$$c = \left(1 - \frac{1}{d}\right) c + \frac{1}{8\pi} \delta_{d=2}. \qquad (2.55)$$

Moreover,

$$\Delta R(0) = - \left(1 - \frac{1}{d}\right) I. \qquad (2.56)$$

*Proof.* From (2.47), we have the following asymptotic expansion for the trace of $G_{\text{per}}$:

$$\mathrm{Tr}(G_{\mathrm{per}}) = \mathrm{Tr}(G) + d\mathfrak{c} + \frac{1}{2}\nabla^2\mathrm{Tr}(R)(0) \cdot x^2 + O(|x|^4) \text{ as } |x| \to 0. \qquad (2.57)$$

By using Euler's identity $x_i \partial_i \Gamma = (2-d)\Gamma + \frac{1}{2\pi}\delta_{d=2}$ and (2.4), we find on the one hand

$$\mathrm{Tr}(G(x)) = \sum_{i=1}^{d} e_i^T G_i(x) = \frac{d}{2}\Gamma(x) - \frac{1}{2}\sum_{i=1}^{d} x_i \partial_i \Gamma(x) = \frac{d}{2}\Gamma(x) - \frac{2-d}{2}\Gamma(x) - \frac{1}{4\pi}\delta_{d=2} = (d-1)\Gamma(x) - \frac{1}{4\pi}\delta_{d=2}.$$

On the other hand, using (2.53),

$$\begin{aligned} \mathrm{Tr}(G_{\mathrm{per}}(x)) &= -\frac{1}{4\pi^2} \sum_{k \in \mathbb{Z}^d \setminus \{0\}} \sum_{i=1}^{d} \left( \frac{1}{|k|^2} - \frac{k_i^2}{|k|^4} \right) e^{2i\pi k \cdot x} \\ &= -\frac{(d-1)}{4\pi^2} \sum_{k \in \mathbb{Z}^d \setminus \{0\}} \frac{e^{2i\pi k \cdot x}}{|k|^2} = (d-1) \left( \Gamma(x) + c - \frac{1}{2d} |x|^2 + O(|x|^4) \right) \text{ as } |x| \to 0. \end{aligned} \qquad (2.58)$$

Equating (2.55) and (2.57) and identifying the constant terms, we obtain $(d-1)c = -\frac{1}{4\pi}\delta_{d=2} + d\mathfrak{c}$ which yields (2.55). Let us then observe that due to the cubic symmetry, only the coefficients $\partial_{ii}R_{jj}(0)$, $\partial_{ij}R_{ij}(0)$ and $\partial_{ii}R_{ii}(0)$ with $1 \le i \ne j \le d$ are not zero, and these do not depend on the chosen pair $\{i, j\}$. Consequently, for $1 \le i_1, i_2 \le d$,

$$\frac{1}{2} \nabla^2 \mathrm{Tr}(R)(0) \cdot x^2 = \frac{1}{2} \sum_{1 \le i_1, i_2 \le d} \sum_{k=1}^{d} \partial_{i_1 i_2} R_{kk}(0) x_{i_1} x_{i_2} = \frac{1}{2} (\partial_{ii} R_{ii}(0) + (d-1) \partial_{jj} R_{jj}(0)) |x|^2.$$

By identifying the quadratic term in (2.57) and (2.58), we thus obtain:

$$\frac{1}{2}(\partial_{ii}R_{ii}(0) + (d-1)\partial_{jj}R_{jj}(0)) = -\frac{d-1}{2d}.$$

Finally, using the cubic symmetries as above, we find that $\Delta R(0)$ is a diagonal matrix, whose components are

$$\Delta R_{lm}(0) = \sum_{i_1=1}^{d} \partial_{i_1 i_1} R_{lm}(0) = \delta_{lm}(\partial_{ii} R_{ii}(0) + (d-1)\partial_{ii} R_{jj}(0)), \quad 1 \le l, m \le d.$$

The two above equalities imply (2.56). □

From (2.54), we deduce the following numerical values for the universal constant c:

$$c \approx \begin{cases} 0.1505273 & \text{if } d = 3, \\ 0.1440777 & \text{if } d = 2. \end{cases}$$

In the remainder of this part, we introduce two family of matrix valued tensors $(\mathfrak{G}^k)_{k \ge 0}$ and $\mathfrak{R}^k_{k \ge 0}$ that play a role in the derivation of the asymptotic expansions of the permeability matrix $\mathcal{X}^*$. For any $k \in \mathbb{N}$, we define the k-th order matrix valued tensor $\mathfrak{G}^k$ by the formula

$$\mathfrak{G}^k := \int_{\partial P} \nabla^k G(x) \boldsymbol{x} \cdot \boldsymbol{n} \, d\sigma(x) = \frac{1}{2} \int_{\partial P} \nabla^k G(x) \, d\sigma(x). \qquad (2.59)$$

The tensors $\mathfrak{R}^k$ are defined as the cell average of the derivatives of the smooth matrix $R$ of (2.51):

$$\mathfrak{R}^k := \int_P \nabla^k R \, dx. \qquad (2.60)$$
Therefore, for any radial symmetry preserving P, performing the change of variable $x = Sx'$ in (2.59) yields

$$ \mathfrak{G}_{i_1 \dots i_k}^k = \frac{1}{2} \int_{\partial P} \partial_{i_1 \dots i_k}^k G \circ S \, d\sigma = \frac{1}{2} S_{i_1 j_1} \dots S_{i_k j_k} \int_{\partial P} S \partial_{j_1 \dots j_k}^k GS \, d\sigma = S_{i_1 j_1} \dots S_{i_k j_k} \mathfrak{G}_{i_1 \dots i_k}^k S, $$

which can be written as

$$ \mathfrak{G}_{i_1 \dots i_k, l_1 l_2}^k = S_{i_1 j_1} \dots S_{i_k j_k} S_{l_1 m_1} S_{l_2 m_2} \mathfrak{G}_{j_1 \dots j_k, m_1 m_2}, $$

with summation over the repeated indices $1 \le j_1 \dots j_k, m_1, m_2 \le d$. Points (i) and (ii) follow by setting $S$ to be the symmetry with respect to respectively the cell axes $e_i$ and the diagonal axes $e_i - e_j$ for $1 \le i \ne j \le d$. $\square$

Remarkably, $\mathfrak{G}^2 = 0$ and the tensors $\mathfrak{G}^k$ and $\mathfrak{R}^k$ are related by simple formulas.

**Proposition 2.11.** The following identities hold:

(i) $\mathfrak{R}^0 = -\frac{1}{2}\mathfrak{G}^0 - cI + \frac{1}{8\pi}\delta_{d=2}I$;

(ii) for any $p \in \mathbb{N} \setminus \{0\}$:

$$ \mathfrak{G}^p = (p-2)\mathfrak{R}^p. \qquad (2.61) $$

(iii) $\mathfrak{R}^2$ is also given by the formula

$$ \mathfrak{R}_{i_1 i_2}^2 = - \int_{\partial P} n_{i_2} \partial_{i_1} G \, d\sigma, \quad 1 \le i_1, i_2 \le d. \qquad (2.62) $$

*Proof.* (i) Noticing that $G_{\text{per}}$ and $G$ are integrable functions and $\int_P G_{\text{per}} dx = 0$, the average of $\mathbf{R}_i$ can be obtained by integrating (2.47) which yields

$$ \int_P \mathbf{R}_i \, dx = - \int_P \mathbf{G}_i \, dx - c e_i. $$

We then use Euler's identity $\text{div}(\mathbf{G}_i \mathbf{x}^T) = 2\mathbf{G}_i + \frac{\delta_{d=2}}{4\pi} e_i$ and the divergence theorem to obtain

$$ \int_P \mathbf{G}_i \, dx = \frac{1}{2} \int_P \left( \text{div}(\mathbf{G}_i \mathbf{x}^T) - \frac{\delta_{d=2}}{4\pi} e_i \right) \, dx = \frac{1}{2} \int_{\partial P} \mathbf{G}_i(x) \mathbf{x} \cdot \mathbf{n} \, d\sigma(x) - \frac{\delta_{d=2}}{8\pi} e_i, $$

which proves (i).

(ii) For $p=1$, the result could be obtained similarly, however the proof does not extend to $p \ge 2$ because the singularity of $\nabla^p G_{\text{per}}$ and $\nabla^p G$ at zero are then not integrable. Alternatively, for a given $p \ge 1$ and indices $1 \le i_1, \dots, i_p \le d$, we use the divergence theorem on $R$ to compute

$$ \mathfrak{R}_{i_1 \dots i_p}^p = \int_P \partial_{i_1 \dots i_p}^p R \, dx = \int_{\partial P} \partial_{i_2 \dots i_p}^{p-1} R n_{i_1} \, d\sigma = \int_{\partial P} \partial_{i_2 \dots i_p}^{p-1} (G_{\text{per}} - G) n_{i_1} \, d\sigma = - \int_{\partial P} \partial_{i_2 \dots i_p}^{p-1} G n_{i_1} \, d\sigma, \quad (2.63) $$

where we used $\int_{\partial P} \partial_{i_2 \dots i_p}^{p-1} G_{\text{per}} n_{i_1} \, d\sigma = 0$ due to the $P$-periodicity of $\nabla^{p-1} G_{\text{per}}$. Introducing a small arbitrary parameter $\epsilon > 0$ and using the divergence theorem on $P \setminus B(0, \epsilon)$, we can rewrite (2.63) as

$$ \mathfrak{R}_{i_1 \dots i_p}^p = - \int_{P \setminus B(0, \epsilon)} \partial_{i_1 \dots i_p}^p G \, dx - \int_{\partial B(0, \epsilon)} \partial_{i_2 \dots i_p}^{p-1} G n_{i_1} \, d\sigma. $$

Using now Euler's identity $\text{div}(\mathbf{x} \nabla^p G) = (2-p)\nabla^p G$ and a change of variable in the integral over $\partial B(0, \epsilon)$, we arrive at

$$ \begin{aligned} (p-2)\mathfrak{R}_{i_1 \dots i_p}^p &= \int_{\partial P} \partial_{i_1 \dots i_p}^p G \mathbf{x} \cdot \mathbf{n} \, d\sigma - \int_{\partial B(0, \epsilon)} \partial_{i_1 \dots i_p}^p G \mathbf{x} \cdot \mathbf{n} \, d\sigma - (p-2) \int_{\partial B(0, \epsilon)} \partial_{i_2 \dots i_p}^{p-1} G n_{i_1} \, d\sigma \\ &= \mathfrak{G}_{i_1 \dots i_p}^p - \epsilon^{2-p} \int_{\partial B(0, 1)} \partial_{i_1 \dots i_p}^p G \, d\sigma - (p-2)\epsilon^{2-p} \int_{\partial B(0, 1)} \partial_{i_2 \dots i_p}^{p-1} G n_{i_1} \, d\sigma. \end{aligned} \qquad (2.64) $$

If $p=1$, (2.61) follows by taking the limit $\epsilon \to 0$. If $p > 2$, then the sum of the two integrals on $\partial B(0, 1)$ must vanish because the left-hand side is a finite quantity. This also implies (2.61).

It remains to prove that $\mathfrak{G}^2 = 0$ for the proof to be complete. Due to the symmetry properties (i) and (ii), only the coefficients $\mathfrak{G}_{ii,jj}^2$, $\mathfrak{G}_{ij,ij}^2$ and $\mathfrak{G}_{ii,ii}^2$ for $1 \le i \ne j \le d$ are not zero, and each of these three coefficients do not depend on the chosen pair $\{i, j\}$. By using $\text{div}(G) = 0$ and $\text{Tr}(\Delta G) = \text{Tr}(\nabla \Psi) = 0$, we find, for fixed $i \ne j$:

$$ \begin{aligned} \mathfrak{G}_{ij,ij}^2 &= \frac{1}{d-1} \sum_{k \ne i} \mathfrak{G}_{ik,ik}^2 = \frac{1}{d-1} \left( \sum_{k=1}^d \mathfrak{G}_{ik,ik}^2 - \mathfrak{G}_{ii,ii}^2 \right) = \frac{1}{d-1} \left( \frac{1}{2} \int_{\partial P} \partial_i (\text{div}(\mathbf{G}_i)) \, d\sigma - \mathfrak{G}_{ii,ii}^2 \right) \\ &= -\frac{1}{d-1} \mathfrak{G}_{ii,ii}^2, \end{aligned} $$
Due to the radial symmetry of G and the cubic symmetry of P, many components of the tensor $\mathfrak{G}^k$ and $\mathfrak{R}^k$ turn out to be zero.

**Lemma 2.2.** The following properties hold for the matrix valued tensors $\mathfrak{G}^k$ and $\mathfrak{R}^k$ defined in (2.59) and (2.60):

(i) $\mathfrak{G}^k_{i_1...i_k,lm} = 0$ and $\mathfrak{R}^k_{i_1...i_k,lm} = 0$ whenever there exists an integer $j \in \{1, 2, ..., k\}$ occurring an odd number of times in the indices $i_1, ..., i_k, l, m$. In particular, the tensors $\mathfrak{G}^k$ and $\mathfrak{R}^k$ of odd order $k$ are zero: for any $k \in \mathbb{N}$, $\mathfrak{G}^{2k+1} = 0$.

(ii) For any permutation $\zeta \in \mathfrak{S}_d$ of $\{1, ..., d\}$, it holds

$$ \mathfrak{G}^k_{\zeta(i_1)...\zeta(i_k),\zeta(l)\zeta(m)} = \mathfrak{G}^k_{i_1...i_k,lm}, \quad \mathfrak{R}^k_{\zeta(i_1)...\zeta(i_k),\zeta(l)\zeta(m)} = \mathfrak{R}^k_{i_1...i_k,lm}. $$

*Proof.* (i) and (ii). For any orthogonal symmetry $S$ satisfying $S^T = S = S^{-1}$, it holds $G \circ S = SGS$ (Lemma 2.1). By differentiating $k$ times and applying the chain rule, we find

$$ \sum_{1 \le j_1...j_k \le d} S_{i_1j_1}...S_{i_kj_k} \partial_{j_1...j_k}^k G \circ S = S \partial_{i_1...i_k}^k GS. $$

Therefore, for any radial symmetry preserving $P$, performing the change of variable $x = Sx'$ in (2.59) yields

$$ \mathfrak{G}^k_{i_1...i_k} = \frac{1}{2} \int_{\partial P} \partial_{i_1...i_k}^k G \circ S \, d\sigma = \frac{1}{2} S_{i_1j_1}...S_{i_kj_k} \int_{\partial P} S \partial_{j_1...j_k}^k GS \, d\sigma = S_{i_1j_1}...S_{i_kj_k} \mathfrak{G}^k_{i_1...i_k} S, $$

which can be written as

$$ \mathfrak{G}^k_{i_1...i_k,l_1l_2} = S_{i_1j_1}...S_{i_kj_k}S_{l_1m_1}S_{l_2m_2}\mathfrak{G}_{j_1...j_k,m_1m_2}, $$

with summation over the repeated indices $1 \le j_1...j_k, m_1, m_2 \le d$. Points (i) and (ii) follow by setting $S$ to be the symmetry with respect to respectively the cell axes $e_i$ and the diagonal axes $e_i - e_j$ for $1 \le i \ne j \le d$. $\square$

Remarkably, $\mathfrak{G}^2 = 0$ and the tensors $\mathfrak{G}^k$ and $\mathfrak{R}^k$ are related by simple formulas.

**Proposition 2.11.** The following identities hold:

(i) $\mathfrak{R}^0 = -\frac{1}{2}\mathfrak{G}^0 - cI + \frac{1}{8\pi}\delta_{d=2}I$;

(ii) for any $p \in \mathbb{N} \setminus \{0\}$:

$$ \mathfrak{G}^p = (p-2)\mathfrak{R}^p. \qquad (2.61) $$

(iii) $\mathfrak{R}^2$ is also given by the formula

$$ \mathfrak{R}_{i_1 i_2}^2 = - \int_{\partial P} n_{i_2} \partial_{i_1} G \, d\sigma, \quad 1 \le i_1, i_2 \le d. \qquad (2.62) $$

*Proof.* (i) Noticing that $G_{\text{per}}$ and $G$ are integrable functions and $\int_P G_{\text{per}} dx = 0$, the average of $\mathbf{R}_i$ can be obtained by integrating (2.47) which yields

$$ \int_P \mathbf{R}_i \, dx = - \int_P \mathbf{G}_i \, dx - c \mathbf{e}_i. $$

We then use Euler's identity $\text{div}(\mathbf{G}_i x^T) = 2\mathbf{G}_i + \frac{\delta_{d=2}}{4\pi} \mathbf{e}_i$ and the divergence theorem to obtain

$$ \int_P \mathbf{G}_i \, dx = \frac{1}{2} \int_P \left( \text{div}(\mathbf{G}_i x^T) - \frac{\delta_{d=2}}{4\pi} \mathbf{e}_i \right) \, dx = \frac{1}{2} \int_{\partial P} \mathbf{G}_i(x) \, x \cdot \mathbf{n} \, d\sigma(x) - \frac{\delta_{d=2}}{8\pi} \mathbf{e}_i, $$

which proves (i).

(ii) For $p = 1$, the result could be obtained similarly, however the proof does not extend to $p \ge 2$ because the singularity of $\nabla^p G_{\text{per}}$ and $\nabla^p G$ at zero are then not integrable. Alternatively, for a given $p \ge 1$ and indices $1 \le i_1, ..., i_p \le d$, we use the divergence theorem on $R$ to compute

$$ \mathfrak{R}_{i_1...i_p}^p = \int_P \partial_{i_1...i_p}^p R \, dx = \int_{\partial P} \partial_{i_2...i_p}^{p-1} R n_{i_1} \, d\sigma = \int_{\partial P} \partial_{i_2...i_p}^{p-1} (G_{\text{per}} - G) n_{i_1} \, d\sigma = - \int_{\partial P} \partial_{i_2...i_p}^{p-1} G n_{i_1} \, d\sigma, \quad (2.63) $$

where we used $\int_{\partial P} \partial_{i_2...i_p}^{p-1} G_{\text{per}} n_{i_1} \, d\sigma = 0$ due to the $P$-periodicity of $\nabla^{p-1} G_{\text{per}}$. Introducing a small arbitrary parameter $\epsilon > 0$ and using the divergence theorem on $P \setminus B(0, \epsilon)$, we can rewrite (2.63) as

$$ \mathfrak{R}_{i_1...i_p}^p = - \int_{P \setminus B(0, \epsilon)} \partial_{i_1...i_p}^p G \, dx - \int_{\partial B(0, \epsilon)} \partial_{i_2...i_p}^{p-1} G n_{i_1} \, d\sigma. $$

Using now Euler's identity $\text{div}(x \nabla^p G) = (2-p)\nabla^p G$ and a change of variable in the integral over $\partial B(0, \epsilon)$, we arrive at

$$ \begin{aligned} (p-2)\mathfrak{R}_{i_1...i_p}^p &= \int_{\partial P} \partial_{i_1...i_p}^p G x \cdot n \, d\sigma - \int_{\partial B(0, \epsilon)} \partial_{i_1...i_p}^p G x \cdot n \, d\sigma - (p-2) \int_{\partial B(0, \epsilon)} \partial_{i_2...i_p}^{p-1} G n_{i_1} \, d\sigma \\ &= \mathfrak{G}_{i_1...i_p}^p - \epsilon^{2-p} \int_{\partial B(0, 1)} \partial_{i_1...i_p}^p G \, d\sigma - (p-2)\epsilon^{2-p} \int_{\partial B(0, 1)} \partial_{i_2...i_p}^{p-1} G n_{i_1} \, d\sigma. \end{aligned} \qquad (2.64) $$
$$
\begin{align*}
\mathfrak{G}_{ii,jj}^2 &= \frac{1}{d(d-1)} \sum_{k=1}^{d} \sum_{l \neq k} \mathfrak{G}_{kk,ll}^2 = \frac{1}{d(d-1)} \sum_{k=1}^{d} \left( \sum_{l=1}^{d} \mathfrak{G}_{kk,ll}^2 - \mathfrak{G}_{kk,kk}^2 \right) = \frac{1}{d(d-1)} \left( \frac{1}{2} \int_{\partial P} \mathrm{Tr}(\Delta G) \, d\sigma - d\mathfrak{G}_{ii,ii}^2 \right) \\
&= -\frac{1}{d-1} \mathfrak{G}_{ii,ii}^2.
\end{align*}
$$

Thus, to show that $\mathfrak{G}^2$ is zero, it is sufficient to show that the coefficient $\mathfrak{G}_{ii,ii}^2$ vanishes for a fixed $1 \le i \le d$. By using the expression (2.4), we find that $\partial_{ii}^2 G_{ii}(x) = -\frac{1}{2}\partial_{ii}\Gamma(x) - \frac{1}{2}x_i \partial_{iii}^3 \Gamma(x)$. Hence by using $\Delta\Gamma = 0$ on $\partial P$ and computing explicitly the third order derivative of $\Gamma$, we infer

$$
\mathfrak{G}_{ii,ii}^2 = \frac{1}{d} \sum_{l=1}^{d} \mathfrak{G}_{ll,ll}^2 = -\frac{1}{4d} \sum_{k=1}^{d} \int_{\partial P} x_i \partial_{ll}^3 \Gamma(x) \, d\sigma(x) = -\frac{1}{4|\partial B(0,1)|} \int_{\partial P} \left( (d+2) \frac{\sum_{l=1}^{d} x_l^4}{|x|^{d+4}} - \frac{3}{|x|^d} \right) \, d\sigma(x).
$$

This integral turns out to be zero due to the identity (B.2) with $k=1$ proved in Appendix B. This implies $\mathfrak{G}_{ii,ii}^2 = 0$ for any $1 \le i \le d$, which concludes the proof that $\mathfrak{G}^2 = 0$. We note that this result could also be obtained from (2.64) for $p=2$: since $\mathfrak{G}^2 = \int_{\partial B(0,1)} \nabla^2 G \, d\sigma$, the same computation can be done on the unit sphere, resorting to the identities of Lemma C.1.

(iii) is (2.63) with $p=2$.

<figure><img src="image_3.png" /></figure>

3. FULL ASYMPTOTIC EXPANSIONS OF THE SOLUTIONS TO THE PERFORATED PERIODIC STOKES SYSTEM

In this section, we derive full asymptotic expansions of the cell solutions ($\mathcal{X}_j, \alpha_j$) to the periodic cell problem (1.3) which are needed to deduce the one of the permeability matrix $\mathcal{X}^*$ in the next Section 4. Our methodology is similar to that of [18] introduced for dealing with the scalar model (1.10) and declines in the following steps:

(i) we find a layer potential representation which allows to prove some analiticity results for the map $\eta \mapsto (\mathcal{X}_j, \alpha_j)$, and hence the existence of an asymptotic expansion with respect to $\eta \to 0^+$;

(ii) the coefficients of the expansion are identified by inserting it as an ansatz in the original cell problem (1.3): this results in a well-defined procedure for computing them by solving a cascade of interior and exterior problems;

(iii) we prove $L^2$ error estimates for the approximation of $(\mathcal{X}_j, \alpha_j)$ by this asymptotic expansion truncated at arbitrary orders.

For an arbitrary constant matrix $\kappa \in \mathbb{R}^{d \times d}$, we introduce the periodic single layer potential $S_{\text{per},\eta T}^{\kappa}$ defined by

$$
S_{\text{per},\eta T}^{\kappa}[\phi](x) := \int_{\partial(\eta T)} G_{\text{per}}(x-x')\phi(x') \, d\sigma(x') + \eta^{2-d} \left(-\frac{\log \eta}{4\pi} \delta_{d=2} \mathbf{I} - \mathbf{c}\mathbf{I} + \kappa\right) \int_{\partial(\eta T)} \phi \, d\sigma, \quad (3.1)
$$

$$
\mathcal{P}_{\text{per},\eta T}[\phi](x) := \int_{\partial(\eta T)} \psi_{\text{per}}(x-x') \cdot \phi(x') \, d\sigma(x'), \quad (3.2)
$$

for any $\phi \in H^{-\frac{1}{2}}(\partial(\eta T), \mathbb{R}^d)$ and $x \in P$. $S_{\text{per},\eta T}^{\kappa}$ is the standard periodic single layer potentials with the kernel $G_{\text{per}}$ (defined in (2.50)) augmented by a constant matrix depending on $\kappa$ and $\eta$. Since $G_{\text{per}}$ and $\psi_{\text{per}}$ are smooth perturbations of $G$ and $\psi$, these potentials satisfy jump relations (2.11) to (2.12) similar to those of the “exterior” potentials $S_T$ and $P_T$. Moreover, it still holds $S_{\text{per},\eta T}^{\kappa}[\mathbf{n}] = 0$ due to $\text{div}(G_{\text{per}}) = 0$. The invariance of the pressure by adding a potential function $\phi$ proportional to the normal $\mathbf{n}$ can be summarized as follows.

$$
\textbf{Lemma 3.1.} \quad (i) \quad S_{\text{per},\eta T}^{\kappa}[\mathbf{n}] = 0 \\
\text{(ii)} \quad \mathcal{P}_{\text{per},\eta T}[\mathbf{n}](y) = \begin{cases} -1 + \eta^d |T| & \text{if } y \in (\eta T), \\ -\frac{1}{2} + \eta^d |T| & \text{if } y \in \partial(\eta T), \\ \eta^d |T| & \text{if } y \in P \setminus (\eta T). \end{cases}
$$

Proof. (i) is a simple consequence of $\text{div}(G_{\text{per}}) = 0$

(ii) can be obtained by writing

$$
\mathcal{P}_{\text{per},\eta T}[\mathbf{n}](x-x') = -\int_{\partial(\eta T)} \nabla_x \Gamma(x-x') \cdot \mathbf{n}(x') \, d\sigma(x') - \int_{\partial(\eta T)} \nabla_x r(x-x') \cdot \mathbf{n}(x') \, d\sigma(x')
$$

where $\Gamma$ is the fundamental solution of the Laplace operator and $r$ is a smooth function in $P$ satisfying $-\Delta r = 1$. The first term is simply $-\mathcal{D}_{\partial(\eta T)}[1]$ where $\mathcal{D}_{\partial(\eta T)}$ is the Laplace double layer potential, which has values $-1, -1/2$ and $0$ respectively outside, on the boundary and inside $\eta T$. Applying Stokes' theorem, we find that the second term is equal to

$$
- \int_{\eta T} \Delta r(x - x') \, dx' = \eta^d |T|.
$$

The result follows.
If $p = 1$, (2.61) follows by taking the limit $\epsilon \to 0$. If $p > 2$, then the sum of the two integrals on $\partial B(0, 1)$ must vanish because the left-hand side is a finite quantity. This also implies (2.61).

It remains to prove that $\mathfrak{G}^2 = 0$ for the proof to be complete. Due to the symmetry properties (i) and (ii), only the coefficients $\mathfrak{G}_{ii,jj}^2$, $\mathfrak{G}_{ij,ii}^2$ and $\mathfrak{G}_{ii,ii}^2$ for $1 \le i \ne j \le d$ are not zero, and each of these three coefficients do not depend on the chosen pair $\{i, j\}$. By using $\mathrm{div}(G) = 0$ and $\mathrm{Tr}(\Delta G) = \mathrm{Tr}(\nabla \Psi) = 0$, we find, for fixed $i \ne j$:

$$ \begin{aligned} \mathfrak{G}_{ij,ij}^2 &= \frac{1}{d-1} \sum_{k \ne i} \mathfrak{G}_{ik,ik}^2 = \frac{1}{d-1} \left( \sum_{k=1}^d \mathfrak{G}_{ik,ik}^2 - \mathfrak{G}_{ii,ii}^2 \right) = \frac{1}{d-1} \left( \frac{1}{2} \int_{\partial P} \partial_i (\mathrm{div}(G_i)) \, d\sigma - \mathfrak{G}_{ii,ii}^2 \right) \\ &= -\frac{1}{d-1} \mathfrak{G}_{ii,ii}^2, \end{aligned} $$

$$ \begin{aligned} \mathfrak{G}_{ii,jj}^2 &= \frac{1}{d(d-1)} \sum_{k=1}^{d} \sum_{l \ne k} \mathfrak{G}_{kk,ll}^2 = \frac{1}{d(d-1)} \sum_{k=1}^{d} \left( \sum_{l=1}^{d} \mathfrak{G}_{kk,ll}^2 - \mathfrak{G}_{kk,kk}^2 \right) = \frac{1}{d(d-1)} \left( \frac{1}{2} \int_{\partial P} \mathrm{Tr}(\Delta G) \, d\sigma - d\mathfrak{G}_{ii,ii}^2 \right) \\ &= -\frac{1}{d-1} \mathfrak{G}_{ii,ii}^2. \end{aligned} $$

Thus, to show that $\mathfrak{G}^2$ is zero, it is sufficient to show that the coefficient $\mathfrak{G}_{ii,ii}^2$ vanishes for a fixed $1 \le i \le d$. By using the expression (2.4), we find that $\partial_{ii}^2 G_{ii}(x) = -\frac{1}{2} \partial_{ii}\Gamma(x) - \frac{1}{2} x_i \partial_{iii}^3 \Gamma(x)$. Hence by using $\Delta\Gamma = 0$ on $\partial P$ and computing explicitly the third order derivative of $\Gamma$, we infer

$$ \mathfrak{G}_{ii,ii}^2 = \frac{1}{d} \sum_{l=1}^{d} \mathfrak{G}_{ll,ll}^2 = -\frac{1}{4d} \sum_{k=1}^{d} \int_{\partial P} x_l \partial_{ll}^3 \Gamma(x) \, d\sigma(x) = -\frac{1}{4|\partial B(0,1)|} \int_{\partial P} \left( (d+2) \frac{\sum_{l=1}^{d} x_l^4}{|x|^{d+4}} - \frac{3}{|x|^d} \right) \, d\sigma(x). $$

This integral turns out to be zero due to the identity (B.2) with $k = 1$ proved in Appendix B. This implies $\mathfrak{G}_{ii,ii}^2 = 0$ for any $1 \le i \le d$, which concludes the proof that $\mathfrak{G}^2 = 0$. We note that this result could also be obtained from (2.64) for $p = 2$: since $\mathfrak{G}^2 = \int_{\partial B(0,1)} \nabla^2 G \, d\sigma$, the same computation can be done on the unit sphere, resorting to the identities of Lemma C.1.

(iii) is (2.63) with $p = 2$.

□

### 3. FULL ASYMPTOTIC EXPANSIONS OF THE SOLUTIONS TO THE PERFORATED PERIODIC STOKES SYSTEM

In this section, we derive full asymptotic expansions of the cell solutions $(\mathcal{X}_j, \alpha_j)$ to the periodic cell problem (1.3) which are needed to deduce the one of the permeability matrix $\mathcal{X}^*$ in the next Section 4. Our methodology is similar to that of [19] introduced for dealing with the scalar model (1.10) and declines in the following steps:

(i) we find a layer potential representation which allows to prove some analiticity results for the map $\eta \mapsto (\mathcal{X}_j, \alpha_j)$, and hence the existence of an asymptotic expansion with respect to $\eta \to 0^+$;

(ii) the coefficients of the expansion are identified by inserting it as an ansatz in the original cell problem (1.3): this results in a well-defined procedure for computing them by solving a cascade of interior and exterior problems;

(iii) we prove $L^2$ error estimates for the approximation of $(\mathcal{X}_j, \alpha_j)$ by this asymptotic expansion truncated at arbitrary orders.

For an arbitrary constant matrix $\kappa \in \mathbb{R}^{d \times d}$, we introduce the periodic single layer potential $S_{\mathrm{per},\eta T}^{\kappa}$ defined by

$$ S_{\mathrm{per},\eta T}^{\kappa}[\phi](x) := \int_{\partial(\eta T)} G_{\mathrm{per}}(x-x')\phi(x') \, d\sigma(x') + \eta^{2-d} \left( -\frac{\log \eta}{4\pi} \delta_{d=2} \mathbf{I} - \mathbf{c}\mathbf{I} + \kappa \right) \int_{\partial(\eta T)} \phi \, d\sigma, \quad (3.1) $$

$$ \mathcal{P}_{\mathrm{per},\eta T}[\phi](x) := \int_{\partial(\eta T)} \psi_{\mathrm{per}}(x-x') \cdot \phi(x') \, d\sigma(x'), \quad (3.2) $$

for any $\phi \in H^{-\frac{1}{2}}(\partial(\eta T), \mathbb{R}^d)$ and $x \in P$. $S_{\mathrm{per},\eta T}^{\kappa}$ is the standard periodic single layer potentials with the kernel $G_{\mathrm{per}}$ (defined in (2.50)) augmented by a constant matrix depending on $\kappa$ and $\eta$. Since $G_{\mathrm{per}}$ and $\psi_{\mathrm{per}}$ are smooth perturbations of $G$ and $\psi$, these potentials satisfy jump relations (2.11) to (2.12) similar to those of the “exterior” potentials $S_T$ and $P_T$. Moreover, it still holds $S_{\mathrm{per},\eta T}^{\kappa}[\mathbf{n}] = 0$ due to $\mathrm{div}(G_{\mathrm{per}}) = 0$. The invariance of the pressure by adding a potential function $\phi$ proportional to the normal $\mathbf{n}$ can be summarized as follows.

$$ \textbf{Lemma 3.1.} \quad (i) S_{\mathrm{per},\eta T}^{\kappa}[\mathbf{n}] = 0 \\ \qquad (ii) \mathcal{P}_{\mathrm{per},\eta T}[\mathbf{n}](y) = \begin{cases} -1 + \eta^d |T| & \text{if } y \in (\eta T), \\ -\frac{1}{2} + \eta^d |T| & \text{if } y \in \partial(\eta T), \\ \eta^d |T| & \text{if } y \in P \setminus (\eta T). \end{cases} $$
The dependency in $\eta$ and $\log \eta$ in (3.1) is motivated by simplifications occurring in the integrals when performing the change of variable $x' = \eta y'$. More precisely, let us introduce the rescaling function $\tau_\eta$ defined by

$$ \tau_\eta(y) := \eta y \text{ for all } y \in \partial T, $$

and the rescaling operator $\mathcal{T}_\eta : H^s(\partial(\eta T), \mathbb{R}^d) \to H^s(\partial T, \mathbb{R}^d)$ defined by

$$ \mathcal{T}_\eta[\phi] := \phi \circ \tau_\eta, \text{ for any } \phi \in H^s(\partial(\eta T), \mathbb{R}^d), \forall s \in \mathbb{R}. $$

The rescaling operator $\mathcal{T}_\eta$ enables to factorize $S_{\text{per},\eta T}^\kappa$ and $\mathcal{P}_{\text{per},\eta T}$ in terms of operators analytic in $\eta$.

**Proposition 3.1.** The following factorization holds:

$$ S_{\text{per},\eta T}^{\kappa} = \eta \mathcal{T}_{\eta}^{-1} S_{T}^{\kappa}(\eta) \mathcal{T}_{\eta}, \qquad (3.3) $$

where $S_T^\kappa : H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d) \to H^{\frac{1}{2}}(\partial T, \mathbb{R}^d)$ is the operator given by

$$ S_T^\kappa(\eta)[\phi](y) = S_T[\phi](y) + \kappa \int_{\partial T} \phi \, d\sigma + \eta^{d-2} \int_{\partial T} R(\eta(y-y'))\phi(y') \, d\sigma(y'), \quad \forall y \in \partial T. \qquad (3.4) $$

*Proof.* We have, for any $\phi \in H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d)$, and $y \in \partial T$:

$$ \int_{\partial(\eta T)} G(\eta y - x') \phi(x'/\eta) \, d\sigma(x') = \eta^{d-1} \left( \eta^{2-d} \int_{\partial T} G(y - y') \phi(y') \, d\sigma(y') + \frac{1}{4\pi} \log \eta \delta_{d=2} \int_{\partial T} \phi(y') \, d\sigma(y') \right) \\ = \eta \left( S_T[\phi](y) + \frac{\delta_{d=2}}{4\pi} \log \eta \int_{\partial T} \phi(y') \, d\sigma(y') \right), $$

and

$$ \int_{\partial(\eta T)} R(\eta y - x') \phi(x'/\eta) \, d\sigma(x') = \eta^{d-1} \int_{\partial T} R(\eta(y - y')) \phi(y') \, d\sigma(y'). $$

The factorization (3.3) follows.

□

For $x \in P\setminus(\eta T)$ and $\phi \in H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d)$, the following identities will be useful later on:

$$ S_{\text{per},\eta T}^{\kappa}[\phi \circ \tau_{\eta}](x) = \eta \left( S_T[\phi] \left( \frac{x}{\eta} \right) + \kappa \int_{\partial T} \phi \, d\sigma + \eta^{d-2} \int_{\partial T} R(x - \eta y') \phi(y') \, d\sigma(y') \right), \quad (3.5) $$

$$ \mathcal{P}_{\text{per},\eta T}[\phi \circ \tau_{\eta}](x) = \mathcal{P}_T[\phi] \left(\frac{x}{\eta}\right) + \eta^{d-1} \int_{\partial T} \mathbf{s}(x - \eta y') \cdot \phi(y') \, d\sigma(y'). \quad (3.6) $$

In what follows, we consider the following extension of the solution to the problem (1.3) to the whole cell $P$, that we still denote by $(\mathbf{X}_i, \alpha_i)$:

$$ \left\{ \begin{array}{ll} -\Delta \mathbf{X}_i + \nabla \alpha_i = \mathbf{e}_i & \text{in } P\setminus\partial(\eta T), \\ \text{div}(\mathbf{X}_i) = 0 & \text{in } P\setminus\partial(\eta T), \\ \mathbf{X}_i = 0 & \text{on } \partial(\eta T), \\ (\mathbf{X}_i, \alpha_i) & \text{is } P\text{-periodic,} \end{array} \right. \quad 1 \le j \le d. \qquad (3.7) $$

Note that $\mathbf{X}_i$ is not extended by zero inside $\eta T$, and that $\alpha_i$ is defined up to some constant inside $\eta T$. This definition is convenient for characterizing $\mathbf{X}_i$ from the kernel of $S_{\text{per},\eta T}^\kappa$.

**Proposition 3.2.** The following inclusion holds for the kernel of $S_{\text{per},\eta T}^\kappa$:

$$ \ker(S_{\text{per},\eta T}^{\kappa}) \subset \operatorname{span}(\mathbf{n}) \oplus \operatorname{span}(\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket)_{1 \le i \le d}. \qquad (3.8) $$

Moreover, if $\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket \in \ker(S_{\text{per},\eta T}^\kappa)$, then $(\mathbf{X}_i, \alpha_i)$ admits the following single layer potential representation:

$$ \mathbf{X}_i(x) = S_{\text{per},\eta T}^{\kappa} \llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket (x), \quad \alpha_i(x) = \mathcal{P}_{\text{per},\eta T} \llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket (x) + C \text{ for } x \in P\setminus(\eta T), \qquad (3.9) $$

for some additive pressure constant $C \in \mathbb{R}$.

*Proof.* Let $\phi \in H^{-\frac{1}{2}}(\partial(\eta T), \mathbb{R}^d)$ be such that $S_{\text{per},\eta T}^\kappa[\phi] = 0$ on $\partial(\eta T)$. The functions $\mathbf{u} := S_{\text{per},\eta T}^\kappa[\phi]$ and $p := \mathcal{P}_{\text{per},\eta T}[\phi]$ solve the following Stokes problem

$$ \left\{ \begin{array}{ll} -\Delta \mathbf{u} + \nabla p = \displaystyle\int_{\partial(\eta T)} \phi \, d\sigma \quad \text{in } P\setminus\partial(\eta T), \\ \text{div}(\mathbf{u}) = 0 \quad \text{in } P\setminus\partial(\eta T), \\ \mathbf{u} = 0 \quad \text{on } \partial(\eta T), \\ \mathbf{u} \text{ is } P\text{-periodic,} \\ \sigma(\mathbf{u}, p)\mathbf{n} \text{ is } P\text{-periodic.} \end{array} \right. \qquad (3.10) $$
Proof. (i) is a simple consequence of $\mathrm{div}(G_{\mathrm{per}}) = 0$

(ii) can be obtained by writing

$$ \mathcal{P}_{\mathrm{per}, \eta T}[\mathbf{n}](x-x') = - \int_{\partial(\eta T)} \nabla_x \Gamma(x-x') \cdot \mathbf{n}(x') \mathrm{d}\sigma(x') - \int_{\partial(\eta T)} \nabla_x r(x-x') \cdot \mathbf{n}(x') \mathrm{d}\sigma(x') $$

where $\Gamma$ is the fundamental solution of the Laplace operator and $r$ is a smooth function in $P$ satisfying $-\Delta r = 1$. The first term is simply $-\mathcal{D}_{\partial(\eta T)}[1]$ where $\mathcal{D}_{\partial(\eta T)}$ is the Laplace double layer potential, which has values $-1, -1/2$ and $0$ respectively outside, on the boundary and inside $\eta T$. Applying Stokes' theorem, we find that the second term is equal to

$$ - \int_{\eta T} \Delta r(x - x') \mathrm{d}x' = \eta^d |T|. $$

The result follows.

<figure><img src="image_4.png" /></figure>

The dependency in $\eta$ and $\log \eta$ in (3.1) is motivated by simplifications occurring in the integrals when performing the change of variable $x' = \eta y'$. More precisely, let us introduce the rescaling function $\tau_{\eta}$ defined by

$$ \tau_{\eta}(y) := \eta y \quad \text{for all } y \in \partial T, $$

and the rescaling operator $\mathcal{T}_{\eta} : H^s(\partial(\eta T), \mathbb{R}^d) \to H^s(\partial T, \mathbb{R}^d)$ defined by

$$ \mathcal{T}_{\eta}[\phi] := \phi \circ \tau_{\eta}, \text{ for any } \phi \in H^s(\partial(\eta T), \mathbb{R}^d), \forall s \in \mathbb{R}. $$

The rescaling operator $\mathcal{T}_{\eta}$ enables to factorize $S_{\mathrm{per}, \eta T}^{\kappa}$ and $\mathcal{P}_{\mathrm{per}, \eta T}$ in terms of operators analytic in $\eta$.

**Proposition 3.1.** *The following factorization holds:*

$$ S_{\mathrm{per}, \eta T}^{\kappa} = \eta \mathcal{T}_{\eta}^{-1} S_{T}^{\kappa}(\eta) \mathcal{T}_{\eta}, \qquad (3.3) $$

where $S_T^\kappa : H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d) \to H^{\frac{1}{2}}(\partial T, \mathbb{R}^d)$ is the operator given by

$$ S_T^\kappa(\eta)[\phi](y) = S_T[\phi](y) + \kappa \int_{\partial T} \phi \mathrm{d}\sigma + \eta^{d-2} \int_{\partial T} R(\eta(y-y'))\phi(y') \mathrm{d}\sigma(y'), \quad \forall y \in \partial T. \qquad (3.4) $$

*Proof.* We have, for any $\phi \in H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d)$, and $y \in \partial T$:

$$ \int_{\partial(\eta T)} G(\eta y - x') \phi(x'/\eta) \mathrm{d}\sigma(x') = \eta^{d-1} \left( \eta^{2-d} \int_{\partial T} G(y-y') \phi(y') \mathrm{d}\sigma(y') + \frac{1}{4\pi} \log \eta \delta_{d=2} \int_{\partial T} \phi(y') \mathrm{d}\sigma(y') \right) \\ = \eta \left( S_T[\phi](y) + \frac{\delta_{d=2}}{4\pi} \log \eta \int_{\partial T} \phi(y') \mathrm{d}\sigma(y') \right), $$

and

$$ \int_{\partial(\eta T)} R(\eta y - x') \phi(x'/\eta) \mathrm{d}\sigma(x') = \eta^{d-1} \int_{\partial T} R(\eta(y - y')) \phi(y') \mathrm{d}\sigma(y'). $$

The factorization (3.3) follows.

<figure><img src="image_5.png" /></figure>

For $x \in P \setminus (\eta T)$ and $\phi \in H^{-\frac{1}{2}}(\partial T, \mathbb{R}^d)$, the following identities will be useful later on:

$$ S_{\mathrm{per}, \eta T}^{\kappa}[\phi \circ \tau_{\eta}](x) = \eta \left( S_T[\phi] \left( \frac{x}{\eta} \right) + \kappa \int_{\partial T} \phi \mathrm{d}\sigma + \eta^{d-2} \int_{\partial T} R(x - \eta y') \phi(y') \mathrm{d}\sigma(y') \right), \qquad (3.5) $$

$$ \mathcal{P}_{\mathrm{per}, \eta T}[\phi \circ \tau_{\eta}](x) = \mathcal{P}_T[\phi] \left(\frac{x}{\eta}\right) + \eta^{d-1} \int_{\partial T} \mathbf{s}(x - \eta y') \cdot \phi(y') \mathrm{d}\sigma(y'). \qquad (3.6) $$

In what follows, we consider the following extension of the solution to the problem (1.3) to the whole cell $P$, that we still denote by $(\mathbf{X}_i, \alpha_i)$:

$$ \left\{ \begin{array}{l} -\Delta \mathbf{X}_i + \nabla \alpha_i = \mathbf{e}_i \text{ in } P \setminus \partial(\eta T), \\ \mathrm{div}(\mathbf{X}_i) = 0 \text{ in } P \setminus \partial(\eta T), \\ \mathbf{X}_i = 0 \text{ on } \partial(\eta T), \\ (\mathbf{X}_i, \alpha_i) \text{ is } P\text{-periodic,} \end{array} \right. \quad 1 \le j \le d. \qquad (3.7) $$

Note that $\mathbf{X}_i$ is not extended by zero inside $\eta T$, and that $\alpha_i$ is defined up to some constant inside $\eta T$. This definition is convenient for characterizing $\mathbf{X}_i$ from the kernel of $S_{\mathrm{per}, \eta T}^{\kappa}$.
By linearity, $\mathbf{u} = \sum_{i=1}^{d} \left( \int_{\partial(\eta T)} \phi \cdot e_i \, d\sigma \right) \mathbf{X}_i$ in $P$ and there may exist constants $C_{\pm}$ such that

$$p = \sum_{i=1}^{d} \left( \int_{\partial(\eta T)} \phi \cdot e_i \, d\sigma \right) \alpha_i + C_{\pm} \mathbf{n},$$

where $C_+$ (resp. $C_-$) stands for the additional pressure constant inside $P \setminus (\eta T)$ (resp. $\eta T$). By using the jump relation, we deduce

$$\phi = \sum_{i=1}^{d} \left[ \int_{\partial(\eta T)} \phi \cdot e_i \, d\sigma \right] \left[ \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \right] - (C_+ - C_-)\mathbf{n}. \quad (3.11)$$

This proves the inclusion (3.8). Now, if $\mathbf{u} := S_{\text{per},\eta T}^{\kappa} [\left[ \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \right]] = 0$ on $\partial(\eta T)$, then integrating $-\Delta \mathbf{X}_i + \nabla \alpha_i = \mathbf{e}_i$ over $P \setminus (\eta T)$ yields:

$$\int_{\partial(\eta T)} [\sigma(\mathbf{X}_i, \alpha_i)\mathbf{n}] \, d\sigma = \mathbf{e}_i. \quad (3.12)$$

Since $(\mathbf{u}, p)$ solves (3.10), we find $\mathbf{u} = \mathbf{X}_i$ and $p = \alpha_i + 1_{\eta T} C_- + 1_{P \setminus (\eta T)} C_+$. From (3.11) we obtain $C_+ = C_-$. $\square$

In the next proposition, we show that it is possible to specify how $\kappa$ should vary with $\eta$ in order to ensure that $[\sigma(\mathbf{X}_i, \alpha_i)\mathbf{n}]$ be in the kernel of $S_{\text{per},\eta T}^{\kappa}$. We recall the definitions (2.28) and (2.29) of the exterior Stokes flows $(\mathbf{\Phi}_i, \sigma_i)$ in dimension $d \ge 3$ and $d = 2$, and we denote (when $d = 2$) by $\mathcal{F}^{\infty}\mathbf{\Phi}$ the matrix

$$\mathcal{F}^{\infty}\mathbf{\Phi} := (\mathcal{F}^{\infty}\mathbf{\Phi}_j \cdot \mathbf{e}_i)_{1 \le i,j \le d}.$$

**Proposition 3.3.** There exists a real analytic matrix valued function $\eta \mapsto \kappa_{\eta}$ such that the operator $S_{\text{per},\eta T}^{\kappa_{\eta}}$ defined according to (3.1) satisfies

$$S_{\text{per},\eta T}^{\kappa_{\eta}} [\left[ \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \right]] = 0 \text{ on } \partial(\eta T), \text{ for all } 1 \le i \le d,$$

for any $\eta$ in a sufficiently small neighborhood of zero. Moreover, the constants determining the pressure $\alpha_i$ can be chosen in a way that $\kappa_{\eta}$ and $[\sigma(\mathbf{X}_i, \alpha_i)\mathbf{n}]$ admit the following analytic series representations:

$$\begin{align}
\kappa_{\eta} &= M^{-1} + \sum_{p \ge d} \eta^p \kappa_p, & \text{if } d \ge 3 \nonumber \\
\kappa_{\eta} &= \mathcal{F}^{\infty}\Phi + \sum_{p \ge 2} \eta^p \kappa_p, & \text{if } d = 2 \tag{3.13}
\end{align}$$

$$[\sigma(\mathbf{X}_i, \alpha_i)\mathbf{n}] = \eta^{1-d}\sigma(\mathbf{\Phi}_i, \tau_i)\mathbf{n} \circ \tau_{\eta}^{-1} + \sum_{p \ge 1} \eta^p \phi_{i,p} \circ \tau_{\eta}^{-1}, \quad \forall d \ge 2,$$

for some constant $d \times d$ matrices $(\kappa_p)_{p \ge d}$ and functions $(\phi_{i,p})_{p \ge 1}$ of $V^{-1/2}(\partial T, \mathbb{R}^d)$ satisfying

$$\int_{\partial T} \phi_{i,p} \, d\sigma = 0 \text{ for all } p \ge 1.$$

*Proof.* For a given $1 \le i \le d$, consider the function $F_i := (V^{-1/2} \times \mathbb{R}^d) \times \mathbb{R} \to V^{1/2} \times \mathbb{R}^d$ defined by

$$F_i((\boldsymbol{\phi}, \mathbf{a}), \eta) := \left( S_T[\boldsymbol{\phi}] + \mathbf{a} + \eta^{d-2} \int_{\partial T} R(\eta(\cdot - y'))\boldsymbol{\phi}(y') \, d\sigma(y'), \int_{\partial T} \boldsymbol{\phi} \, d\sigma - \mathbf{e}_i \right).$$

The partial Fréchet derivative of $F_i$ with respect to the variable $(\boldsymbol{\phi}, \mathbf{a})$ at $\eta = 0$ is precisely the map $\mathcal{A}$ of (2.36). Moreover $F$ is analytic with respect to $\eta$ and

$$\left\{
\begin{array}{ll}
F_i((\sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n} - C\mathbf{n}, M^{-1}\mathbf{e}_i), 0) = (0,0), & \text{if } d \ge 3, \\
F_i((\sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n} - C\mathbf{n}, \mathcal{F}^{\infty}\mathbf{\Phi}_i), 0) = (0,0), & \text{if } d = 2,
\end{array}
\right.
\quad (3.14)$$

where $C$ is a constant set to ensure $\sigma(\mathbf{\Phi}_i, \tau_i)\mathbf{n} - C\mathbf{n} \in V^{-1/2}$. According to Proposition 2.6, the assumptions of the analytic implicit function theorem are fulfilled and guarantee the existence of analytic mappings $\eta \mapsto \mathbf{a}_{i,\eta} \in \mathbb{R}^d$ and $\eta \mapsto \phi_{i,\eta} \in V^{-1/2}$ such that $F_i((\mathbf{a}_{i,\eta}, \phi_{i,\eta}), \eta) = 0$ for all $\eta$ in a neighborhood of zero $\mathcal{V}$. We define then the matrix $\kappa_{\eta}$ by gathering the vectors $(\mathbf{a}_{j,\eta})$ into columns, namely $\kappa_{\eta,ij} := \mathbf{a}_{j,\eta} \cdot \mathbf{e}_i$ for all $1 \le i,j \le d$. Due to $R(0) = 0$ and $\nabla R(0) = 0$ (Proposition 2.9), it is readily verified that the first order variation of these mappings are of order $O(\eta^d)$, leading to the existence of coefficients $\kappa_p \in \mathbb{R}^{d \times d}$ and functions $\phi_{i,p} \in V^{-1/2}$ such that

$$\begin{align}
\kappa_{\eta} &= M^{-1} + \sum_{p \ge d} \eta^p \kappa_p \quad \text{and} \quad \phi_i = \sigma(\mathbf{\Psi}_i, \tau_i)\mathbf{n} + \sum_{p \ge d} \eta^p \phi_{i,p}, & \text{if } d \ge 3, \nonumber \\
\kappa_{\eta} &= \mathcal{F}^{\infty}\Phi + \sum_{p \ge 2} \eta^p \kappa_p \quad \text{and} \quad \phi_i = \sigma(\mathbf{\Phi}_i, \sigma_i)\mathbf{n} + \sum_{p \ge 2} \eta^p \phi_{i,p}, & \text{if } d = 2. \tag{3.15}
\end{align}$$

Then, by construction, for all $\eta \in \mathcal{V}$:

$$S_T^{\kappa_{\eta}}(\eta)[\phi_{i,\eta}] = 0 \text{ on } \partial T \text{ and } \int_{\partial T} \phi_{i,\eta} \, d\sigma = \mathbf{e}_i,$$
**Proposition 3.2.** The following inclusion holds for the kernel of $S_{\text{per},\eta T}^{\kappa}$:

$$ \ker(S_{\text{per},\eta T}^{\kappa}) \subset \operatorname{span}(\mathbf{n}) \oplus \operatorname{span}(\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket)_{1\le i\le d}. \qquad (3.8) $$

Moreover, if $\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket \in \ker(S_{\text{per},\eta T}^{\kappa})$, then $(\mathbf{X}_i, \alpha_i)$ admits the following single layer potential representation:

$$ \mathbf{X}_i(x) = S_{\text{per},\eta T}^{\kappa} [\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket](x), \quad \alpha_i(x) = \mathcal{P}_{\text{per},\eta T} [\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket](x) + C \quad \text{for } x \in P\setminus(\eta T), \qquad (3.9) $$

for some additive pressure constant $C \in \mathbb{R}$.

**Proof.** Let $\phi \in H^{-\frac{1}{2}}(\partial(\eta T), \mathbb{R}^d)$ be such that $S_{\text{per},\eta T}^{\kappa}[\phi] = 0$ on $\partial(\eta T)$. The functions $\mathbf{u} := S_{\text{per},\eta T}^{\kappa}[\phi]$ and $p := \mathcal{P}_{\text{per},\eta T}[\phi]$ solve the following Stokes problem

$$ \left\{ \begin{array}{l} -\Delta \mathbf{u} + \nabla p = \displaystyle\int_{\partial(\eta T)} \phi \, d\sigma \quad \text{in } P \setminus \partial(\eta T), \\ \operatorname{div}(\mathbf{u}) = 0 \quad \text{in } P \setminus \partial(\eta T), \\ \mathbf{u} = 0 \quad \text{on } \partial(\eta T), \\ \mathbf{u} \text{ is P-periodic}, \\ \sigma(\mathbf{u}, p)\mathbf{n} \text{ is P-periodic}. \end{array} \right. \qquad (3.10) $$

By linearity, $\mathbf{u} = \sum_{i=1}^{d} \left( \int_{\partial(\eta T)} \phi \cdot e_i \, d\sigma \right) \mathbf{X}_i$ in $P$ and there may exist constants $C_{\pm}$ such that

$$ p = \sum_{i=1}^{d} \left( \int_{\partial(\eta T)} \phi \cdot e_i \, d\sigma \right) \alpha_i + C_{\pm} \mathbf{n}, $$

where $C_+$ (resp. $C_-$) stands for the additional pressure constant inside $P\setminus(\eta T)$ (resp. $\eta T$). By using the jump relation, we deduce

$$ \phi = \sum_{i=1}^{d} \left[ \int_{\partial(\eta T)} \phi \cdot e_i \, d\sigma \right] \llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket - (C_+ - C_-)\mathbf{n}. \qquad (3.11) $$

This proves the inclusion (3.8). Now, if $\mathbf{u} := S_{\text{per},\eta T}^{\kappa} [\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket] = 0$ on $\partial(\eta T)$, then integrating $-\Delta \mathbf{X}_i + \nabla \alpha_i = e_i$ over $P\setminus(\eta T)$ yields:

$$ \int_{\partial(\eta T)} \llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket \, d\sigma = e_i. \qquad (3.12) $$

Since $(\mathbf{u}, p)$ solves (3.10), we find $\mathbf{u} = \mathbf{X}_i$ and $p = \alpha_i + 1_{\eta T} C_- + 1_{P\setminus(\eta T)} C_+$. From (3.11) we obtain $C_+ = C_-$. $\square$

In the next proposition, we show that it is possible to specify how $\kappa$ should vary with $\eta$ in order to ensure that $\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket$ be in the kernel of $S_{\text{per},\eta T}^{\kappa}$. We recall the definitions (2.28) and (2.29) of the exterior Stokes flows $(\mathbf{\Phi}_i, \sigma_i)$ in dimension $d \ge 3$ and $d = 2$, and we denote (when $d = 2$) by $\mathcal{F}^{\infty}\mathbf{\Phi}$ the matrix

$$ \mathcal{F}^{\infty}\mathbf{\Phi} := (\mathcal{F}^{\infty}\mathbf{\Phi}_j \cdot e_i)_{1\le i,j\le d}. $$

**Proposition 3.3.** There exists a real analytic matrix valued function $\eta \mapsto \kappa_{\eta}$ such that the operator $S_{\text{per},\eta T}^{\kappa_{\eta}}$ defined according to (3.1) satisfies

$$ S_{\text{per},\eta T}^{\kappa_{\eta}} [\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket] = 0 \text{ on } \partial(\eta T), \text{ for all } 1 \le i \le d, $$

for any $\eta$ in a sufficiently small neighborhood of zero. Moreover, the constants determining the pressure $\alpha_i$ can be chosen in a way that $\kappa_{\eta}$ and $\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket$ admit the following analytic series representations:

$$ \begin{align} \kappa_{\eta} &= M^{-1} + \sum_{p \ge d} \eta^p \kappa_p, & \text{if } d \ge 3 \notag \\ \kappa_{\eta} &= \mathcal{F}^{\infty}\mathbf{\Phi} + \sum_{p \ge 2} \eta^p \kappa_p, & \text{if } d = 2 \tag{3.13} \end{align} $$

$$ \llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket = \eta^{1-d} \sigma(\mathbf{\Phi}_i, \tau_i)\mathbf{n} \circ \tau_{\eta}^{-1} + \sum_{p \ge 1} \eta^p \phi_{i,p} \circ \tau_{\eta}^{-1}, \quad \forall d \ge 2, $$

for some constant $d \times d$ matrices $(\kappa_p)_{p \ge d}$ and functions $(\phi_{i,p})_{p \ge 1}$ of $V^{-1/2}(\partial T, \mathbb{R}^d)$ satisfying

$$ \int_{\partial T} \phi_{i,p} \, d\sigma = 0 \text{ for all } p \ge 1. $$

**Proof.** For a given $1 \le i \le d$, consider the function $F_i : (V^{-\frac{1}{2}} \times \mathbb{R}^d) \times \mathbb{R} \to V^{\frac{1}{2}} \times \mathbb{R}^d$ defined by

$$ F_i((\phi, \mathbf{a}), \eta) := \left( S_T[\phi] + \mathbf{a} + \eta^{d-2} \int_{\partial T} R(\eta(\cdot - y'))\phi(y') \, d\sigma(y'), \int_{\partial T} \phi \, d\sigma - e_i \right). $$
where the second equality imply $\int_{\partial T} \phi_{i,p} d\sigma = 0$ for all $p \ge d$. Let us consider then the functions $u_i := S_{\text{per}}^{\kappa_\eta}[\phi_{i,\eta} \circ \tau_\eta^{-1}]$ and $p_i := \mathcal{P}_{\text{per}}[\phi_{i,\eta} \circ \tau_\eta^{-1}]$. Due to the factorization (3.3), $S_{\text{per}}^{\kappa_\eta}[\phi_{i,\eta} \circ \tau_\eta^{-1}] = 0$ on $\partial(\eta T)$ and the functions $u_i$ and $p_i$ solve the problem (3.10) with right-hand side

$$ \int_{\partial(\eta T)} \phi_{i,\eta} \circ \tau_{\eta}^{-1} d\sigma = \eta^{d-1} \int_{\partial T} \phi_{i,\eta} d\sigma = \eta^{d-1} e_i. $$

It follows that $u_i = \eta^{d-1} \mathcal{X}_i$ and that, up to redefining the additive constants determining the pressure, $p_i = \eta^{d-1} \alpha_i$. By using the jump relation we obtain

$$ \phi_{i,\eta} = \eta^{d-1} \llbracket \sigma(\mathcal{X}_i, \alpha_i) \mathbf{n} \rrbracket, $$

and the result is obtained after redefining the functions $\phi_{i,p}$ to be the coefficients of $\eta^{1-d}\phi_{i,\eta}$, as well as the constants of $\alpha_i$ in order to absorb the constant C of (3.14). $\square$

**Remark 3.1.** The precedent proposition requires the dependence of the additive pressure constants determining $\alpha_i$ with respect to $\eta$ to be carefully selected. The adjustment of the pressure constant for the representation (3.13) corresponds to setting the pressure constants arising in $\alpha_i$ such that for any $\eta > 0$ sufficiently small:

$$ \int_{\partial T} \llbracket \sigma(\mathcal{X}_i, \alpha_i) \mathbf{n} \rrbracket \cdot \mathbf{n} \, d\sigma = \int_{\partial T} \mathbf{n} \cdot \sigma(\mathbf{\Phi}_i, \sigma_i) \mathbf{n}|_{+} \, d\sigma. $$

Combining the two previous propositions, we obtain the following result.

**Proposition 3.4.** The solution $(\mathcal{X}_i, \alpha_i)$ to (1.3) admits the following asymptotic series representation:

$$ \left\{ \begin{array}{l} \mathcal{X}_i(x) = \eta^{2-d} \mathbf{\Phi}_i(x/\eta) + \sum_{p=2}^{+\infty} \eta^p \mathbf{v}_{i,p}(x/\eta) + \sum_{p=0}^{+\infty} \eta^p \mathbf{w}_{i,p}(x), \\ \\ \alpha_i(x) = \eta^{1-d} \sigma_i(x/\eta) + \sum_{p=2}^{+\infty} \eta^{p-1} \xi_{i,p}(x/\eta) + \sum_{p=0}^{+\infty} \eta^p \zeta_{i,p}(x), \end{array} \right. \qquad (3.16) $$

where

(i) $(\mathbf{\Phi}_i, \sigma_i)_{1 \le i \le d}$ are the exterior solutions defined in (2.28) and (2.29) for respectively $d \ge 3$ and $d = 2$,

(ii) the expansion converges point-wise for any fixed $x \in P \setminus \{0\}$,

(iii) there exists a constant $C > 0$ independent of $\eta$ such that the following series convergences hold:

$$ \sum_{p=2}^{+\infty} \eta^p ||\nabla \mathbf{v}_{i,p}||_{L^2(\mathbb{R}^d, \mathbb{R}^{d \times d})} \le C, \quad \sum_{p=0}^{+\infty} \eta^p ||\mathbf{w}_{i,p}||_{L^{\infty}(P, \mathbb{R}^d)} \le C, \qquad (3.17) $$

$$ \sum_{p=2}^{+\infty} \eta^{p-1} ||\xi_{i,p}||_{L^2(\mathbb{R}^d \setminus T)} \le C, \quad \sum_{p=0}^{+\infty} \eta^p ||\zeta_{i,p}||_{L^{\infty}(P)/\mathbb{R}} \le C; \qquad (3.18) $$

(iv) $((\mathbf{v}_{i,p}, \xi_{i,p}))_{p \ge 2}$ and $((\mathbf{w}_{i,p}, \zeta_{i,p}))_{p \ge 0}$ are the solutions to the following recursive systems of interior Stokes problems in $P$ and exterior Stokes problems in $\mathbb{R}^d$:

$$ \left\{ \begin{array}{l} -\Delta \mathbf{w}_{i,p} + \nabla \zeta_{i,p} = \begin{cases} e_i & \text{if } p=0, \\ 0 & \text{if } p \ge 1, \end{cases} \quad \text{in } P, \\ \\ \text{div}(\mathbf{w}_{i,p}) = 0 \quad \text{in } P, \\ \\ \mathbf{w}_{i,p} + \nabla^p G \cdot \mathcal{F}^p \mathbf{\Phi}_i \text{ is } P\text{-periodic,} \quad \text{for } 0 \le p \le d, \\ \\ \left\{ \begin{array}{l} \sigma(\mathbf{w}_{i,p}, \zeta_{i,p}) \mathbf{n} \mathbf{n}^T + \sigma(\nabla^p G \cdot \mathcal{F}^p \mathbf{\Phi}_i, \nabla^p \psi \cdot \mathcal{F}^p \mathbf{\Phi}_i) \mathbf{n} \mathbf{n}^T \text{ is } P\text{-periodic,} \\ \\ \mathbf{w}_{i,p} + \nabla^p G \cdot \mathcal{F}^p \mathbf{\Phi}_i + \sum_{k=d-1}^{p-2} \nabla^{k-d+2} G \cdot \mathcal{F}^{k-d+2} \mathbf{v}_{i,p-k} \text{ is } P\text{-periodic,} \\ \\ \sigma(\mathbf{w}_{i,p}, \zeta_{i,p}) \mathbf{n} \mathbf{n}^T + \sigma(\nabla^p G \cdot \mathcal{F}^p \mathbf{\Phi}_i, \nabla^p \psi \cdot \mathcal{F}^p \mathbf{\Phi}_i) \mathbf{n} \mathbf{n}^T \quad \text{for } p > d, \\ \\ + \sum_{k=d-1}^{p-2} \sigma(\nabla^{k-d+2} G \cdot \mathcal{F}^{k-d+2} \mathbf{v}_{i,p-k}, \nabla^{k-d+2} \psi \cdot \mathcal{F}^{k-d+2} \mathbf{v}_{i,p-k}) \mathbf{n} \mathbf{n}^T \text{ is } P\text{-periodic,} \end{array} \right. \end{array} \right. \qquad (3.19) $$
The partial Fréchet derivative of $F_i$ with respect to the variable $(\Phi, a)$ at $\eta = 0$ is precisely the map $\mathcal{A}$ of (2.36). Moreover $F$ is analytic with respect to $\eta$ and

$$ \left\{ \begin{array}{ll} F_i((\sigma(\Phi_i, \sigma_i)n - Cn, M^{-1}e_i), 0) = (0,0), & \text{if } d \ge 3, \\ F_i((\sigma(\Phi_i, \sigma_i)n - Cn, F^{\infty}\Phi_i), 0) = (0,0), & \text{if } d = 2, \end{array} \right. \qquad (3.14) $$

where $C$ is a constant set to ensure $\sigma(\Phi_i, \tau_i)n - Cn \in V^{-\frac{1}{2}}$. According to Proposition 2.6, the assumptions of the analytic implicit function theorem are fulfilled and guarantee the existence of analytic mappings $\eta \mapsto a_{i,\eta} \in \mathbb{R}^d$ and $\eta \mapsto \phi_{i,\eta} \in V^{-\frac{1}{2}}$ such that $F_i((a_{i,\eta}, \phi_{i,\eta}), \eta) = 0$ for all $\eta$ in a neighborhood of zero $\mathcal{V}$. We define then the matrix $\kappa_{\eta}$ by gathering the vectors $(a_{j,\eta})$ into columns, namely $\kappa_{\eta,ij} := a_{j,\eta} \cdot e_i$ for all $1 \le i,j \le d$. Due to $R(0) = 0$ and $\nabla R(0) = 0$ (Proposition 2.9), it is readily verified that the first order variation of these mappings are of order $O(\eta^d)$, leading to the existence of coefficients $\kappa_p \in \mathbb{R}^{d \times d}$ and functions $\phi_{i,p} \in V^{-\frac{1}{2}}$ such that

$$ \begin{aligned} \kappa_{\eta} &= M^{-1} + \sum_{p \ge d} \eta^p \kappa_p \text{ and } \phi_i = \sigma(\Psi_i, \tau_i)n + \sum_{p \ge d} \eta^p \phi_{i,p}, & \text{if } d \ge 3, \\ \kappa_{\eta} &= F^{\infty}\Phi + \sum_{p \ge 2} \eta^p \kappa_p \text{ and } \phi_i = \sigma(\Phi_i, \sigma_i)n + \sum_{p \ge 2} \eta^p \phi_{i,p}, & \text{if } d = 2. \end{aligned} \qquad (3.15) $$

Then, by construction, for all $\eta \in \mathcal{V}$:

$$ S_T^{\kappa_{\eta}}(\eta)[\phi_{i,\eta}] = 0 \text{ on } \partial T \text{ and } \int_{\partial T} \phi_{i,\eta} \, d\sigma = e_i, $$

where the second equality imply $\int_{\partial T} \phi_{i,p} \, d\sigma = 0$ for all $p \ge d$. Let us consider then the functions $u_i := S_{\text{per}}^{\kappa_{\eta}}[\phi_{i,\eta} \circ \tau_{\eta}^{-1}]$ and $p_i := P_{\text{per}}[\phi_{i,\eta} \circ \tau_{\eta}^{-1}]$. Due to the factorization (3.3), $S_{\text{per}}^{\kappa_{\eta}}[\phi_{i,\eta} \circ \tau_{\eta}^{-1}] = 0$ on $\partial(\eta T)$ and the functions $u_i$ and $p_i$ solve the problem (3.10) with right-hand side

$$ \int_{\partial(\eta T)} \phi_{i,\eta} \circ \tau_{\eta}^{-1} \, d\sigma = \eta^{d-1} \int_{\partial T} \phi_{i,\eta} \, d\sigma = \eta^{d-1} e_i. $$

It follows that $u_i = \eta^{d-1} \mathcal{X}_i$ and that, up to redefining the additive constants determining the pressure, $p_i = \eta^{d-1} \alpha_i$. By using the jump relation we obtain

$$ \phi_{i,\eta} = \eta^{d-1} \llbracket \sigma(\mathcal{X}_i, \alpha_i)n \rrbracket, $$

and the result is obtained after redefining the functions $\phi_{i,p}$ to be the coefficients of $\eta^{1-d}\phi_{i,\eta}$, as well as the constants of $\alpha_i$ in order to absorb the constant $C$ of (3.14). $\square$

**Remark 3.1.** The precedent proposition requires the dependence of the additive pressure constants determining $\alpha_i$ with respect to $\eta$ to be carefully selected. The adjustment of the pressure constant for the representation (3.13) corresponds to setting the pressure constants arising in $\alpha_i$ such that for any $\eta > 0$ sufficiently small:

$$ \int_{\partial T} \llbracket \sigma(\mathcal{X}_i, \alpha_i)n \rrbracket \cdot n \, d\sigma = \int_{\partial T} n \cdot \sigma(\Phi_i, \sigma_i)n|_+ \, d\sigma. $$

Combining the two previous propositions, we obtain the following result.

**Proposition 3.4.** The solution $(\mathcal{X}_i, \alpha_i)$ to (1.3) admits the following asymptotic series representation:

$$ \left\{ \begin{array}{l} \mathcal{X}_i(x) = \eta^{2-d} \Phi_i(x/\eta) + \sum_{p=2}^{+\infty} \eta^p v_{i,p}(x/\eta) + \sum_{p=0}^{+\infty} \eta^p w_{i,p}(x), \\ \\ \alpha_i(x) = \eta^{1-d} \sigma_i(x/\eta) + \sum_{p=2}^{+\infty} \eta^{p-1} \xi_{i,p}(x/\eta) + \sum_{p=0}^{+\infty} \eta^p \zeta_{i,p}(x), \end{array} \right. \qquad (3.16) $$

where

(i) $(\Phi_i, \sigma_i)_{1 \le i \le d}$ are the exterior solutions defined in (2.28) and (2.29) for respectively $d \ge 3$ and $d = 2$,

(ii) the expansion converges point-wise for any fixed $x \in P \setminus \{0\}$,

(iii) there exists a constant $C > 0$ independent of $\eta$ such that the following series convergences hold:

$$ \sum_{p=2}^{+\infty} \eta^p ||\nabla v_{i,p}||_{L^2(\mathbb{R}^d, \mathbb{R}^{d \times d})} \le C, \quad \sum_{p=0}^{+\infty} \eta^p ||w_{i,p}||_{L^\infty(P, \mathbb{R}^d)} \le C, \qquad (3.17) $$

$$ \sum_{p=2}^{+\infty} \eta^{p-1} ||\xi_{i,p}||_{L^2(\mathbb{R}^d \setminus T)} \le C, \quad \sum_{p=0}^{+\infty} \eta^p ||\zeta_{i,p}||_{L^\infty(P)/\mathbb{R}} \le C; \qquad (3.18) $$
and

$$
\left\{
\begin{array}{l}
-\Delta \mathbf{v}_{i,p} + \nabla \xi_{i,p} = 0 \quad \text{in } \mathbb{R}^d \setminus \partial T, \\
\text{div}(\mathbf{v}_{i,p}) = 0 \quad \text{in } \mathbb{R}^d \setminus \partial T, \\
\mathbf{v}_{i,p}(y) = -\mathbf{w}_{i,p}(0) - \sum_{k=1}^{p} \frac{1}{k!} \nabla^k \mathbf{w}_{i,p-k}(0) \cdot y^k \quad \text{for } y \in \partial T, \quad p \ge 2, \\
\mathbf{v}_{i,p}(y) = O(|y|^{1-d}) \quad \text{as } |y| \to +\infty,
\end{array}
\right.
\tag{3.20}
$$

where $\mathbf{w}_{i,0}(0) = \mathbf{w}_{i,1}(0) = 0$, $\mathbf{w}_{i,p}(0)$ is determined by the condition that $\mathbf{v}_{i,p}(y) = O(|y|^{1-d})$ as $|y| \to +\infty$ for $p \ge 2$ (namely $\mathcal{F}^0\mathbf{v}_{i,p} = 0$ if $d \ge 3$ and $\mathcal{F}^\infty\mathbf{v}_{i,p} = 0$ if $d = 2$).

Proof. We insert the expansion (3.13) for $[\sigma(\mathcal{X}_i, \alpha_i)\mathbf{n}]$ into (3.5) and (3.6), to read, thanks to (3.9):

$$
\begin{align*}
\mathcal{X}_i(x) = S_{\text{per},\eta T}^{\kappa,\eta} [[\sigma(\mathcal{X}_i, \alpha_i)\mathbf{n}]](x) = & -\eta^{2-d} \sum_{j=1}^{d} M_{ij}^{-1} \mathbf{\Phi}_j \left(\frac{x}{\eta}\right) + \sum_{p \ge 1} \eta^{p+1} S_T[\phi_{i,p}](x/\eta) \\
& + \eta^{2-d} \left( M^{-1} + \sum_{p \ge d} \eta^p \kappa_p \right) e_i + \int_{\partial T} R(x - \eta y') \hat{\phi}_{i,\eta}(y') \, d\sigma(y'),
\end{align*}
$$

where $\hat{\phi}_{i,\eta} := \eta^{d-1}[\sigma(\mathcal{X}_i, \alpha_i)\mathbf{n}]$ is analytic in $\eta$. The decomposition (3.16) follows by setting $\mathbf{v}_{i,p} := S_T[\phi_{i,p-1}]$ and defining the functions $\mathbf{w}_{i,p}$ by identifying equal powers of $\eta$ in the analytic series representation of

$$
\sum_{p=0}^{+\infty} \eta^p \mathbf{w}_{i,p}(x) = \sum_{p \ge 2} \eta^p \kappa_{p+d-2} e_i + \int_{\partial T} R(x - \eta y') \hat{\phi}_{i,\eta}(y') \, d\sigma(y'),
$$

where we note that $\mathbf{w}_{i,0}(0) = \mathbf{w}_{i,1}(0) = 0$ owing to $R(0) = \nabla R(0) = 0$ (equation (2.48)). The summability properties (3.17) and (3.18) are consequences of the analiticity of $\hat{\phi}_{i,\eta}$ and the continuity bound (2.15) for the series in $\mathbf{v}_{i,p}$ and $\xi_{i,p}$, and of the analiticity of the mapping $\eta \mapsto \int_{\partial T} R(x - \eta y') \hat{\phi}_{i,\eta}(y') d\sigma(y')$ in $H^1(P, \mathbb{R}^d)$.
For the pressure $\alpha_i$, we find similarly

$$
\alpha_i(x) = \mathcal{P}_{\text{per},\eta T} [[\sigma(\mathcal{X}_i, \alpha_i)\mathbf{n}]](x) = -\eta^{1-d} \sum_{j=1}^{d} M_{ij}^{-1} \sigma_i(x/\eta) + \sum_{p \ge 1} \eta^p \mathcal{P}_T[\phi_{i,p}] (x/\eta) + \int_{\partial T} s(x - \eta y') \cdot \hat{\phi}_{i,\eta}(y') \, d\sigma(y'),
$$

whence $\xi_{i,p} := \mathcal{P}_T[\Phi_{i,p-1}]$ and $\zeta_{i,p}(x)$ is obtained by identifying the powers in $\eta$ in the last term. For $d=2$, the derivation is identical, noticing that the constant term disappears by using $S_T[\Phi_i](y) = \Phi_i(y) - \mathcal{F}^\infty\Phi_i$.

It remains to derive the set of equations satisfied by the coefficients of the expansions (3.16). Inserting these as an ansatz in $-\Delta \mathcal{X}_i + \nabla \alpha_i = e_i$, we read $-\Delta \mathbf{w}_{i,p} + \nabla \zeta_{i,p} = \delta_{p=0} e_i$ by identifying identical powers of $\eta$. Considering now $x = \eta y$ with $y \in \partial T$, we obtain

$$
\begin{align*}
0 = \mathcal{X}_i(\eta y) &= \sum_{p=2}^{+\infty} \eta^p \mathbf{v}_{i,p}(y) + \sum_{p=0}^{+\infty} \eta^p \mathbf{w}_{i,p}(\eta y) = \sum_{p=2}^{+\infty} \eta^p \mathbf{v}_{i,p}(y) + \sum_{p=0}^{+\infty} \sum_{k=0}^{+\infty} \frac{\eta^{p+k}}{k!} \nabla^k \mathbf{w}_{i,p}(0) \cdot y^k \\
&= \sum_{p=2}^{+\infty} \eta^p \mathbf{v}_{i,p}(y) + \sum_{p=0}^{+\infty} \eta^p \sum_{k=0}^{p} \frac{1}{k!} \nabla^k \mathbf{w}_{i,p-k}(0) \cdot y^k.
\end{align*}
$$

Identifying equal powers of $\eta$ implies $\mathbf{w}_{i,0}(0) = \mathbf{w}_{i,1}(0) = 0$ and (3.20). Considering now $x \in \partial P$ and using the far field expansions of $\mathbf{\Phi}_i$ and $\mathbf{v}_{i,p}$, we find

$$
\begin{align*}
\mathcal{X}_i(x) &= \eta^{2-d} \left( \mathcal{F}^\infty \mathbf{\Phi}_i - \delta_{d=2} \frac{\log \eta}{4\pi} + \sum_{p=0}^{+\infty} \eta^{d-2+p} \nabla^p G(x) \cdot \mathcal{F}^p \mathbf{\Phi}_i \right) + \sum_{p=2}^{+\infty} \sum_{k=1}^{+\infty} \eta^{d-2+p+k} \nabla^k G(x) \mathcal{F}^k \mathbf{v}_{i,p} + \sum_{p=0}^{+\infty} \eta^p \mathbf{w}_{i,p}(x) \\
&= \eta^{2-d} \left( \mathcal{F}^\infty \mathbf{\Phi}_i - \delta_{d=2} \frac{\log \eta}{4\pi} \right) + \sum_{p=0}^{+\infty} \eta^p \nabla^p G(x) \cdot \mathcal{F}^p \mathbf{\Phi}_i + \sum_{p=d+1}^{+\infty} \sum_{k=d-1}^{p-2} \eta^p \nabla^{k-d+2} G(x) \cdot \mathcal{F}^{k-d+2} \mathbf{v}_{i,p-k} \\
&\quad + \sum_{p=0}^{+\infty} \eta^p \mathbf{w}_{i,p}.
\end{align*}
$$

Identifying again equal powers of $\eta$ yields the boundary conditions of (3.19) (the constant being absorbed by the periodicity condition). Finally, let us verify that the systems (3.19) and (3.20) are well-posed by checking that their respective compatibility conditions are satisfied. The compatibility condition for (3.19) reads

$$
\int_{\partial P} \mathbf{n} \cdot \left( \nabla^p G \cdot \mathcal{F}^p \mathbf{\Phi}_i + \sum_{k=d-1}^{p-2} \nabla^{k-d+2} G \cdot \mathcal{F}^{k-d+2} \mathbf{v}_{i,p-k} \right) \, d\sigma = 0,
$$
(iv) $((\mathbf{v}_{i,p}, \xi_{i,p}))_{p \ge 2}$ and $((\mathbf{w}_{i,p}, \zeta_{i,p}))_{p \ge 0}$ are the solutions to the following recursive systems of interior Stokes problems in $\mathbb{P}$ and exterior Stokes problems in $\mathbb{R}^d$:

$$ egin{cases}
-\Delta \mathbf{w}_{i,p} + \nabla \zeta_{i,p} = \begin{cases} \mathbf{e}_i & \text{if } p=0, \\ 0 & \text{if } p \ge 1, \end{cases} & \text{in } \mathbb{P}, \\
\operatorname{div}(\mathbf{w}_{i,p}) = 0 & \text{in } \mathbb{P}, \\
\begin{cases} \mathbf{w}_{i,p} + \nabla^p G \cdot \mathcal{F}^p \mathbf{\Phi}_i \text{ is } \mathbb{P}\text{-periodic}, \\ \sigma(\mathbf{w}_{i,p}, \zeta_{i,p}) \mathbf{n}\mathbf{n}^T + \sigma(\nabla^p G \cdot \mathcal{F}^p \mathbf{\Phi}_i, \nabla^p \psi \cdot \mathcal{F}^p \mathbf{\Phi}_i) \mathbf{n}\mathbf{n}^T \text{ is } \mathbb{P}\text{-periodic}, \end{cases} & \text{for } 0 \le p \le d, \\
\mathbf{w}_{i,p} + \nabla^p G \cdot \mathcal{F}^p \mathbf{\Phi}_i + \sum_{k=d-1}^{p-2} \nabla^{k-d+2} G \cdot \mathcal{F}^{k-d+2} \mathbf{v}_{i,p-k} \text{ is } \mathbb{P}\text{-periodic}, \\
\sigma(\mathbf{w}_{i,p}, \zeta_{i,p}) \mathbf{n}\mathbf{n}^T + \sigma(\nabla^p G \cdot \mathcal{F}^p \mathbf{\Phi}_i, \nabla^p \psi \cdot \mathcal{F}^p \mathbf{\Phi}_i) \mathbf{n}\mathbf{n}^T & \text{for } p > d, \\
+ \sum_{k=d-1}^{p-2} \sigma(\nabla^{k-d+2} G \cdot \mathcal{F}^{k-d+2} \mathbf{v}_{i,p-k}, \nabla^{k-d+2} \psi \cdot \mathcal{F}^{k-d+2} \mathbf{v}_{i,p-k}) \mathbf{n}\mathbf{n}^T \text{ is } \mathbb{P}\text{-periodic},
\end{cases} \tag{3.19} $$

and

$$ egin{cases}
-\Delta \mathbf{v}_{i,p} + \nabla \xi_{i,p} = 0 \text{ in } \mathbb{R}^d \setminus \partial T, \\
\operatorname{div}(\mathbf{v}_{i,p}) = 0 \text{ in } \mathbb{R}^d \setminus \partial T, \\
\mathbf{v}_{i,p}(y) = -\mathbf{w}_{i,p}(0) - \sum_{k=1}^{p} \frac{1}{k!} \nabla^k \mathbf{w}_{i,p-k}(0) \cdot y^k \text{ for } y \in \partial T, \quad p \ge 2, \\
\mathbf{v}_{i,p}(y) = O(|y|^{1-d}) \text{ as } |y| \to +\infty,
\end{cases} \tag{3.20} $$

where $\mathbf{w}_{i,0}(0) = \mathbf{w}_{i,1}(0) = 0$, $\mathbf{w}_{i,p}(0)$ is determined by the condition that $\mathbf{v}_{i,p}(y) = O(|y|^{1-d})$ as $|y| \to +\infty$ for $p \ge 2$ (namely $\mathcal{F}^0 \mathbf{v}_{i,p} = 0$ if $d \ge 3$ and $\mathcal{F}^\infty \mathbf{v}_{i,p} = 0$ if $d = 2$).

*Proof.* We insert the expansion (3.13) for $\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket$ into (3.5) and (3.6), to read, thanks to (3.9):

$$ \begin{aligned} \mathbf{X}_i(x) = S_{\text{per},\eta T} [\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket](x) = & -\eta^{2-d} \sum_{j=1}^{d} M_{ij}^{-1} \mathbf{\Phi}_j \left(\frac{x}{\eta}\right) + \sum_{p \ge 1} \eta^{p+1} S_T[\phi_{i,p}](x/\eta) \\ & + \eta^{2-d} \left( M^{-1} + \sum_{p \ge d} \eta^p \kappa_p \right) e_i + \int_{\partial T} R(x - \eta y') \hat{\phi}_{i,\eta}(y') \, d\sigma(y'), \end{aligned} $$

where $\hat{\phi}_{i,\eta} := \eta^{d-1} \llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket$ is analytic in $\eta$. The decomposition (3.16) follows by setting $\mathbf{v}_{i,p} := S_T[\phi_{i,p-1}]$ and defining the functions $\mathbf{w}_{i,p}$ by identifying equal powers of $\eta$ in the analytic series representation of

$$ \sum_{p=0}^{+\infty} \eta^p \mathbf{w}_{i,p}(x) = \sum_{p \ge 2} \eta^p \kappa_{p+d-2} e_i + \int_{\partial T} R(x - \eta y') \hat{\phi}_{i,\eta}(y') \, d\sigma(y'), $$

where we note that $\mathbf{w}_{i,0}(0) = \mathbf{w}_{i,1}(0) = 0$ owing to $R(0) = \nabla R(0) = 0$ (equation (2.48)). The summability properties (3.17) and (3.18) are consequences of the analiticity of $\hat{\phi}_{i,\eta}$ and the continuity bound (2.15) for the series in $\mathbf{v}_{i,p}$ and $\xi_{i,p}$, and of the analiticity of the mapping $\eta \mapsto \int_{\partial T} R(x - \eta y') \hat{\phi}_{i,\eta}(y') \, d\sigma(y')$ in $H^1(\mathbb{P}, \mathbb{R}^d)$.
For the pressure $\alpha_i$, we find similarly

$$ \alpha_i(x) = \mathcal{P}_{\text{per},\eta T} [\llbracket \sigma(\mathbf{X}_i, \alpha_i)\mathbf{n} \rrbracket](x) = -\eta^{1-d} \sum_{j=1}^{d} M_{ij}^{-1} \sigma_i(x/\eta) + \sum_{p \ge 1} \eta^p \mathcal{P}_T[\phi_{i,p}] (x/\eta) + \int_{\partial T} s(x - \eta y') \hat{\phi}_{i,\eta}(y') \, d\sigma(y'), $$

whence $\xi_{i,p} := \mathcal{P}_T[\phi_{i,p-1}]$ and $\zeta_{i,p}(x)$ is obtained by identifying the powers in $\eta$ in the last term. For $d=2$, the derivation is identical, noticing that the constant term disappears by using $S_T[\mathbf{\Phi}_i](y) = \mathbf{\Phi}_i(y) - \mathcal{F}^\infty \mathbf{\Phi}_i$.

It remains to derive the set of equations satisfied by the coefficients of the expansions (3.16). Inserting these as an ansatz in $-\Delta \mathbf{X}_i + \nabla \alpha_i = \mathbf{e}_i$, we read $-\Delta \mathbf{w}_{i,p} + \nabla \zeta_{i,p} = \delta_{p=0} \mathbf{e}_i$ by identifying identical powers of $\eta$. Considering now $x = \eta y$ with $y \in \partial T$, we obtain

$$ \begin{aligned} 0 = \mathbf{X}_i(\eta y) &= \sum_{p=2}^{+\infty} \eta^p \mathbf{v}_{i,p}(y) + \sum_{p=0}^{+\infty} \eta^p \mathbf{w}_{i,p}(\eta y) = \sum_{p=2}^{+\infty} \eta^p \mathbf{v}_{i,p}(y) + \sum_{p=0}^{+\infty} \sum_{k=0}^{+\infty} \frac{\eta^{p+k}}{k!} \nabla^k \mathbf{w}_{i,p}(0) \cdot y^k \\ &= \sum_{p=2}^{+\infty} \eta^p \mathbf{v}_{i,p}(y) + \sum_{p=0}^{+\infty} \eta^p \sum_{k=0}^{p} \frac{1}{k!} \nabla^k \mathbf{w}_{i,p-k}(0) \cdot y^k. \end{aligned} $$
which is satisfied due to $\mathrm{div}(G) = 0$ in $P$. The compatibility condition for (3.20) reads

$$ - \int_{\partial T} \left( \sum_{k=0}^{p} \nabla^k \boldsymbol{w}_{i,p}(0) \cdot y^k \right) \cdot \boldsymbol{n} \, d\sigma(y) = 0. $$

This condition is also satisfied because $\mathrm{div}(\boldsymbol{w}_{i,p}) = 0$ implies

$$ \mathrm{div}(\nabla^k \boldsymbol{w}_{i,p}(0) \cdot y^k) = \sum_{1 \le i_1, \dots, i_k \le d} \sum_{j=1}^{d} \partial_{y_{i_1} \dots y_{i_k}} w_{i,p,j}(0) \partial_{y_j}(y_{i_1} \dots y_{i_k}) = k(\nabla^{k-1} \mathrm{div}(\boldsymbol{w}_{i,p}))(0) \cdot y^{k-1} = 0. $$

□

The next proposition proves that the asymptotic expansions (3.16) hold in the energy norm.

**Proposition 3.5.** For any $N \in \mathbb{N}$, denote by $(\boldsymbol{X}_i^N, \alpha_i^N)$ the approximations of $(\boldsymbol{X}_i, \alpha_i)$ obtained by truncating the expansions (3.16) at rank $N$:

$$ \left\{ \begin{array}{l} \boldsymbol{X}_i^N(x) := \eta^{2-d} \Phi_i(x/\eta) + \sum_{p=2}^{N} \eta^p \boldsymbol{v}_{i,p}(x/\eta) + \sum_{p=0}^{N} \eta^p \boldsymbol{w}_{i,p}(x), \\ \\ \alpha_i^N(x) := \eta^{1-d} \sigma_i(x/\eta) + \sum_{p=2}^{N} \eta^{p-1} \xi_{i,p}(x/\eta) + \sum_{p=0}^{N} \eta^p \zeta_{i,p}(x). \end{array} \right. \qquad (3.21) $$

There exists a uniform constant $C > 0$ such that for any $N \ge 0$, $(\boldsymbol{X}_i^N, \alpha_i^N)$ is an approximation of the solution $(\boldsymbol{X}_i, \alpha_i)$ to (1.3) at order $O(\eta^{N+1})$ according to the following error bounds:

$$ \begin{aligned} & \|\boldsymbol{X}_i - \boldsymbol{X}_i^N\|_{L^2(P\setminus(\eta T), \mathbb{R}^{d})} + \|\nabla\boldsymbol{X}_i - \nabla\boldsymbol{X}_i^N\|_{L^2(P\setminus(\eta T), \mathbb{R}^{d \times d})} \\ & \qquad + \|\alpha_i - \alpha_i^N\|_{L^2(P\setminus(\eta T))/\mathbb{R}} \le C\eta^{N+1} \quad \text{if } d \ge 3 \text{ or } (d=2 \text{ and } N=0), \end{aligned} \qquad (3.22) $$

$$ \begin{aligned} & |\log \eta|^{-\frac{1}{2}} \|\boldsymbol{X}_i - \boldsymbol{X}_i^N\|_{L^2(P\setminus(\eta T), \mathbb{R}^d)} + \|\nabla\boldsymbol{X}_i - \nabla\boldsymbol{X}_i^N\|_{L^2(P\setminus(\eta T), \mathbb{R}^{d \times d})} \\ & \qquad + \|\alpha_i - \alpha_i^N\|_{L^2(P\setminus(\eta T))/\mathbb{R}} \le C\eta^{N+1} \quad \text{if } d=2 \text{ and } N \ge 1. \end{aligned} \qquad (3.23) $$

*Proof.* Let us assume without loss of generality that $N \ge 1$. By utilizing the summability properties (3.17) and (3.18), we find the existence of a constant $C$ uniform in $N$ and $\eta$ such that

$$ \begin{aligned} \sum_{p=N+1}^{+\infty} \eta^{p-N-1} \|\nabla(\boldsymbol{v}_{i,p}(\cdot/\eta))\|_{L^2(P\setminus(\eta T), \mathbb{R}^{d \times d})} &= \eta^{\frac{d}{2}-1} \sum_{p=N+1}^{+\infty} \eta^{p-N-1} \|\nabla\boldsymbol{v}_{i,p}\|_{L^2(\eta^{-1}P\setminus T, \mathbb{R}^{d \times d})} \\ &\le \eta^{\frac{d}{2}-1} \sum_{p=N+1}^{+\infty} \eta^{p-N-1} \|\nabla\boldsymbol{v}_{i,p}\|_{L^2(\mathbb{R}^d \setminus T, \mathbb{R}^{d \times d})} \le C\eta^{\frac{d}{2}-1}, \end{aligned} \qquad (3.24) $$

$$ \begin{aligned} \sum_{p=N+1}^{+\infty} \eta^{p-N-1} ( & \|\boldsymbol{w}_{i,p}\|_{L^2(P\setminus(\eta T), \mathbb{R}^d)} + \|\nabla\boldsymbol{w}_{i,p}\|_{L^2(P\setminus(\eta T), \mathbb{R}^{d \times d})} + \|\zeta_{i,p}\|_{L^2(P\setminus(\eta T))/\mathbb{R}} \\ & \le \sum_{p=N+1}^{+\infty} \eta^{p-N-1} ( \|\boldsymbol{w}_{i,p}\|_{H^1(P, \mathbb{R}^d)} + \|\zeta_{i,p}\|_{L^2(P)/\mathbb{R}} ) \le C. \end{aligned} \qquad (3.25) $$

$$ \begin{aligned} \sum_{p=N+1}^{+\infty} \eta^{p-N-2} \|\xi_{i,p}(\cdot/\eta)\|_{L^2(P\setminus(\eta T))/\mathbb{R}} &= \sum_{p=N+1}^{+\infty} \eta^{p-N-1} \eta^{\frac{d}{2}-1} \|\xi_{i,p}\|_{L^2(\eta^{-1}P\setminus T)} \\ &\le \eta^{\frac{d}{2}-1} \sum_{p=N+1}^{+\infty} \eta^{p-N-1} \|\xi_{i,p}\|_{L^2(\mathbb{R}^d\setminus T)} \le C\eta^{\frac{d}{2}-1}. \end{aligned} $$

Combining (3.24) and the Poincaré inequality (A.3), we also obtain

$$ \sum_{p=N+1}^{+\infty} \eta^{p-N-1} \|\boldsymbol{v}_{i,p}(\cdot/\eta)\|_{L^2(P\setminus(\eta T), \mathbb{R}^{d \times d})} \le \begin{cases} C |\log \eta|^{\frac{1}{2}} & \text{if } d=2, \\ C & \text{if } d \ge 3. \end{cases} $$

By using these bounds, we can directly estimate the norm of $\boldsymbol{X}_i - \boldsymbol{X}_i^N$ and $\alpha_i - \alpha_i^N$:

$$ \begin{aligned} & \|\nabla \boldsymbol{X}_i - \nabla \boldsymbol{X}_i^N\|_{L^2(P\setminus(\eta T), \mathbb{R}^{d \times d})} \le \eta^{N+1} \sum_{p=N+1}^{+\infty} \eta^{p-N-1} (\|\nabla(\boldsymbol{v}_{i,p}(\cdot/\eta))\|_{L^2(P\setminus(\eta T), \mathbb{R}^{d \times d})} + \|\nabla \boldsymbol{w}_{i,p}\|_{L^2(P\setminus(\eta T), \mathbb{R}^{d \times d})}) \\ & \le C\eta^{N+1}, \end{aligned} $$