format: text

text:
# Weak convergence results for inhomogeneous rotating fluid equations

## 1. Introduction

The aim of this article is to study the asymptotics of solutions of rotating fluid equations, in the case when the rotation vector is non homogeneous. We consider a domain $`\mathrm{\Omega }=\mathrm{\Omega }_h\times \mathrm{\Omega }_3`$, where $`\mathrm{\Omega }_h`$ denotes either the whole space $`𝐑^2`$ or any periodic domain of $`𝐑^2`$, and similarly $`\mathrm{\Omega }_3`$ denotes $`𝐑`$ or $`𝐓`$. We are interested in the following system:

(1.1) $`_tu+uu\nu \mathrm{\Delta }u+{\displaystyle \frac{1}{\epsilon }}uB+p=0\text{ on }𝐑^+\times \mathrm{\Omega },`$
$`u=0\text{ on }𝐑^+\times \mathrm{\Omega },`$
$`u_{|t=0}=u^0\text{ on }\mathrm{\Omega }`$
where $`B=be_3`$ is the rotation vector, and $`b`$ is a smooth vector field defined in $`\mathrm{\Omega }_h`$. We shall suppose throughout this paper that $`b`$ does not vanish, and is equal to a positive constant $`b_0`$ perturbated by a smooth $`C_c^{\mathrm{}}`$ function; more assumptions on $`b`$ will be made as we go along. Before stating the result we shall prove here, let us recall some well-known facts in the constant case ($`b=1`$). The rotating fluid equations, with $`b`$ constant and homogeneous, modelize the movement of the atmosphere or the oceans at mid-latitudes (see for instance  or ). The fluid is supposed to be incompressible (which corresponds to the hydrostatic approximation), and its viscosity is $`\nu >0`$. The vector field $`u`$ is the velocity and the scalar $`p`$ is the pressure, both are unknown. The parameter $`\epsilon `$ is the Rossby number, and its inverse stands for the speed of rotation of the Earth. Taking the limit $`\epsilon 0`$ means that the scale of motion of the fluid is much smaller than that of the Earth. Note that one can also see $`B`$ as a magnetic field, in which case it makes sense to understand what happens when $`b`$ is not homogeneous; that also holds if one wants to study the movement of the atmosphere in other regions than mid-latitudes.

In the constant case, those equations have been studied by a number of authors. We refer for instance to the works of A. Babin, A. Mahalov and B. Nicolaenko -, I. Gallagher , E. Grenier  for the periodic case, and J.-Y. Chemin, B. Desjardins, I. Gallagher and E. Grenier  for the whole space case as well as  for the case of horizontal plates with Dirichlet boundary conditions (for such boundary conditions we refer also to the work of E. Grenier and N. Masmoudi  as well as N. Masmoudi ). The results in those papers concern both weak and strong solutions; in this article we shall only be concerned with Leray–type weak solutions (): we will see in Section 2 below that their existence is an easy adaptation of the proof of Leray’s existence theorem . In the constant case, it is known that weak solutions converge towards the solution of the two–dimensional Navier-Stokes equations. Such a result in the whole space case is due to Strichartz-type estimates (which are obtained by writing the solution of the linearized problem in Fourier space), whereas in the periodic case it follows from the study of the (discrete) spectrum of the rotating fluid operator (following methods introduced by S. Schochet in ).

The problem here if we want to follow those methods is that it does not seem a good idea to take the Fourier transform of the operator

$$Lu\stackrel{\text{def}}{=}P(uB),u=0,$$

where $`P`$ denotes the Leray projection onto divergence free vector fields, when $`B`$ is not homogeneous; moreover the study of the spectrum of $`L`$ is not an easy matter. So our strategy to study this problem is first to try and recover the well-known results of the constant case without using any information on the spectrum of $`L`$ (other than the determination of its kernel), and without using the Fourier transform. This will be achieved in Section 3. Then the study of the variable case will be an adaptation of the constant case, in Section 4.

Before stating the results we shall prove in this paper, let us comment on the difficulties compared with the constant case: as stated above, it is easy to construct a bounded family of weak solutions to our problem, whether $`b`$ is constant or not. Hence one can construct a weak limit point $`\overline{u}`$, and the question we want to address is to find the equation satisfied by $`\overline{u}`$. Of course the problem consists in taking the limit in the non linear part of the equation. As noted above, we do not wish to study the spectrum of the operator $`L`$ since that seems to be a difficult issue. So we cannot apply the usual, constant $`b`$ methods, as to our knowledge they all involve spectral properties of $`L`$. The idea therefore is to turn to what is known as “weak compactness methods”, in the spirit of Lions and Masmoudi - (for the incompressible limit). We shall recall briefly below what those methods are, and then we shall state the main results of this paper.

### 1.1.  Weak compactness methods

Let us explain what weak compactness methods are all about. The idea is as follows: as usual the trouble to find the limit of the equation comes from the bilinear terms. They can be separated into three categories:

* products involving only elements of the kernel of the penalization $`L`$, which can be shown to be compact (see Corollary 2.5);
* products of elements of the kernel against elements of $`(\text{Ker}L)^{}`$, for which one can take the limit since elements of $`(\text{Ker}L)^{}`$ converge weakly to zero (see Corollary 2.5);
* products involving only elements of $`(\text{Ker}L)^{}`$, which are the problem.

The idea now is to prove that in the last situation, the limit is in fact zero for algebraic reasons: in previous works on rotating fluids, that result was proved essentially by writing the product of two elements of $`(\text{Ker}L)^{}`$ by projection onto eigenvectors of $`L`$. In the periodic case, a “miracle” in the formulation yielded the result (see - or ), whereas in the whole space case, Strichartz estimates did the job (and the convergence was strong), see . In this paper we will show that the result has in fact not much to do with spectral properties of $`L`$, but is due to simple algebraic properties. Let us recall the result in the case of the incompressible limit, where such properties were first used (see ).

###### Proposition 1.1.

 Let $`(\rho _\epsilon ),(u_\epsilon ),(\theta _\epsilon )`$ be bounded families of $`L^2([0,T],H^1(\mathrm{\Omega }))`$ such that

$$\rho _\epsilon \rho ,u_\epsilon u,\theta _\epsilon \theta \text{ as }\epsilon 0.$$

Assume that

$$\begin{array}{c}\epsilon _t\rho _\epsilon +u_\epsilon =0,\hfill \\ \epsilon _tu_\epsilon +(\rho _\epsilon +\theta _\epsilon )=s_\epsilon ,\hfill \\ \epsilon _t\theta _\epsilon +\frac{2}{3}u_\epsilon =s_\epsilon ^{},\hfill \end{array}$$

where $`s_\epsilon ,s_\epsilon ^{}0`$ in $`L^1([0,T],H^s(\mathrm{\Omega }))`$ for some $`s>0`$. Then

$$P(u_\epsilon u_\epsilon )P(uu)\text{ and }(u_\epsilon \theta _\epsilon )(u\theta )$$

in the sense of distributions.

###### Proof.

This result has to be compared with the so-called “compensated compactness” theorems, in the sense that the convergences of some quadratic quantities in $`\rho _\epsilon ,u_\epsilon ,\theta _\epsilon `$ are established under the assumption that some combinations of the derivatives of these functions converge strongly in time to 0. The proof consists in checking that the acoustic oscillations do not bring any contribution to the limiting terms. We introduce the following decompositions:

$$u_\epsilon =Pu_\epsilon +\psi _\epsilon ,\text{ and }\theta _\epsilon =\frac{3\theta _\epsilon 2\rho _\epsilon }{5}+\pi _\epsilon ,$$

so that

$`Pu_\epsilon \text{ and }{\displaystyle \frac{3\theta _\epsilon 2\rho _\epsilon }{5}}\text{ are bounded in }W^{1,1}([0,T],H^s(\mathrm{\Omega })),\text{and}`$
$`\epsilon _t\psi _\epsilon +\pi _\epsilon =(IdP)s_\epsilon 0,\epsilon _t\pi _\epsilon +{\displaystyle \frac{2}{3}}\mathrm{\Delta }\psi _\epsilon ={\displaystyle \frac{2}{5}}s_\epsilon ^{}0\text{ in }L^1([0,T],H^s(\mathrm{\Omega })).`$
We shall note in the following $`S_\epsilon \stackrel{\text{def}}{=}(IdP)s_\epsilon `$ and $`S_\epsilon ^{}\stackrel{\text{def}}{=}{\displaystyle \frac{2}{5}}s_\epsilon ^{}`$. The incompressibility and Boussinesq relations

$$u=0,(\rho +\theta )=0$$

allow to identify the limits

$`Pu_\epsilon u,{\displaystyle \frac{3\theta _\epsilon 2\rho _\epsilon }{5}}\theta `$ $`\text{ in }L^2([0,T]\times \mathrm{\Omega }),`$
$`\psi _\epsilon 0,\pi _\epsilon 0`$ $`\text{ in }wL^2([0,T]\times \mathrm{\Omega }),`$
from which we deduce that, in the sense of distributions

$`u_\epsilon u_\epsilon \psi _\epsilon \psi _\epsilon `$ $`uu,`$
$`\theta _\epsilon u_\epsilon \pi _\epsilon \psi _\epsilon `$ $`\theta u.`$
The key argument is therefore the following formal computation (which can be made rigorous by introducing regularizations with respect to the space variable $`x`$)

$`P.(\psi _\epsilon \psi _\epsilon )`$ $`={\displaystyle \frac{1}{2}}P|\psi _\epsilon |^2+P(\mathrm{\Delta }\psi _\epsilon \psi _\epsilon )`$
$`={\displaystyle \frac{3}{2}}P(_t(\epsilon \pi _\epsilon \psi _\epsilon )\pi _\epsilon \pi _\epsilon +\pi _\epsilon S_\epsilon +S_{}^{}{}_{\epsilon }{}^{}\psi _\epsilon )`$
$`={\displaystyle \frac{3}{2}}P(_t(\epsilon \pi _\epsilon \psi _\epsilon )+\pi _\epsilon S_\epsilon +S_{}^{}{}_{\epsilon }{}^{}\psi _\epsilon ),`$
$`.(\pi _\epsilon \psi _\epsilon )`$ $`=\pi _\epsilon \mathrm{\Delta }\psi _\epsilon +\psi _\epsilon .\pi _\epsilon `$
$`={\displaystyle \frac{3}{2}}\pi _\epsilon (S_{}^{}{}_{\epsilon }{}^{}\epsilon _t\pi _\epsilon )+\psi _\epsilon .(S_\epsilon \epsilon _t\psi _\epsilon )`$
$`={\displaystyle \frac{3}{2}}\pi _\epsilon S_{}^{}{}_{\epsilon }{}^{}+\psi _\epsilon .S_\epsilon {\displaystyle \frac{3\epsilon }{4}}_t|\pi _\epsilon |^2{\displaystyle \frac{\epsilon }{2}}_t|\psi _\epsilon |^2,`$
which shows that the contribution of the acoustic oscillations is negligible. ∎

Inspired by the previous computation, we shall in this article try to use a similar method in the case of rotating fluids: we refer to the proofs of Propositions 3.4 and 4.4 for precise computations.

### 1.2.  Main results

Since we consider incompressible flows, we introduce the following subspaces of $`L^2(\mathrm{\Omega })`$ and $`H^1(\mathrm{\Omega })`$

$$\text{H}=\{uL^2(\mathrm{\Omega })/u=0\},\text{V}=\{uH^1(\mathrm{\Omega })/u=0\}.$$

We will also use the following notation for the inhomogeneous Sobolev spaces

$$H^s(\mathrm{\Omega })=\{u𝒟^{}(\mathrm{\Omega })/(Id\mathrm{\Delta })^{s/2}uL^2(\mathrm{\Omega })\}.$$

Similarly homogeneous Sobolev spaces will be defined by

$$\dot{H}^s(\mathrm{\Omega })=\{u𝒟^{}(\mathrm{\Omega })/(\mathrm{\Delta })^{s/2}uL^2(\mathrm{\Omega })\}.$$

It will appear clearly in the following that the horizontal variables play a special role in this problem. Consequently we shall use the following notation: if $`x`$ is a point in $`\mathrm{\Omega }`$, then we shall note its cartesian coordinates by $`(x_1,x_2,x_3)`$, and the horizontal part of $`x`$ will be denoted $`x_h\stackrel{\text{def}}{=}(x_1,x_2)\mathrm{\Omega }_h`$. Similarly we will denote the horizontal part of any vector field $`f`$ by $`f_h`$, the horizontal gradient by $`_h\stackrel{\text{def}}{=}(_1,_2)`$ and its orthogonal by $`_h^{}=(_2,_1)`$, and the horizontal divergence and Laplacian respectively by $`\text{div}_hf\stackrel{\text{def}}{=}_1f_1+_2f_2=_hf_h`$ and $`\mathrm{\Delta }_h\stackrel{\text{def}}{=}_1^2+_2^2`$.

Finally as usual, $`C`$ will denote a constant which can change from line to line, and $`p`$ will denote the gradient of a function which can also change from line to line.

Before stating the main theorems of this paper, let us give some additional definitions. We will note by

$$S\stackrel{\text{def}}{=}\{x\mathrm{\Omega }/b(x)=0\}\text{and}𝒪\stackrel{\text{def}}{=}\{x\mathrm{\Omega }/b(x)0\}.$$

Finally $`\mathrm{SS}`$ will be the interior of the singular set $`S`$. We will assume in the sequel that

$`(H0)`$
$$\text{ The set }\mathrm{\Omega }\left(\mathrm{SS}𝒪\right)\text{ is of Lebesgue measure 0 ;}$$

$`(H1)`$
$$\mathrm{SS}\text{ is a smooth domain,}$$

$`(H2)`$
$$\begin{array}{c}\text{ On each connected component }𝒪_j\text{ of }𝒪,\text{ there is a smooth function }\sigma _j\\ \text{ such that }(b,\sigma _j,x_3)\text{ is a global smooth coordinate system }\\ \text{ and }𝒪_j=\{x_h𝐑^2|(b(x_h),\sigma _j(x_h))B_j\times \mathrm{\Sigma }_j)\}\times \mathrm{\Omega }_3.\end{array}$$

Now we are ready to state the main theorems of this paper. The first result, rather standard, shows that there are weak solutions to the system (1.1).

###### Theorem 1.

Let $`u^0`$ be any vector field in H. Then for all $`\epsilon >0`$, Equation (1.1) has at least one weak solution $`u_\epsilon L^{\mathrm{}}(𝐑^+,\text{H})L^2(𝐑^+,\dot{H}^1)`$. Moreover, for all $`t>0`$, the following energy estimate holds:

(1.2)
$$u_\epsilon (t)_{L^2}^2+2\nu _0^tu_\epsilon (s)_{L^2}^2𝑑su^0_{L^2}^2.$$

Now the aim of the paper is to describe the limit of $`u_\epsilon `$ as $`\epsilon `$ goes to zero. We will first concentrate on the constant case.

###### Theorem 2.

Suppose that $`B=be_3`$ where $`b`$ is constant and homogeneous. Let $`u^0`$ be any vector field in H, and let $`u_\epsilon `$ be any weak solution of (1.1) in the sense of Theorem 1. Then $`u_\epsilon `$ converges weakly in $`L_{loc}^2(𝐑^+\times \mathrm{\Omega })`$ to a limit $`\overline{u}`$ which if $`\mathrm{\Omega }_3=𝐑`$ is zero, and if $`\mathrm{\Omega }_3=𝐓`$ is the solution of the two dimensional Navier–Stokes equations

$$(NS2D)_t\overline{u}\nu \mathrm{\Delta }_h\overline{u}+\overline{u}_h_h\overline{u}=(_hp,0),\text{div}_h\overline{u}_h=0,\overline{u}_{|t=0}=_𝐓u^0(x_h,x_3)𝑑x_3.$$

###### Remark 1.2.

This theorem is by no means a novelty, it is even rather less precise than other such results one can find in the literature (-, , , ). As we will see in Section 3, the interest of this result lies in its proof, which contrary to the references above, does not depend on the boundary conditions (which can be the whole space or periodic, in each direction).

Now let us state the new result of this paper, concerning the case when $`b`$ is not homogeneous.

###### Theorem 3.

Suppose that $`B=be_3`$ where $`b=b(x_h)`$ is a nonnegative smooth function, say a $`C_c^{\mathrm{}}(\mathrm{\Omega }_h)`$ perturbation of a constant, and where assumptions $`(H0)`$ to $`(H2)`$ are satisfied. Let $`u^0`$ be any vector field in H, and let $`u_\epsilon `$ be any weak solution of (1.1) in the sense of Theorem 1. Then $`u_\epsilon `$ converges weakly in $`L_{loc}^2(𝐑^+\times \mathrm{\Omega })`$ to a limit $`\overline{u}`$ which if $`\mathrm{\Omega }_3=𝐑`$ is zero, and if $`\mathrm{\Omega }_3=𝐓`$ is defined as follows: the third component $`\overline{u}^3`$ satisfies the transport equation

$$_t\overline{u}^3\nu \mathrm{\Delta }_h\overline{u}^3+\overline{u}_h_h\overline{u}^3=0,_3\overline{u}^3=0,\overline{u}_{|t=0}^3=_𝐓u^{0,3}(x_h,x_3)𝑑x_3\text{in}𝐑^+\times (𝒪\mathrm{SS}),$$

while the horizontal component $`\overline{u}_h`$ depends on the region of space considered:

* in $`\mathrm{SS}`$, $`\overline{u}_h`$ satisfies the two dimensional Navier–Stokes equations $`(NS2D)`$ with Dirichlet boundary conditions.
* in $`𝒪`$, $`\overline{u}_h`$ satisfies the following heat equation:

$$_t\overline{u}_h\nu \mathrm{\Pi }\mathrm{\Delta }_h\overline{u}_h=0,$$

where $`\mathrm{\Pi }`$ is the $`L^2`$ orthogonal projection onto the kernel of $`L`$ (which can be extended to $`𝒟^{}(𝒪\mathrm{SS})`$) which satisfies in particular

$$\left(\mathrm{\Pi }\mathrm{\Delta }_h\overline{u}_h|\overline{u}_h\right)_{L^2(𝒪)}=_h\overline{u}_h_{L^2(𝒪)}^2.$$

###### Remark 1.3.

In the regions where $`b`$ is homogeneous, we recover at the limit the 2D Navier–Stokes equations as usual. The Dirichlet boundary conditions appear quite naturally, considering that on the other side of the boundary one finds that $`\overline{u}_h`$ is proportional to $`_h^{}b`$ which vanishes on the boundary of $`\mathrm{SS}`$. More surprising is certainly the result in the region where $`b`$ is not homogeneous. This can be understood as some sort of turbulent behaviour, where all scales are mixed due to the variation of $`b`$. Technically the result is due to the fact that the kernel of $`L`$ is very small as soon as $`b`$ is not a constant, which induces a lot of rigidity in the limit equation.

The structure of the paper is as follows. In the next section, we present the operator $`L`$ and study its main properties (proof of Theorem 1, kernel, projections onto $`\text{Ker}L`$). The following section is devoted to the proof of Theorem 2. Although the result is not new, we present an alternative proof which holds regardless of the domain (with no boundary). This serves as a warm–up to the final section, in which the general variable case is presented, with the proof of Theorem 3.

###### Remark 1.4.

One can wonder about what remains of those results under more general assumptions on the rotation vector $`B`$. First let us consider the case when $`B=b(t,x_h)e_3`$ depends also on time. The results of Section 2 are identical, and if the sets $`\mathrm{SS}`$ and $`𝒪`$ are independent of time one recovers the same type of theorem as in the constant case. In particular the equation on $`\overline{u}`$ in $`𝒪`$ is derived in an identical way to Section 4.2. If the sets $`\mathrm{SS}`$ and $`𝒪`$ do depend on time, then one has to be a little bit more careful and this issue will not be treated here. A more physical problem is the case when the direction of $`B`$ is not fixed, in other words when $`B`$ is a three component vector, depending on all three variables. Then geometrical problems appear, simply to determine the kernel of $`L`$ ; this will be dealt with in a forecoming paper.

## 2. Study of the singular perturbation

### 2.1. Energy estimate

In this section we shall prove Theorem 1 stated in the introduction.

###### Proof.

The structure of the equation (1.1) governing the rotating fluids is very similar to the one of the usual Navier-Stokes equation, since the singular perturbation is just a linear skew-symmetric operator. Therefore weak solutions “ à la Leray ” can be constructed by the approximation scheme of Friedrichs : approximate solutions are obtained by a standard truncation $`J_n`$ of high frequencies. In order to obtain uniform bounds on these approximate solutions, we have just to check that the energy inequality is still satisfied. Computing formally the $`L^2`$ scalar product of (1.1) by $`u`$ leads to

$$\frac{1}{2}\frac{d}{dt}u_{L^2}^2=(\frac{1}{2}(u)|u|^2\nu u\mathrm{\Delta }u+\frac{1}{\epsilon }uuB+up)𝑑x.$$

Integrating by parts (without boundary) and using the incompressibility constraint, we get

$$\frac{1}{2}\frac{d}{dt}u_{L^2}^2=\nu u_{L^2}^2,$$

which holds for any smooth solution of (1.1).

The energy inequality for weak solutions is obtained by taking limits in the approximation scheme. ∎

In particular, the energy estimate provides uniform bounds in $`L^{\mathrm{}}(𝐑^+,\text{H})L^2(𝐑^+,\dot{H}^1)`$ on any family $`(u_\epsilon )_{\epsilon >0}`$ of weak solutions of (1.1) provided that the initial data $`u^0`$ belongs to H.

###### Corollary 2.1.

Let $`u^0`$ be any vector field in H. For all $`\epsilon >0`$, denote by $`u_\epsilon `$ a weak solution of (1.1). Then there exists $`\overline{u}L^{\mathrm{}}(𝐑^+,\text{H})L^2(𝐑^+,\dot{H}^1)`$, such that, up to extraction of a subsequence,

$$u_\epsilon \overline{u}\text{ in }w\text{-}L_{loc}^2(𝐑^+\times \mathrm{\Omega })\text{ as }\epsilon 0.$$

### 2.2.  Characterization of the kernel

We are interested in describing the asymptotic behaviour of $`(u_\epsilon )`$, i.e. in characterizing its limit points. Of course, the equations satisfied by such a limit point $`\overline{u}`$ depend strongly on the structure of the singular perturbation

(2.1)
$$L:u\text{H}P(uB)\text{H}$$

where $`P`$ denotes the Leray projection from $`L^2(\mathrm{\Omega })`$ onto its subspace H of divergence-free vector fields. In particular, we will prove that $`\overline{u}`$ belongs to the kernel $`\text{Ker}(L)`$ of $`L`$, which is characterized in the following proposition.

###### Proposition 2.2.

Define the linear operator $`L`$ by (2.1). Then $`u\text{H}`$ belongs to $`\text{Ker}(L)`$ if and only if there exist $`_h\phi L^2(\mathrm{\Omega }_h)`$ and $`\alpha L^2(\mathrm{\Omega }_h)`$ with

$$_hb_h^{}\phi =0,$$

such that

$$u=_h^{}\phi +\alpha e_3.$$

###### Proof.

We have

$$P(uB)=0.$$

Then, in the sense of distributions,

$$\text{rot}(uB)=0,$$

which can be rewritten

$$(B)u+(B)u(u)B(u)B=0.$$

As $`B=u=0`$ and $`B=be_3`$, we get

(2.2)
$$b_3u(u)be_3=0.$$

In particular, $`_3u_1=_3u_2=0`$ from which we deduce that

(2.3)
$$u_1,u_2L^2(\mathrm{\Omega }_h).$$

Note that in the case where $`\mathrm{\Omega }_3=𝐑`$, the invariance with respect to $`x_3`$ and the fact that $`uL^2(\mathrm{\Omega })`$ imply that $`u_1=u_2=0`$ (and therefore $`u_3=0`$ by the divergence free condition).

Differentiating the incompressibility constraint with respect to $`x_3`$ leads then to

$$_{33}^2u_3=_{13}^2u_1_{23}^2u_2=0$$

in the sense of distributions. The function $`_3u_3`$ depends only on $`x_1`$ and $`x_2`$, and satisfies $`{\displaystyle _3u_3dx_3}=0`$. So $`_3u_3=0`$ and

(2.4)
$$u_3L^2(\mathrm{\Omega }_h),_1u_1+_2u_2=0.$$

Combining (2.3) and (2.4) provides the existence of $`_h\phi L^2(\mathrm{\Omega }_h)`$ such that

$$u_1=_2\phi ,u_2=_1\phi .$$

Replacing in (2.2) leads to

$$_h^{}\phi _hb=b_3u_3=0,$$

which concludes the proof. ∎

Before applying this result to the characterization of the weak limit $`\overline{u}`$, let us just specify it in two important cases. If $`b=0`$ almost everywhere, $`u\text{H}`$ belongs to $`\text{Ker}(L)`$ if and only if

$$u=_h^{}\phi +\alpha e_3,$$

for some $`_h\phi L^2(\mathrm{\Omega }_h)`$ and $`\alpha L^2(\mathrm{\Omega }_h)`$. If $`b0`$ almost everywhere (in other words, if $`\mathrm{\Omega }𝒪`$ is of Lebesgue measure zero), then the condition arising on $`u`$ is much more restrictive : $`u\text{H}`$ belongs to $`\text{Ker}(L)`$ if and only if on each connected component of $`𝒪`$,

$$u=F(b)_h^{}b+\alpha e_3,$$

for some square integrable function $`F(b)`$ and some $`\alpha L^2(\mathrm{\Omega }_h)`$.

From this characterization of $`\text{Ker}(L)`$, we deduce some constraints on the weak limit $`\overline{u}`$.

###### Corollary 2.3.

Let $`u^0`$ be any vector field in H. Denote by $`(u_\epsilon )_{\epsilon >0}`$ a family of weak solutions of (1.1), and by $`\overline{u}`$ any of its limit points. Then, there exist $`\phi L^{\mathrm{}}(𝐑^+,\dot{H}^1(\mathrm{\Omega }_h))L^2(𝐑^+,\dot{H}^2(\mathrm{\Omega }_h))`$ and $`\alpha L^{\mathrm{}}(𝐑^+,L^2(\mathrm{\Omega }_h))L^2(𝐑^+,\dot{H}^1(\mathrm{\Omega }_h))`$ with

$$_hb_h^{}\phi =0,$$

such that

$$\overline{u}=_h^{}\phi +\alpha e_3.$$

###### Proof.

Let $`\chi 𝒟(𝐑^+\times \mathrm{\Omega })`$ be any divergence-free test function. Multiplying (1.1) by $`\epsilon \chi `$ and integrating with respect to all variables leads to

$$u_\epsilon (\epsilon _t\chi +\epsilon u_\epsilon \chi +\epsilon \nu \mathrm{\Delta }\chi +\chi B)𝑑x𝑑t=0.$$

Because of the bounds coming from the energy estimate, we can take limits in the previous identity as $`\epsilon 0`$ to get

$$\overline{u}B\chi dxdt=0.$$

This means that there exists some $`p`$ such that

$$\overline{u}B=p.$$

As $`u_\epsilon `$ satisfies the incompressibility relation for all $`\epsilon >0`$,

$$\overline{u}=0.$$

Then $`\overline{u}(t)\text{Ker}(L)`$ for almost all $`t𝐑^+`$, and we conclude by Proposition 2.2 that there exist $`\phi L^{\mathrm{}}(𝐑^+,\dot{H}^1(\mathrm{\Omega }_h))L^2(𝐑^+,\dot{H}^2(\mathrm{\Omega }_h))`$ and $`\alpha L^{\mathrm{}}(𝐑^+,L^2(\mathrm{\Omega }_h))L^2(𝐑^+,\dot{H}^1(\mathrm{\Omega }_h))`$ with

$$_hb_h^{}\phi =0,$$

such that

$$u=_h^{}\phi +\alpha e_3.$$

### 2.3.  Decomposition by projection on the kernel

To go further in the description of the asymptotic behaviour (i.e. in the characterization of $`\overline{u}`$), we have to isolate the fast oscillations generated by the singular perturbation $`L`$, which produce “big” terms in (1.1), but converge weakly to 0.

Therefore we introduce the following decomposition

$$u_\epsilon =\overline{u}_\epsilon +w_\epsilon ,$$

where $`\overline{u}_\epsilon =\mathrm{\Pi }u_\epsilon `$ is the $`L^2`$ orthogonal projection of $`u_\epsilon `$ on $`\text{Ker}(L)`$ and $`w_\epsilon =\mathrm{\Pi }_{}u_\epsilon `$ is the projection of $`u_\epsilon `$ on $`\text{Ker}(L)^{}`$.

We have seen in the previous paragraph that the characterization of the kernel $`\text{Ker}(L)`$ is strongly linked to the geometry of the vector field $`b`$. In order to obtain further regularity properties on $`\mathrm{\Pi }`$ and $`\mathrm{\Pi }_{}`$, we then need a precise description of the singular set

$$S=\{x\mathrm{\Omega }/b(x)=0\},$$

which justifies Assumptions $`(H0)`$ to $`(H2)`$ given in the introduction.

Before stating the main properties of $`\mathrm{\Pi }`$ and $`\mathrm{\Pi }_{}`$, let us give the following definitions: by Assumption $`(H0)`$ it is enough to describe the limiting function $`\overline{u}`$ on $`𝒪\mathrm{SS}`$. So it is natural to define the following function spaces: for all $`s0`$, $`H^s(𝒪\mathrm{SS})`$ is the closure of $`C_c^{\mathrm{}}(𝒪\mathrm{SS})`$ for the $`H^s`$ norm, and we will note, for $`s0`$, $`H^s(𝒪\mathrm{SS})`$ the dual space of $`H^s(𝒪\mathrm{SS})`$.

It will be useful in the following to note that

(2.5)
$$s0,H^s(𝒪\mathrm{SS})H^s(\mathrm{\Omega })\text{and}s0,H^s(\mathrm{\Omega })H^s(𝒪\mathrm{SS}).$$

###### Proposition 2.4.

Define the linear operator $`L`$ by (2.1). Denote by $`\mathrm{\Pi }`$ the orthogonal projection of H onto $`\text{Ker}(L)`$ and by $`\mathrm{\Pi }_{}=Id\mathrm{\Pi }`$ the orthogonal projection of H onto $`\text{Ker}(L)^{}`$. The operators $`\mathrm{\Pi }`$ and $`\mathrm{\Pi }_{}`$ so defined have natural extensions to tempered distributions on $`𝒪\mathrm{SS}`$, and for all $`s𝐑`$, there exists some $`C_s>0`$ such that for any function $`uH^s(𝒪\mathrm{SS})`$,

$$\mathrm{\Pi }u_{H^s(𝒪\mathrm{SS})}C_su_{H^s(𝒪\mathrm{SS})}\text{and}\mathrm{\Pi }_{}u_{H^s(𝒪\mathrm{SS})}C_su_{H^s(𝒪\mathrm{SS})}.$$

###### Proof.

By Proposition 2.2, for all $`u\text{H}`$

$$\mathrm{\Pi }u=_h^{}\phi +\alpha e_3,$$

for some $`_h\phi L^2(\mathrm{\Omega }_h)`$ and $`\alpha L^2(\mathrm{\Omega }_h)`$ with

$$_hb_h^{}\phi =0.$$

By definition $`\mathrm{\Pi }_{}u=u\mathrm{\Pi }u`$ is orthogonal to any element of $`\text{Ker}(L)`$. So for all $`\beta L^2(\mathrm{\Omega }_h)`$,

$$(u\mathrm{\Pi }u)\beta e_3𝑑x=0,$$

which implies that

(2.6)
$$\alpha =\frac{1}{|\mathrm{\Omega }_3|}u_3𝑑x_3.$$

In order to determine $`\phi `$, we consider separately the domains $`𝒪`$ and $`S`$.

On $`S`$, $`\phi `$ is defined as

$$_h^{}\phi =\left(\frac{1}{|\mathrm{\Omega }_3|}_{\mathrm{\Omega }_3}u_h𝑑x_3\right),$$

and the smoothness properties stated in Proposition 2.4 are obvious.

On $`𝒪`$, we use Assumption $`(H2)`$ which implies that on each connected component $`𝒪_j`$ of $`𝒪`$, $`u`$ can be written $`u(x)=f_j(b,\sigma _j,x_3)`$ where $`f_j`$ has the same smoothness as $`u`$ since the change of coordinate is in $`C^{\mathrm{}}(\mathrm{\Omega })`$. Then clearly the projection $`\mathrm{\Pi }`$ is simply defined by

$$\mathrm{\Pi }u\stackrel{\text{def}}{=}\frac{1}{|\mathrm{\Omega }_3||\mathrm{\Sigma }_j|}_{\mathrm{\Sigma }_j}_{\mathrm{\Omega }_3}f_j(b,\sigma _j,x_3)𝑑\sigma _j𝑑x_3,$$

and the result follows. ∎

###### Corollary 2.5.

Let $`u^0`$ be any vector field in H. Denote by $`(u_\epsilon )_{\epsilon >0}`$ a family of weak solutions of (1.1), and by $`\overline{u}`$ any of its limit points. Consider a subsequence of $`(u_\epsilon )`$ (abusively denoted $`(u_\epsilon )`$) such that

$$u_\epsilon \overline{u}\text{ in }w\text{-}L_{loc}^2(𝐑^+\times \mathrm{\Omega })\text{ as }\epsilon 0.$$

Then, if $`\overline{u}_\epsilon =\mathrm{\Pi }u_\epsilon `$ is the projection of $`u_\epsilon `$ on $`\text{Ker}(L)`$,

$$\overline{u}_\epsilon \overline{u}\text{ strongly in }L_{loc}^2(𝐑^+\times (𝒪\mathrm{SS}))\text{ as }\epsilon 0.$$

###### Proof.

By Proposition 2.4, the projection $`\mathrm{\Pi }`$ is continuous in $`L^2`$. Then, by the energy estimate,

$$\overline{u}_\epsilon \text{ is uniformly bounded in }L_{loc}^2(𝐑^+,𝒦)$$

for some compact subset $`𝒦`$ of $`L^2(\mathrm{\Omega })`$, which provides regularity with respect to space variables.

The second step consists in getting regularity with respect to time. Apply $`\mathrm{\Pi }`$ to the convection equation in (1.1) :

(2.7)
$$_t\overline{u}_\epsilon +\mathrm{\Pi }(u_\epsilon u_\epsilon \nu \mathrm{\Delta }u_\epsilon )=0.$$

As $`u_\epsilon `$ is divergence-free,

$$u_\epsilon u_\epsilon =(u_\epsilon u_\epsilon )\text{ is uniformly bounded in }L^{\mathrm{}}(𝐑^+,W^{1,1}(\mathrm{\Omega })).$$

From the energy estimate, we also deduce that

$$\mathrm{\Delta }u_\epsilon \text{ is uniformly bounded in }L^2(𝐑^+,H^1(\mathrm{\Omega })),$$

from which we infer that

$$u_\epsilon u_\epsilon \nu \mathrm{\Delta }u_\epsilon \text{ is uniformly bounded in }L_{loc}^2(𝐑^+,H_{loc}^{5/2}(\mathrm{\Omega })).$$

Combining this with (2.5), we get by Proposition 2.4

$$_t\overline{u}_\epsilon =\mathrm{\Pi }(u_\epsilon u_\epsilon \nu \mathrm{\Delta }u_\epsilon )\text{ is uniformly bounded in }L_{loc}^2(𝐑^+,H_{loc}^{5/2}(𝒪\mathrm{SS})),$$

which provides the expected regularity in $`t`$.

Aubin’s lemma  then gives the following interpolation result

$$(\overline{u}_\epsilon )\text{ is strongly compact in }L_{loc}^2(𝐑^+\times (𝒪\mathrm{SS})).$$

By Proposition 2.4, $`\mathrm{\Pi }`$ belongs to $`C(L^2(\mathrm{\Omega }),L^2(\mathrm{\Omega }))`$ and therefore to $`C(wL^2(\mathrm{\Omega }),wL^2(\mathrm{\Omega }))`$, from which we deduce that

$$\overline{u}_\epsilon =\mathrm{\Pi }u_\epsilon \mathrm{\Pi }\overline{u}=\overline{u}.$$

Combining both results shows that $`\overline{u}_\epsilon `$ converges strongly to $`\overline{u}`$ in $`L_{loc}^2(𝐑^+\times (𝒪\mathrm{SS})).`$

### 2.4.  Remarks concerning the regularity

#### 2.4.1. Comparison with the gyrokinetic approximation

As mentionned in the introduction, the study of the asymptotics for an inhomogeneous penalization is a natural question in the magnetohydrodynamic framework, when $`B`$ represents the magnetic field. Such a study has been performed for the gyrokinetic approximation , that is for a kinetic model perturbed by a singular magnetic constraint :

* in the case where $`B=b(x_h)e_3`$, the singular limit is exactly the same as in the constant case : the fast rotation has an averaging effect in the plane orthogonal to the magnetic lines;
* in the case where $`B`$ has constant modulus but variable direction, extra drift terms are obtained due to the curvature of the field.

A simplified version of this result can be written as follows.

###### Theorem 4.

 Let $`f^0`$ be a function of $`L^{\mathrm{}}(\mathrm{\Omega }\times 𝐑^3)`$, and $`(f_\epsilon )`$ be a family of solutions of

$$_tf_\epsilon +v_xf_\epsilon +\frac{1}{\epsilon }vB_vf_\epsilon =0,t𝐑_+^{},(x,v)\mathrm{\Omega }\times 𝐑^3,$$

with initial condition

$$f_\epsilon (|t=0)=f^0.$$

Then the family $`(f_\epsilon )`$ is relatively compact in $`wL^{\mathrm{}}(𝐑_+\times \mathrm{\Omega }\times 𝐑^3)`$, as well as the family $`(g_\epsilon )`$ defined by

$$g_\epsilon (t,x,w)=f_\epsilon (t,x,R(x,\frac{t}{\epsilon })w)$$

where $`R(x,\theta )`$ denotes the rotation of angle $`\theta `$ around the oriented axis of direction $`B(x)`$. Moreover,

* if $`B=be_3`$ with $`bC^1(\mathrm{\Omega }_h,𝐑_+^{})`$, any limit point of $`(g_\epsilon )`$ satisfies

$$_tg+v_3_{x_3}g=0;$$
* if $`BC^1(\mathrm{\Omega })`$ with $`_xB=0`$ and $`|B|1`$, any of its limit points satisfies

$$_tg+(wB)B_xg=\frac{1}{2}w(3(w.B)(B_BB)B_wB_{Bw}B)_wg$$

with the notation $`_V\mathrm{\Phi }\stackrel{\text{def}}{=}V\mathrm{\Phi }`$.

The result obtained in this paper is very different because of the incompressibility constraint, which imposes a lot of rigidity to the system. In particular, the kernel of the penalization is much smaller and the limiting system has less degrees of freedom.

#### 2.4.2. A remark in the inviscid case

The weak compactness method used here allows to study the singular limit without regularity with respect to the time variable. However it uses crucially the strong compactness in $`x`$ given by the energy estimate (1.2). Implicitly we have actually considered the penalization

$$L_\epsilon :u\text{V}P(uB)\epsilon \mathrm{\Delta }uH^1(\mathrm{\Omega }).$$

That rules out the possibility to manage an analogous study for inviscid rotating fluids, the first obstacle being to prove the existence of solutions for

$$_tu_\epsilon +(u_\epsilon )u_\epsilon +\frac{1}{\epsilon }u_\epsilon B+p=0,u_\epsilon =0$$

on a uniform time interval $`[0,T]`$. Indeed the operator $`\mathrm{exp}(tL/\epsilon )`$ is not even uniformly bounded on $`H^s(\mathrm{\Omega })`$ for $`s\frac{1}{2}`$

###### Lemma 2.6.

Define the linear operator $`L`$ by (2.1). Then, the group $`(\mathrm{exp}(tL))_{t𝐑}`$ generated by $`L`$ is not uniformly bounded on $`H^s(\mathrm{\Omega })`$ for $`s\frac{1}{2}`$

###### Proof.

The proof of that result is simply due to the fact that by definition of $`\mathrm{\Pi }`$ seen above, the trace of $`\mathrm{\Pi }u`$ is not defined on $`𝒮`$ even if $`u\text{V}`$. So $`\mathrm{\Pi }`$ is not continuous on $`H^s(\mathrm{\Omega })`$ for $`s\frac{1}{2}`$ Then the formula

$$\mathrm{\Pi }=\underset{T+\mathrm{}}{lim}\frac{1}{T}_0^Te^{tL}𝑑t$$

implies that $`e^{tL}`$ cannot be uniformly bounded in $`H^s(\mathrm{\Omega })`$ for $`s\frac{1}{2}`$ which proves the lemma. ∎

## 3. The case of a constant vector field $`B`$ : <br>the 2D Navier-Stokes limit

In the previous section, we have obtained a constraint equation on the limiting velocity field, which expresses that $`\overline{u}`$ belongs to the kernel of the singular perturbation $`L`$. This comes from the fact that the component $`w_\epsilon =\mathrm{\Pi }_{}u_\epsilon `$ has fast oscillations with respect to time, and consequently converges weakly to 0. In the case where $`\mathrm{\Omega }_3=𝐑`$, this characterizes completely the weak limit $`\overline{u}=0`$.

Then it remains to get an evolution equation for $`\overline{u}`$ in the case where $`\mathrm{\Omega }_3=𝐓`$. A natural idea consists in projecting the evolution equation (1.1) for $`u_\epsilon `$ on the kernel of $`L`$, and to study its limit as $`\epsilon 0`$. The difficulty is to take limits in the nonlinear terms : as Corollary 2.5 provides strong compactness on the non-oscillating component $`\overline{u}_\epsilon `$, the problem comes actually to prove that the oscillating terms $`w_\epsilon `$ do not product any constructive interference.

In order to have a good understanding of this phenomena and of its mathematical formulation, we propose to consider first the case where the vector field $`B`$ is constant and homogeneous. The convergence result established here is not so precise as the ones given in - or , since it does not describe the oscillating component and consequently does not provide any strong convergence. Nevertheless the proof is less technical (in particular it does not require any knowledge on the spectral structure of $`L`$), which allows to consider more general cases in the sequel.

### 3.1. Projection on the kernel

In order to obtain the evolution of the limiting velocity field $`\overline{u}`$, the idea is to use the strong convergence

$$\overline{u}_\epsilon =\mathrm{\Pi }u_\epsilon \overline{u}\text{ in }L_{loc}^2(𝐑^+\times (𝒪\mathrm{SS})),$$

rather than the weak convergence

$$u_\epsilon \overline{u}\text{ in }w\text{-}L_{loc}^2(𝐑^+\times \mathrm{\Omega }).$$

Having this idea in view, we first determine the evolution equation for $`\overline{u}_\epsilon `$.

###### Proposition 3.1.

Let $`u^0`$ be any vector field in H. For all $`\epsilon >0`$, denote by $`u_\epsilon `$ a weak solution of (1.1), and by $`\overline{u}_\epsilon =\mathrm{\Pi }u_\epsilon `$ its projection on $`\text{Ker}(L)`$. Then,

(3.1)
$$_t\overline{u}_\epsilon \mathrm{\Pi }(u_\epsilon \text{rot}u_\epsilon )\nu \mathrm{\Delta }\overline{u}_\epsilon =0.$$

###### Proof.

Identity (3.1) is essentially a variant of (2.7). Indeed, in the case of a constant $`B`$, the projection $`\mathrm{\Pi }`$ commutes with any partial derivative, in particular with the Laplacian $`\mathrm{\Delta }`$ :

$$\mathrm{\Pi }(\nu \mathrm{\Delta }u_\epsilon )=\nu \mathrm{\Delta }\mathrm{\Pi }u_\epsilon =\nu \mathrm{\Delta }\overline{u}_\epsilon .$$

Then the key argument is the following identity :

(3.2)
$$u\text{rot}u=\frac{|u|^2}{2}(uu)+uu$$

As $`u_\epsilon `$ is divergence-free,

$$u_\epsilon u_\epsilon =0.$$

As $`\text{Ker}(L)`$ is embedded in the space of divergence-free vector fields,

$$\mathrm{\Pi }\left(\frac{|u_\epsilon |^2}{2}\right)=0.$$

Replacing in (2.7) leads to the expected result. ∎

### 3.2.  Brief description of the oscillations

###### Proposition 3.2.

Let $`u^0`$ be any vector field in H. For all $`\epsilon >0`$, denote by $`u_\epsilon `$ a weak solution of (1.1), and by $`w_\epsilon =\mathrm{\Pi }_{}u_\epsilon `$ its projection on $`\text{Ker}(L)^{}`$. Then there exists $`W_\epsilon L_{loc}^2(𝐑^+\times \mathrm{\Omega })`$ such that $`w_\epsilon =_3W_\epsilon `$. Moreover,

(3.3) $`\epsilon `$ $`_t(_3W_\epsilon W_{\epsilon ,3})+_3W_\epsilon B=r_\epsilon ,`$
$`\epsilon `$ $`_t(_2W_{\epsilon ,1}_1W_{\epsilon ,2})+b(_1W_{\epsilon ,1}+_2W_{\epsilon ,2})=s_\epsilon ,`$
where $`r_\epsilon ,s_\epsilon `$ converge to 0 in $`L_{loc}^2(𝐑^+,H_{loc}^{5/2}(𝒪\mathrm{SS}))`$ as $`\epsilon 0`$.

###### Proof.

As we have supposed in this section that $`\mathrm{SS}=\mathrm{\Omega }`$, the projection of $`u_\epsilon `$ on $`\text{Ker}(L)`$ satisfies, due to Proposition 2.4,

$$\mathrm{\Pi }u_\epsilon =\frac{1}{|\mathrm{\Omega }_3|}u_\epsilon 𝑑x_3.$$

Then

$$w_\epsilon 𝑑x_3=(u_\epsilon \mathrm{\Pi }u_\epsilon )𝑑x_3=0,$$

and there exists $`W_\epsilon `$ such that

$$w_\epsilon =_3W_\epsilon .$$

Moreover, as $`\mathrm{\Omega }_3=𝐓`$, we can always choose $`W_\epsilon `$ so that

$$W_\epsilon 𝑑x_3=0.$$

It remains then to determine the equations for $`W_\epsilon `$. Equation (1.1) implies

$$\epsilon _tw_\epsilon +w_\epsilon be_3+p=\epsilon \nu \mathrm{\Pi }_{}\mathrm{\Delta }u_\epsilon \epsilon \mathrm{\Pi }_{}(u_\epsilon u_\epsilon ),$$

which can be rewritten in terms of $`W_\epsilon `$

(3.4)
$$\epsilon _t_3W_\epsilon +_3W_\epsilon be_3+p=S_\epsilon ,$$

with

$$S_\epsilon =\epsilon (\nu \mathrm{\Pi }_{}\mathrm{\Delta }u_\epsilon \mathrm{\Pi }_{}(u_\epsilon u_\epsilon )).$$

From the energy bound, we deduce that the right hand side in (3.4) is of order $`\epsilon `$ in the space $`L_{loc}^2(𝐑^+,H_{loc}^{3/2}(𝒪\mathrm{SS}))`$. Indeed, using the continuity properties of $`\mathrm{\Pi }_{}`$ as well as (2.5) we get

$$\mathrm{\Pi }_{}\mathrm{\Delta }u_\epsilon _{L^2([0,T],H^1(𝒪\mathrm{SS}))}\mathrm{\Delta }u_\epsilon _{L^2([0,T],H^1(\mathrm{\Omega }))}u_\epsilon _{L^2([0,T],\dot{H}^1(\mathrm{\Omega }))}$$

and

$`\mathrm{\Pi }_{}(u_\epsilon u_\epsilon )_{L^2([0,T],H^{3/2}(K))}`$ $``$ $`u_\epsilon u_\epsilon _{L^2([0,T],H^{3/2}(K))}`$
$``$ $`Cu_\epsilon _{L^{\mathrm{}}(𝐑^+,L^2(\mathrm{\Omega }))}u_\epsilon _{L^2([0,T],L^2(\mathrm{\Omega }))}`$
for all compact subsets $`K`$ of $`𝒪\mathrm{SS}`$. Moreover, the right-hand side in (3.4) belongs to the image of $`\mathrm{\Pi }_{}`$ and can therefore be written as a partial derivative with respect to $`x_3`$,

$$S_\epsilon =_3R_\epsilon \text{with}R_\epsilon 𝑑x_3=0.$$

The equation on the third component provides then

$$\epsilon _t_3W_{\epsilon ,3}+_3p=_3R_\epsilon $$

where $`R_\epsilon `$ converges to 0 in $`L_{loc}^2(𝐑^+,H_{loc}^{3/2}(𝒪\mathrm{SS}))`$ as $`\epsilon 0`$. Integrating with respect to $`x_3`$ provides, since $`p𝑑x_3=0`$,

$$\epsilon _tW_{\epsilon ,3}+p=R_\epsilon .$$

Replacing in (3.4) leads to

$$\epsilon _t_3W_\epsilon +_3W_\epsilon be_3+(R_\epsilon \epsilon _tW_{\epsilon ,3})=S_\epsilon $$

which is the first identity in Proposition 3.2.

In order to establish the second identity, we compute the rotational of (3.4) and write its last component

$$\epsilon _t_3(_2W_{\epsilon ,1}_1W_{\epsilon ,2})+_3(_1W_{\epsilon ,1}+_2W_{\epsilon ,2})=(_2S_{\epsilon ,1}_1S_{\epsilon ,2}).$$

Moreover, the right hand side in (3.4) belongs to the image of $`\mathrm{\Pi }_{}`$ (recall that in this section $`\mathrm{\Pi }_{}`$ commutes with all derivatives) and can therefore be written as a partial derivative with respect to $`x_3`$. Then integrating with respect to $`x_3`$ leads to the expected equality, where the right-hand side converges to 0 in $`L_{loc}^2(𝐑^+,H_{loc}^{5/2}(𝒪\mathrm{SS}))`$ as $`\epsilon `$ goes to $`0`$. ∎

### 3.3. Study of the coupling

The algebraic structure of the propagator (3.3) implies that the oscillating terms cannot interact and produce some contribution in the limiting equation governing $`\overline{u}`$. Indeed the nonlinear term $`w_\epsilon \text{rot}w_\epsilon `$ can be rewritten as the sum of a total derivative with respect to $`x_3`$ and a total derivative with respect to $`t/\epsilon `$, modulo a remainder which converges formally to 0. Then, in order to prove a rigorous convergence result, the first step is to introduce a regularization of the equations (3.3) and to get a control of the source terms in some strong norm.

###### Lemma 3.3.

Let $`u^0`$ be any vector field in H. For all $`\epsilon >0`$, denote by $`u_\epsilon `$ a weak solution of (1.1) in $`L^{\mathrm{}}(𝐑^+,\text{H})L^2(𝐑^+,\dot{H}^1)`$, and by $`w_\epsilon =\mathrm{\Pi }_{}u_\epsilon `$ its projection onto $`\text{Ker}(L)^{}`$. Then, for all $`\delta >0`$, there exists $`W_\epsilon ^\delta L_{loc}^2(𝐑^+,_sH^s(\mathrm{\Omega }))`$ such that

(3.5)
$$w_\epsilon _3W_\epsilon ^\delta 0\text{ in }L_{loc}^2(𝐑^+\times \mathrm{\Omega })\text{ as }\delta 0\text{ uniformly in }\epsilon >0,$$

and

(3.6) $`\epsilon `$ $`_t(_3W_\epsilon ^\delta W_{\epsilon ,3}^\delta )+_3W_\epsilon ^\delta B=r_\epsilon ^\delta ,`$
$`\epsilon `$ $`_t(_2W_{\epsilon ,1}^\delta _1W_{\epsilon ,2}^\delta )+b(_1W_{\epsilon ,1}^\delta +_2W_{\epsilon ,2}^\delta )=s_\epsilon ^\delta ,`$
where for all $`\delta `$, $`r_\epsilon ^\delta `$ and $`s_\epsilon ^\delta `$ converge to 0 in $`L_{loc}^2(𝐑^+,L_{loc}^2(\mathrm{\Omega }))`$ as $`\epsilon 0`$.

###### Proof.

We introduce the following regularization: let $`\kappa C_c^{\mathrm{}}(𝐑^3,𝐑^+)`$ such that $`\kappa (x)=0`$ if $`|x|1`$ and $`\kappa 𝑑x=1`$, we define

$$\kappa _\delta :x\frac{1}{\delta ^3}\kappa \left(\frac{.}{\delta }\right)$$

and

$$w_\epsilon ^\delta =w_\epsilon \kappa _\delta ,W_\epsilon ^\delta =W_\epsilon \kappa _\delta ,$$

so that $`w_\epsilon ^\delta =_3W_\epsilon ^\delta .`$

By the energy estimate, for all $`T>0`$,

$$u_\epsilon \text{ is uniformly bounded in }L^2([0,T],H^1(\mathrm{\Omega })).$$

By Proposition 2.4, we infer that

$$w_\epsilon =\mathrm{\Pi }_{}u_\epsilon \text{ is uniformly bounded in }L^2([0,T],𝒦)$$

for some compact $`𝒦`$ of $`L^2(\mathrm{\Omega })`$. The result (3.5) then follows from the following fact:

$$w_\epsilon ^\delta (t,x)w_\epsilon (t,x)=(w_\epsilon (t,xy)w_\epsilon (t,x))\kappa _\delta (y)𝑑y$$

hence there is some continuity modulus $`\omega `$ such that

$$\epsilon ,|w_\epsilon ^\delta (t,x)w_\epsilon (t,x)|\omega (y)\kappa _\delta (y)𝑑y,$$

and the result follows.

Regularizing (3.3) leads to

$`\epsilon `$ $`_t(_3W_\epsilon ^\delta W_{\epsilon ,3}^\delta )+_3W_\epsilon ^\delta B=r_\epsilon \kappa _\delta ,`$
$`\epsilon `$ $`_t(_2W_{\epsilon ,1}^\delta _1W_{\epsilon ,2}^\delta )+b(_1W_{\epsilon ,1}^\delta +_2W_{\epsilon ,2}^\delta )=s_\epsilon \kappa _\delta ,`$
because $`b`$ is homogeneous. Then, for all $`T>0`$ and all compact subsets $`K`$ of $`𝒪\mathrm{SS}`$ and for $`\delta `$ small enough,

$`r_\epsilon ^\delta _{L^2([0,T],L^2(K))}`$ $`=r_\epsilon \kappa _\delta _{L^2([0,T],L^2(K))}`$
$`C\kappa _\delta _{W^{5/2,1}(𝐑^3)}r_\epsilon _{L^2([0,T],H^{5/2}(K^{}))}`$
$`{\displaystyle \frac{C}{\delta ^{5/2}}}r_\epsilon _{L^2([0,T],H^{5/2}(K^{}))},`$
where $`K^{}`$ is a compact subset of $`𝒪\mathrm{SS}`$ such that

$$\{x\mathrm{\Omega }/d(x,K)\delta \}K^{}.$$

And, in the same way,

$$r_\epsilon ^\delta _{L^2([0,T],L^2(K))}\frac{C}{\delta ^{5/2}}r_\epsilon _{L^2([0,T],H^{5/2}(K^{}))}.$$

For a fixed $`\delta `$, Proposition 3.2 gives the expected convergences. ∎

Equipped with this preliminary result, we are now able to study the coupling between the oscillating terms and to prove the following proposition.

###### Proposition 3.4.

Let $`u^0`$ be any vector field in H. For all $`\epsilon >0`$, denote by $`u_\epsilon `$ a weak solution of (1.1), and by $`w_\epsilon =\mathrm{\Pi }_{}u_\epsilon `$ its projection onto $`\text{Ker}(L)^{}`$. Then,

$$Pw_\epsilon \text{rot}w_\epsilon dx_30\text{ in }𝒟^{}(𝐑^+\times \mathrm{\Omega })\text{ as }\epsilon 0,$$

where $`P`$ denotes the Leray projection.

###### Proof.

We start by proving that

(3.7)
$$P\left(w_\epsilon \text{rot}w_\epsilon dx_3w_\epsilon ^\delta \text{rot}w_\epsilon ^\delta dx_3\right)0\text{in}𝒟^{}(𝐑^+\times \mathrm{\Omega })\text{ as }\delta 0\text{ uniformly in }\epsilon .$$

From the identity (3.2) and the relations $`w_\epsilon =w_\epsilon ^\delta =0`$, we deduce that

$`P`$ $`\left(w_\epsilon \text{rot}w_\epsilon w_\epsilon ^\delta \text{rot}w_\epsilon ^\delta \right)`$
$`=P\left({\displaystyle \frac{|w_\epsilon |^2|w_\epsilon ^\delta |^2}{2}}\right)P(w_\epsilon w_\epsilon w_\epsilon ^\delta w_\epsilon ^\delta )`$
$`=P((w_\epsilon w_\epsilon ^\delta )w_\epsilon )P(w_\epsilon ^\delta (w_\epsilon w_\epsilon ^\delta ))`$
By Lemma 3.3, we get (3.7).

It remains to prove that, for any fixed $`\delta >0`$,

(3.8)
$$w_\epsilon ^\delta \text{rot}w_\epsilon ^\delta dx_30\text{ as }\epsilon 0.$$

Define $`\rho _\epsilon ^\delta =\text{rot}W_\epsilon ^\delta `$ :

$$\left(\begin{array}{c}\rho _{\epsilon ,1}^\delta \hfill \\ \rho _{\epsilon ,2}^\delta \hfill \\ \rho _{\epsilon ,3}^\delta \hfill \end{array}\right)=\left(\begin{array}{c}_2W_{\epsilon ,3}^\delta _3W_{\epsilon ,2}^\delta \hfill \\ _1W_{\epsilon ,3}^\delta +_3W_{\epsilon ,1}^\delta \hfill \\ _1W_{\epsilon ,2}^\delta _2W_{\epsilon ,1}^\delta \hfill \end{array}\right)$$

which is uniformly bounded with respect to $`\epsilon `$ in $`L^2([0,T],_sH^s(\mathrm{\Omega }))`$ by Lemma 3.3. We have

$$w_\epsilon ^\delta \text{rot}w_\epsilon ^\delta =\left(\begin{array}{c}_3W_{\epsilon ,1}^\delta \hfill \\ _3W_{\epsilon ,2}^\delta \hfill \\ _3W_{\epsilon ,3}^\delta \hfill \end{array}\right)\left(\begin{array}{c}_3\rho _{\epsilon ,1}^\delta \hfill \\ _3\rho _{\epsilon ,2}^\delta \hfill \\ _3\rho _{\epsilon ,3}^\delta \hfill \end{array}\right).$$

From (3.6) and the divergence-free relation $`_3W_{\epsilon ,3}^\delta =_1W_{\epsilon ,1}^\delta _2W_{\epsilon ,2}^\delta `$, we deduce that the previous term can be rewritten

$$\left(\begin{array}{c}\epsilon b^1_t(_3W_{\epsilon ,2}^\delta _2W_{\epsilon ,3}^\delta )b^1r_{\epsilon ,2}^\delta \hfill \\ \epsilon b^1_t(_3W_{\epsilon ,1}^\delta +_1W_{\epsilon ,3}^\delta )+b^1r_{\epsilon ,1}^\delta \hfill \\ \epsilon b^1_t(_2W_{\epsilon ,1}^\delta _1W_{\epsilon ,2}^\delta )b^1s_\epsilon ^\delta \hfill \end{array}\right)\left(\begin{array}{c}_3\rho _{\epsilon ,1}^\delta \hfill \\ _3\rho _{\epsilon ,2}^\delta \hfill \\ _3\rho _{\epsilon ,3}^\delta \hfill \end{array}\right),$$

or equivalently

$$\left(\begin{array}{c}\epsilon b^1_t\rho _{\epsilon ,1}^\delta b^1r_{\epsilon ,2}^\delta \hfill \\ \epsilon b^1_t\rho _{\epsilon ,2}^\delta +b^1r_{\epsilon ,1}^\delta \hfill \\ \epsilon b^1_t\rho _{\epsilon ,3}^\delta b^1s_\epsilon ^\delta \hfill \end{array}\right)\left(\begin{array}{c}_3\rho _{\epsilon ,1}^\delta \hfill \\ _3\rho _{\epsilon ,2}^\delta \hfill \\ _3\rho _{\epsilon ,3}^\delta \hfill \end{array}\right).$$

Integrating by parts with respect to $`x_3`$ leads then to

$$b^1\left(\begin{array}{c}\epsilon _t(\rho _{\epsilon ,2}^\delta _3\rho _{\epsilon ,3}^\delta )+_3(\rho _{\epsilon ,2}^\delta \epsilon _t\rho _{\epsilon ,3}^\delta )+(r_{\epsilon ,1}^\delta _3\rho _{\epsilon ,3}^\delta +s_\epsilon ^\delta _3\rho _{\epsilon ,2}^\delta )\hfill \\ \epsilon _t(\rho _{\epsilon ,3}^\delta _3\rho _{\epsilon ,1}^\delta )+_3(\rho _{\epsilon ,3}^\delta \epsilon _t\rho _{\epsilon ,1}^\delta )+(r_{\epsilon ,2}^\delta _3\rho _{\epsilon ,3}^\delta s_\epsilon ^\delta _3\rho _{\epsilon ,1}^\delta )\hfill \\ \epsilon _t(\rho _{\epsilon ,1}^\delta _3\rho _{\epsilon ,2}^\delta )+_3(\rho _{\epsilon ,1}^\delta \epsilon _t\rho _{\epsilon ,2}^\delta )(r_{\epsilon ,2}^\delta _3\rho _{\epsilon ,2}^\delta +r_{\epsilon ,1}^\delta _3\rho _{\epsilon ,1}^\delta )\hfill \end{array}\right),$$

from which we deduce that, for any fixed $`\delta >0`$,

$$Pw_\epsilon ^\delta \text{rot}w_\epsilon ^\delta dx_30$$

in the sense of distributions, as $`\epsilon 0`$.

Combining (3.7) and (3.8) gives the expected convergence. ∎

### 3.4. Passage to the limit

In order to determine the limiting velocity field $`\overline{u}`$, we have now to take limits in (3.1) which can be rewritten

$$_t\overline{u}_\epsilon \mathrm{\Pi }(\overline{u}_\epsilon \text{rot}\overline{u}_\epsilon +w_\epsilon \text{rot}\overline{u}_\epsilon +\overline{u}_\epsilon \text{rot}w_\epsilon +w_\epsilon \text{rot}w_\epsilon )\nu \mathrm{\Delta }\overline{u}_\epsilon =0,$$

using the decomposition

$$u_\epsilon =\overline{u}_\epsilon +w_\epsilon ,$$

where we recall that

$$\overline{u}_\epsilon \overline{u}\text{ strongly in }L_{loc}^2(𝐑^+\times \mathrm{\Omega }),$$

$$w_\epsilon 0\text{ weakly in }L_{loc}^2(𝐑^+\times \mathrm{\Omega }).$$

Then standard arguments using (3.2) show that

$$_t\overline{u}\mathrm{\Pi }(\overline{u}\text{rot}\overline{u})\nu \mathrm{\Delta }\overline{u}=\underset{\epsilon 0}{lim}\mathrm{\Pi }(w_\epsilon \text{rot}w_\epsilon ),$$

and by Proposition 3.4 we get

$$_t\overline{u}\mathrm{\Pi }(\overline{u}\text{rot}\overline{u})\nu \mathrm{\Delta }\overline{u}=0.$$

In the case where $`B`$ is constant, the projection $`\mathrm{\Pi }`$ reduces to a simple averaging with respect to $`x_3`$, and we recover the usual convergence result towards the 2D1/2 Navier-Stokes equation.

## 4. The case of a variable vector field $`B`$ : <br>a turbulent behaviour

In this section we shall prove Theorem 3 stated in the introduction, concerning the general case when the rotation vector $`B=b(x_h)e_3`$ is inhomogeneous. We suppose that Assumptions $`(H0)`$ to $`(H2)`$ are satisfied. If $`\mathrm{\Omega }_3=𝐑`$, then $`\overline{u}=0`$ simply because it is in $`L^2(\mathrm{\Omega })`$ but only depends on the horizontal variables. So from now on we can suppose that $`\mathrm{\Omega }_3=𝐓`$.

The strategy of proof is quite similar to the constant case, so we shall often be referring to the results of the previous section. The first remark to be made is that if $`b`$ is constant in some positive measured region of $`\mathrm{\Omega }_h`$, then in that region the results of the previous section should apply and one should recover at the limit the usual two–dimensional behaviour.

Moreover, the results of Section 2 hold for any $`B`$, so in particular any weak limit point $`\overline{u}`$ of a sequence of weak solutions $`u_\epsilon `$ to (1.1) is in the kernel of $`L`$ according to Corollary 2.3. That means in particular that the third component does not see the difference between $`𝒪`$ and $`\mathrm{SS}`$ since the elements of the kernel of $`L`$ have the same third component whether $`b`$ is homogeneous or not. So in the following, we shall restrict the study of the limit system to the horizontal components only. As in the previous section, the proof of Theorem 3 consists in finding the equation satisfied by $`\overline{u}`$ (at least its horizontal part $`\overline{u}_h`$), by taking the limit of the equation satisfied by the horizontal part of $`\overline{u}_\epsilon =\mathrm{\Pi }u_\epsilon `$. The first result we shall establish is that in the general, variable $`b`$ case, there is no coupling between oscillating vector fields yielding extra terms in the averaged equation. This will be a generalization of Proposition 3.4 to the variable case, and the analysis will follow closely the proof of Proposition 3.4. Then we shall write the averaged equation on $`𝒪`$. Finally we shall concentrate on the $`\mathrm{SS}`$ case and show the limit $`\overline{u}_h`$ satisfies a two–dimensional Navier–Stokes equation with homogeneous boundary conditions on the boundary of $`\mathrm{SS}`$.

### 4.1. The averaged equation

Let $`u_\epsilon `$ be a family of weak solutions to (1.1), and define $`\overline{u}_\epsilon =\mathrm{\Pi }u_\epsilon `$. Recalling that the elements of $`\text{Ker}(L)`$ are divergence free, we have as in the constant case

$$_t\overline{u}_\epsilon \nu \mathrm{\Pi }\mathrm{\Delta }u_\epsilon \mathrm{\Pi }(u_\epsilon \text{rot}u_\epsilon )=0.$$

Of course the projector $`\mathrm{\Pi }`$ no longer commutes with (horizontal) derivatives. However as $`\mathrm{\Pi }`$ belongs to $`C(wH^s(𝒪\mathrm{SS}),wH^s(𝒪\mathrm{SS}))`$ for $`s0`$, we clearly have as $`\epsilon `$ goes to zero,

$$\mathrm{\Pi }\mathrm{\Delta }u_\epsilon \mathrm{\Pi }\mathrm{\Delta }\overline{u}.$$

Let us now take the limit in the quadratic term.

###### Lemma 4.1.

Let $`u^0L^2(\mathrm{\Omega })`$ be a divergence-free vector field. For all $`\epsilon >0`$, denote by $`u_\epsilon `$ a weak solution of (1.1), and by $`\overline{u}_\epsilon =\mathrm{\Pi }u_\epsilon `$ (resp. $`w_\epsilon =\mathrm{\Pi }_{}u_\epsilon `$) its projection onto $`\text{Ker}(L)`$ (resp. $`\text{Ker}(L)^{}`$). Then the following results hold in $`𝒟^{}(𝐑^+\times (𝒪\mathrm{SS}))`$, as $`\epsilon `$ goes to zero:

(4.1)
$$\mathrm{\Pi }(\overline{u}_\epsilon \overline{u}_\epsilon )\mathrm{\Pi }(\overline{u}\overline{u}),$$

(4.2)
$$\mathrm{\Pi }(w_\epsilon \overline{u}_\epsilon +\overline{u}_\epsilon w_\epsilon )0,$$

(4.3)
$$\mathrm{\Pi }(w_\epsilon w_\epsilon )0.$$

###### Proof.

The results (4.1) and (4.2) are simply due to the compactness of $`\overline{u}_\epsilon `$ as well as the fact that $`w_\epsilon `$ goes weakly to zero, results given in Corollary 2.5. We also use the continuity properties of $`\mathrm{\Pi }`$ stated in Proposition 2.4. The more difficult result to prove is of course (4.3). The method will follow the proof of Proposition 3.4, and will be achieved in two steps. First we show that one can smoothen out the equation satisfied by $`w_\epsilon `$, and then we perform some algebra on the bilinear term in the equation, as in the constant case. We shall therefore continuously be referring to the methods of Section 3.

Let us start by proving the following result, analogous to Proposition 3.2.

###### Proposition 4.2.

Let $`u^0L^2(\mathrm{\Omega })`$ be a divergence-free vector field. For all $`\epsilon >0`$, denote by $`u_\epsilon `$ a weak solution of (1.1), and by $`w_\epsilon =\mathrm{\Pi }_{}u_\epsilon `$ its projection onto $`\text{Ker}(L)^{}`$. Then there exists $`W_{\epsilon ,3}L_{loc}^2(𝐑^+\times \mathrm{\Omega })`$ such that $`w_{\epsilon ,3}=_3W_{\epsilon ,3}`$. Moreover,

(4.4) $`\epsilon `$ $`_t(w_{\epsilon ,h}_hW_{\epsilon ,3})+w_{\epsilon ,h}B=\stackrel{~}{r}_\epsilon `$
$`\epsilon `$ $`_t(_2w_{\epsilon ,1}_1w_{\epsilon ,2})+\text{div}_h(bw_{\epsilon ,h})=\stackrel{~}{s}_\epsilon ,`$
where $`\stackrel{~}{r}_\epsilon ,\stackrel{~}{s}_\epsilon `$ converge to 0 in $`L_{loc}^2(𝐑^+,H_{loc}^{5/2}(𝒪\mathrm{SS}))`$ as $`\epsilon 0`$.

###### Proof.

We shall omit the proof of that result here, as it is identical to the proof of Proposition 3.2: we just have to notice that the third component $`w_{\epsilon ,3}`$ is of vertical mean zero, so can as in the constant case be replaced by $`_3W_{\epsilon ,3}`$. The other components, contrary to the constant case, cannot be transformed in that way, so remain as they are. The rest of the proof is identical to the constant case. ∎

Now as in the constant case, let us smoothen out Equation (4.4).

###### Lemma 4.3.

Let $`u^0L^2(\mathrm{\Omega })`$ be a divergence-free vector field. For all $`\epsilon >0`$, denote by $`u_\epsilon `$ a weak solution of (1.1), and by $`w_\epsilon =\mathrm{\Pi }_{}u_\epsilon `$ its projection on $`\text{Ker}(L)^{}`$. Then, for all $`\delta >0`$, there exists $`w_{\epsilon ,h}^\delta L_{loc}^2(𝐑^+,_sH^s(\mathrm{\Omega }))`$ and $`W_{\epsilon ,3}^\delta L_{loc}^2(𝐑^+,_sH^s(\mathrm{\Omega }))`$ such that

(4.5)
$$w_{\epsilon ,h}w_{\epsilon ,h}^\delta 0\text{ in }L_{loc}^2(𝐑^+\times \mathrm{\Omega })\text{ as }\delta 0\text{ uniformly in }\epsilon >0,$$

(4.6)
$$w_{\epsilon ,3}_3W_{\epsilon ,3}^\delta 0\text{ in }L_{loc}^2(𝐑^+\times \mathrm{\Omega })\text{ as }\delta 0\text{ uniformly in }\epsilon >0,$$

(4.7)
$$\delta _hw_{\epsilon ,h}^\delta _{L^2(\mathrm{\Omega })}0\text{ as }\delta 0\text{ uniformly in }\epsilon >0$$

and

(4.8) $`\epsilon `$ $`_t(w_{\epsilon ,h}^\delta _hW_{\epsilon ,3}^\delta )+w_{\epsilon ,h}^\delta B=\stackrel{~}{r}_\epsilon ^\delta `$
$`\epsilon `$ $`_t(_2w_{\epsilon ,1}^\delta _1w_{\epsilon ,2}^\delta )+_1(w_{\epsilon ,1}^\delta b)+_2(w_{\epsilon ,2}^\delta b)=\stackrel{~}{s}_\epsilon ^\delta ,`$
where $`\stackrel{~}{r}_\epsilon ^\delta `$ and $`\stackrel{~}{s}_\epsilon ^\delta `$ converge to 0 in $`L_{loc}^2(𝐑^+,L_{loc}^2(𝒪\mathrm{SS}))`$ as $`\epsilon `$ and $`\delta `$ go to $`0`$. More precisely for any subset $`K`$ of $`𝒪\mathrm{SS}`$ there is a constant $`C`$ (independent of $`\epsilon `$ and $`\delta `$) and a constant $`c_\delta `$ depending only on $`\delta `$ such that

(4.9)
$$\stackrel{~}{r}_\epsilon ^\delta +\stackrel{~}{s}_\epsilon ^\delta _{L^2(K)}c_\delta o_\epsilon (1)+C\delta .$$

###### Proof.

We shall not write all the details of the proof here, since it is very similar to the constant case (Lemma 3.3); let us simply point out where the fact that $`b`$ is not constant appears — note that in (4.9), the part $`c_\delta o_\epsilon (1)`$ is precisely due to the terms of the constant case, and we will see here that the fact that $`b`$ is no longer constant yields terms which are estimated by $`C\delta `$. In the approximation of the equation, the only difference with the constant case is that of course $`(w_\epsilon B)\kappa _\delta `$ is not equal to $`(w_\epsilon \kappa _\delta )B`$. Moreover of course (4.7) is obvious in the constant case, since $`\mathrm{\Pi }_{}`$ commutes with partial derivatives. So we need to deal with those two problems due to the fact that $`b`$ is not constant.

First of all, the difference between $`(w_\epsilon B)\kappa _\delta `$ and $`(w_\epsilon \kappa _\delta )B`$ is small when $`\delta `$ goes to zero, due to the following computation: we have

$$(w_\epsilon B)\kappa _\delta (t,x)(w_\epsilon \kappa _\delta )B(t,x)=w_\epsilon (t,xy)\kappa _\delta (y)(B(xy)B(x))dy$$

hence

$$(w_\epsilon B)\kappa _\delta (t,x)(w_\epsilon \kappa _\delta )B(t,x)=_0^1yB(x\sigma y)w_{t,\epsilon }(xy)\kappa _\delta (y)dyd\sigma .$$

To conclude we need to take the $`L^2`$ norm in $`x`$ of that quantity, and Young’s inequality yields

$$(w_\epsilon B)\kappa _\delta (t,)(w_\epsilon \kappa _\delta )B(t,)_{L^2(\mathrm{\Omega })}C_B||\kappa _\delta _{L^1(\mathrm{\Omega })}w_\epsilon (t)_{L^2(\mathrm{\Omega })}C\delta u_0_{L^2}$$

uniformly in time. The result follows for the first equation in (4.8). The second one is of the same type, since

$$_1(w_{\epsilon ,1}b)+_2(w_{\epsilon ,2}b)=b\text{div}_hw_{\epsilon ,h}+w_{\epsilon ,h}_hb.$$

The term $`(w_{\epsilon ,h}_hb)\kappa _\delta `$ is approximated by $`(w_{\epsilon ,h}\kappa _\delta )_hb`$ exactly as above; to replace the term $`(b\text{div}_hw_{\epsilon ,h})\kappa _\delta `$ by $`b(\text{div}_hw_{\epsilon ,h}\kappa _\delta )`$ we write the same type of computation, with

$$(b\text{div}_hw_{\epsilon ,h})\kappa _\delta (t,x)b(\text{div}_hw_{\epsilon ,h}\kappa _\delta )(t,x)=_3w_{\epsilon ,3}(t,xy)\kappa _\delta (y)(b(xy)b(x))dy$$

hence, since $`\mathrm{\Pi }_{}`$ commutes with $`_3`$,

$$(b\text{div}_hw_{\epsilon ,h})\kappa _\delta b(\text{div}_hw_{\epsilon ,h}\kappa _\delta )_{L^2(𝐑^+\times \mathrm{\Omega })}C_B_3w_{\epsilon ,3}_{L^2(𝐑^+\times \mathrm{\Omega })}\delta C_B\delta u_0_{L^2}.$$

That ends the proof of (4.8).

Now to end the proof of the proposition, we still need to check (4.7). The idea is to use the following estimate, due to the fact that $`w_\epsilon `$ is bounded in $`L_t^2(𝒦)`$ for some compact subspace of $`L^2(\mathrm{\Omega })`$: there is a continuity modulus $`\omega `$ such that

(4.10)
$$y\mathrm{\Omega },w_\epsilon (t,+y)w_\epsilon (t,)_{L^2(\mathrm{\Omega })}\omega (y),\text{uniformly in }t\text{and }\epsilon .$$

Now since $`w_{\epsilon ,h}^\delta =\kappa _\delta w_{\epsilon ,h}`$, we have

$$_hw_{\epsilon ,h}^\delta (x)=_\mathrm{\Omega }\frac{1}{\delta ^4}(_h\kappa )(\frac{y}{\delta })w_{\epsilon ,h}(xy)𝑑y.$$

Since $`_\mathrm{\Omega }_h\kappa (y)𝑑y=0`$, it follows that

$$_hw_{\epsilon ,h}^\delta (x)=_\mathrm{\Omega }\frac{1}{\delta ^4}(_h\kappa )(\frac{y}{\delta })(w_{\epsilon ,h}(xy)w_{\epsilon ,h}(x))𝑑y.$$

Then by (4.10) we find that

$$_hw_{\epsilon ,h}^\delta _{L^2(\mathrm{\Omega })}C\frac{\eta (\delta )}{\delta }_h\kappa _{L^1(\mathrm{\Omega })}$$

where $`\eta (\delta )`$ goes to zero as $`\delta `$ goes to zero. The result is proved.

Now we are ready to prove the following result.

###### Proposition 4.4.

Let $`u^0L^2(\mathrm{\Omega })`$ be a divergence-free vector field. For all $`\epsilon >0`$, denote by $`u_\epsilon `$ a weak solution of (1.1), and by $`w_\epsilon =\mathrm{\Pi }_{}u_\epsilon `$ its projection onto $`\text{Ker}(L)^{}`$. Then

$$\mathrm{\Pi }\left(w_\epsilon \text{rot}w_\epsilon \right)0\text{ in }𝒟^{}(𝐑^+\times (𝒪\mathrm{SS}))\text{ as }\epsilon 0.$$

###### Proof.

Since the result we are looking for is a weak limit, we can restrict our attention to the set $`𝒪`$ where $`b`$ does not vanish, as in $`\mathrm{SS}`$ the result is due to Proposition 3.4. Moreover, to prove the result one can restrict our attention to $`w_\epsilon \text{rot}w_\epsilon dx_3`$ and it is enough to prove that it is proportionnal to $`_hb`$ (up to a small remainder term): taking the scalar product with a function in $`\text{Ker}(L)`$ will then yield automatically zero by definition of $`\text{Ker}(L)`$.

In order to simplify the computations, we shall directly prove the result replacing $`w_{\epsilon ,h}`$ by $`w_{\epsilon ,h}^\delta `$ and $`w_{\epsilon ,3}`$ by $`_3W_{\epsilon ,3}^\delta `$. The difference in the two computations is indeed small when $`\delta `$ is small, uniformly in $`\epsilon `$, exactly as in the proof of (3.7) in the constant case. So writing $`_3W_{\epsilon ,3}^\delta =w_{\epsilon ,3}^\delta `$, and dropping the index $`\delta `$ to simplify, we can perform the following algebraic computations, which will prove the result.

Let us start by recalling that, due to Proposition 4.2, we have

(4.11)
$$w_{\epsilon ,h}=\frac{1}{b}(\epsilon _t\rho _{\epsilon ,h}+r_\epsilon ),$$

where similarly to Section 3 we have noted

$$\rho _{\epsilon ,h}\stackrel{\text{def}}{=}_h^{}W_{\epsilon ,3}w_{\epsilon ,h}^{}.$$

We shall also define

$$\rho _{\epsilon ,3}\stackrel{\text{def}}{=}_1w_{\epsilon ,2}_2w_{\epsilon ,1}.$$

In (4.11) and in the following, the function $`r_\epsilon `$ denotes a remainder term, arbitrarily small in the space $`L_{loc}^2(𝐑^+\times (𝒪\mathrm{SS}))`$. It follows that

$$w_\epsilon \text{rot}w_\epsilon =\left(\begin{array}{c}\frac{1}{b}(\epsilon _t\rho _{\epsilon ,1}+r_\epsilon )\\ \frac{1}{b}(\epsilon _t\rho _{\epsilon ,2}+r_\epsilon )\\ _3W_{\epsilon ,3}\end{array}\right)\left(\begin{array}{c}_3\rho _{\epsilon ,1}\\ _3\rho _{\epsilon ,2}\\ \rho _{\epsilon ,3}\end{array}\right),$$

which implies that

$$w_\epsilon \text{rot}w_\epsilon =\left(\begin{array}{c}\frac{1}{b}(\epsilon _t\rho _{\epsilon ,2}+r_\epsilon )\rho _{\epsilon ,3}_3W_{\epsilon ,3}_3\rho _{\epsilon ,2}\hfill \\ \frac{1}{b}(\epsilon _t\rho _{\epsilon ,1}+r_\epsilon )\rho _{\epsilon ,3}+_3W_{\epsilon ,3}_3\rho _{\epsilon ,1}\hfill \\ \frac{1}{b}(\epsilon _t\rho _{\epsilon ,2}+r_\epsilon )_3\rho _{\epsilon ,1}\frac{1}{b}(\epsilon _t\rho _{\epsilon ,1}+r_\epsilon )_3\rho _{\epsilon ,2}\hfill \end{array}\right).$$

Since the vertical component can be treated exactly as in the constant case, we shall now restrict our attention to the first two components. So calling $`\alpha _h`$ the horizontal components of that vector field, we have after an integration by parts and using the divergence free condition $`_{33}W_{\epsilon ,3}=\text{div}_hw_{\epsilon ,h}`$:

$$\alpha _h𝑑x_3=\frac{1}{b}(\epsilon _t\rho _{\epsilon ,h}^{}+r_\epsilon )\rho _{\epsilon ,3}𝑑x_3\rho _{\epsilon ,h}^{}\text{div}_hw_{\epsilon ,h}𝑑x_3.$$

Now we recall (calling once again generically $`r_\epsilon `$ the small remainder terms) that

$$\epsilon _t\rho _{\epsilon ,3}+\text{div}_h(bw_{\epsilon ,h})=r_\epsilon $$

so

$$\text{div}_hw_{\epsilon ,h}=\frac{1}{b}(\epsilon _t\rho _{\epsilon ,3}w_{\epsilon ,h}_hb+r_\epsilon ).$$

It follows that

$`{\displaystyle \alpha _h𝑑x_3}`$ $`=`$ $`{\displaystyle \frac{1}{b}}{\displaystyle \epsilon _t\rho _{\epsilon ,h}^{}\rho _{\epsilon ,3}dx_3}{\displaystyle \frac{1}{b}}{\displaystyle \rho _{\epsilon ,h}^{}\epsilon _t\rho _{\epsilon ,3}dx_3}+{\displaystyle \frac{1}{b}}{\displaystyle \rho _{\epsilon ,h}^{}w_{\epsilon ,h}_hbdx_3}+{\displaystyle \stackrel{~}{r}_\epsilon 𝑑x_3}`$
$`=`$ $`{\displaystyle \frac{1}{b}}{\displaystyle \epsilon _t(\rho _{\epsilon ,h}^{}\rho _{\epsilon ,3})dx_3}+{\displaystyle \frac{1}{b}}{\displaystyle \rho _{\epsilon ,h}^{}w_{\epsilon ,h}_hbdx_3}+{\displaystyle \stackrel{~}{r}_\epsilon 𝑑x_3},`$
where now $`\stackrel{~}{r}_\epsilon `$ denotes generically the product of $`r_\epsilon `$ by a component of $`\rho _\epsilon `$. But $`\rho _{\epsilon ,3}`$ is a combination of derivatives of $`w_\epsilon `$ whereas $`\rho _{\epsilon ,h}`$ is a combination of components of $`w_\epsilon `$. A product of the type $`\rho _{\epsilon ,h}r_\epsilon `$ clearly goes to zero in $`𝒟^{}(𝐑^+\times (𝒪\mathrm{SS}))`$. For the term $`r_\epsilon \rho _{\epsilon ,3}`$, one uses result (4.7) and (4.9) to infer that for any subset $`K`$ of $`𝒪\mathrm{SS}`$,

$`r_\epsilon \rho _{\epsilon ,3}_{L^1(K)}`$ $``$ $`r_\epsilon _{L^2(K)}\rho _{\epsilon ,3}_{L^2(K)}`$
$``$ $`(c_\delta o_\epsilon (1)+C\delta )_hw_{\epsilon ,h}^\delta _{L^2(K)}0,`$
as $`\epsilon `$ followed by $`\delta `$ go to zero. So from now on $`\stackrel{~}{r}_\epsilon `$ will denote generically a term going to zero in $`𝒟^{}(𝐑^+\times (𝒪\mathrm{SS}))`$.

Recalling (4.11) we get

$$\alpha _h𝑑x_3=\frac{1}{b}\epsilon _t(\rho _{\epsilon ,h}^{}\rho _{\epsilon ,3})dx_3\frac{1}{b^2}\rho _{\epsilon ,h}^{}\epsilon _t\rho _{\epsilon ,h}_hbdx_3+\stackrel{~}{r}_\epsilon 𝑑x_3.$$

In particular we have

$$\alpha _1𝑑x_3=\frac{1}{b}\epsilon _t(\rho _{\epsilon ,2}\rho _{\epsilon ,3})dx_3\frac{1}{b^2}\rho _{\epsilon ,2}\epsilon _t\rho _{\epsilon ,1}_1bdx_3\frac{1}{2b^2}\epsilon _t(\rho _{\epsilon ,2}^2_2b)dx_3+\stackrel{~}{r}_\epsilon 𝑑x_3.$$

Similarly

$$\alpha _2𝑑x_3=\frac{1}{b}\epsilon _t(\rho _{\epsilon ,1}\rho _{\epsilon ,3})dx_3+\frac{1}{b^2}\rho _{\epsilon ,1}\epsilon _t\rho _{\epsilon ,2}_2bdx_3+\frac{1}{2b^2}\epsilon _t(\rho _{\epsilon ,1}^2_1b)dx_3+\stackrel{~}{r}_\epsilon 𝑑x_3.$$

But one can also write

$$\alpha _2𝑑x_3=\frac{1}{b}\epsilon _t(\rho _{\epsilon ,1}\rho _{\epsilon ,3}+\frac{1}{2b}\rho _{\epsilon ,1}^2_1b+\frac{1}{b}\rho _{\epsilon ,1}\rho _{\epsilon ,2}_2b)dx_3\frac{1}{b^2}\rho _{\epsilon ,2}\epsilon _t\rho _{\epsilon ,1}_2bdx_3+\stackrel{~}{r}_\epsilon 𝑑x_3.$$

It follows that up to full derivatives of the type $`\epsilon _t`$, $`{\displaystyle \alpha _h𝑑x_3}`$ is equal to

(4.12)
$$\left(\frac{1}{b^2}\rho _{\epsilon ,2}\epsilon _t\rho _{\epsilon ,1}dx_3\right)_hb+\stackrel{~}{r}_\epsilon 𝑑x_3.$$

Now the proof is almost finished: we recall that we want to take the projector of the term $`{\displaystyle w_\epsilon }\text{rot}w_\epsilon dx_3`$ onto the kernel of $`L`$ restricted to the set $`𝒪`$. Recalling that $`\text{Ker}L`$ is made of vector fields of the type $`_h^{}\phi `$ with $`_hb_h^{}\phi =0`$, we have obviously

$$\mathrm{\Pi }\left(\frac{1}{b^2}\rho _{\epsilon ,2}\epsilon _t\rho _{\epsilon ,1}dx_3\right)_hb=0.$$

The result is therefore proved for the horizontal component of $`{\displaystyle w_\epsilon }\text{rot}w_\epsilon dx_3`$. The third component is identical to the constant case, so the proposition is proved. ∎

This ends the proof of (4.3), hence of Lemma 4.1. ∎

In the following we shall denote by $`()`$ the limiting system:

$$()\{\begin{array}{c}_t\overline{u}\nu \mathrm{\Pi }\mathrm{\Delta }_h\overline{u}+\mathrm{\Pi }(\overline{u}_h\overline{u})=0\text{in}𝒟^{}(𝐑^+\times (𝒪\mathrm{SS})),\mathrm{\Pi }_{}\overline{u}=0\\ \overline{u}_{|t=0}=\mathrm{\Pi }u^0.\end{array}$$

The existence of solutions to this system is easy to prove, as it is of the same form as a 2D Navier–Stokes equation: the only point which might be a problem is that $`\mathrm{\Pi }`$ does not in general commute with the Friedrich frequency truncation $`J_n`$, recalled in Section 2: however approximating the system by for instance

$$_tJ_n\overline{u}_n\nu J_n\mathrm{\Pi }J_n\mathrm{\Delta }_hJ_n\overline{u}_n+J_n\mathrm{\Pi }(J_n\overline{u}_n_hJ_n\overline{u}_n)=0,J_n\mathrm{\Pi }_{}J_n\overline{u}_n=0$$

will do the job. In any case in the next two sections we shall give precise formulations for the solution of $`()`$: in $`\mathrm{SS}`$ $`\overline{u}_h`$ is the unique solution of $`(NS2D)`$, and outside $`\mathrm{SS}`$ it is the solution of a heat equation. The uniqueness of $`\overline{u}`$ implies in particular that the convergence holds for the whole sequence $`u_\epsilon `$ and not only for a subsequence.

As noted earlier in this section, the third component of $`\mathrm{\Pi }u`$ for any vector field $`u`$ is the same whether $`b`$ is constant or not and is simply the vertical average of $`u_3`$. It follows that the third component of this equation is simply

$$_t\overline{u}_3\nu \mathrm{\Delta }_h\overline{u}_3+\overline{u}_h_h\overline{u}_3=0\text{in}𝐑^+\times (𝒪\mathrm{SS}),\overline{u}_{3|t=0}=u_{3|t=0}𝑑x_3.$$

Now all the work consists in determining $`\overline{u}_h`$. We shall consider separately the vector field on $`𝒪`$ and on $`\mathrm{SS}`$, which is the object of the two following sections; so in those sections, our attention will be restricted to the horizontal component $`\overline{u}_h`$.

### 4.2. The averaged equation on $`𝒪`$

We shall prove the following result, which yields the part of Theorem 3 which lies in $`𝒪`$.

###### Proposition 4.5.

Let $`\overline{u}`$ be a vector field satisfying $`()`$ with

$$\overline{u}L^{\mathrm{}}(𝐑^+,\text{H})L^2(𝐑^+,\dot{H}^1)\text{and}\overline{u}\text{Ker}(L).$$

Then the vector field $`\overline{u}_h`$ satisfies the following heat equation:

$`_t\overline{u}_h\nu \mathrm{\Pi }\mathrm{\Delta }_h\overline{u}_h=0\text{in}𝐑^+\times 𝒪`$
$`\overline{u}_{h|t=0}=\mathrm{\Pi }u_{h,0}.`$
###### Proof.

The function $`\overline{u}_h`$ satisfies

$$_t\overline{u}_h\nu \mathrm{\Pi }\mathrm{\Delta }_h\overline{u}_h+\mathrm{\Pi }(\overline{u}_h_h\overline{u}_h)=0,$$

and since $`\overline{u}`$ is in $`\text{Ker}(L)`$, we have

$$\left(\mathrm{\Pi }\mathrm{\Delta }_h\overline{u}_h|\overline{u}_h\right)_{L^2(𝒪)}=\left(\mathrm{\Delta }_h\overline{u}_h|\overline{u}_h\right)_{L^2(𝒪)}=_h\overline{u}_h_{L^2(𝒪)}^2.$$

We note that $`\overline{u}_h`$ is equal to zero on the boundary of $`𝒪`$, since it is a multiple of $`_h^{}b`$.

So to prove the proposition, the only point we need to check is that

(4.13)
$$\mathrm{\Phi }\text{Ker}(L)\dot{H}^s(𝒪),s>1,\left(\overline{u}_h_h\overline{u}_h|\mathrm{\Phi }_h\right)_{L^2(𝒪)}=0.$$

By definition of $`\text{Ker}(L)`$, we have

$$\mathrm{\Phi }_h_hb=\overline{u}_h_hb=0,$$

from which we infer that

(4.14)
$$\mathrm{\Phi }_h\overline{u}_h=0.$$

Now as in (3.2) we can write

$$\mathrm{\Pi }(\overline{u}_h_h\overline{u}_h)=\mathrm{\Pi }(\overline{u}_h\text{rot}\overline{u}_h)$$

hence

$$\left(\overline{u}_h_h\overline{u}_h|\mathrm{\Phi }_h\right)_{L^2(𝒪)}=\left(\overline{u}_h\mathrm{\Phi }_h|\text{rot}\overline{u}_h\right)_{L^2(𝒪)}.$$

Finally Identity (4.14) yields (4.13), and Proposition 4.5 is proved. ∎

### 4.3. The 2D Navier-Stokes limit on $`\mathrm{SS}`$

In this section we shall analyse the equation satisfied by the limit on $`\mathrm{SS}`$, that is to say in the regions where $`b`$ is a constant.

###### Proposition 4.6.

Let $`\overline{u}`$ be a vector field satisfying $`()`$ with

$$\overline{u}L^{\mathrm{}}(𝐑^+,\text{H})L^2(𝐑^+,\dot{H}^1)\text{and}\overline{u}\text{Ker}(L).$$

Then the vector field $`\overline{u}_h`$ satisfies the two–dimensional Navier–Stokes equations in $`\mathrm{SS}`$, with homogeneous Dirichlet boundary conditions:

$$\begin{array}{c}_t\overline{u}_h\nu \mathrm{\Delta }_h\overline{u}_h+\overline{u}_h_h\overline{u}_h=_hp\text{in}𝐑^+\times \mathrm{SS}\\ \overline{u}_{h|t=0}=u_{h,0}𝑑x_3\\ \overline{u}_{h|\mathrm{SS}}=0.\end{array}$$

###### Proof.

First let us recall why the equation on $`\overline{u}_h`$ is the two–dimensional Navier–Stokes equation: we simply consider the weak formulation of the original rotating fluid equations and take its limit, by integrating against a test function $`\mathrm{\Phi }`$, divergence free and compactly supported in $`𝐑^+\times \mathrm{SS}`$. The weak formulation is as follows:

(4.15)
$$_{𝐑^+}_\mathrm{\Omega }\left(u_\epsilon _t\mathrm{\Phi }+\nu u_\epsilon \mathrm{\Phi }u_\epsilon u_\epsilon \mathrm{\Phi }\right)𝑑x𝑑t=_\mathrm{\Omega }u^0\mathrm{\Phi }_{|t=0}𝑑x.$$

Then taking the limit as $`\epsilon `$ goes to zero in (4.15) yields, due to Lemma 4.1,

$$\left(\overline{u}_t\mathrm{\Phi }+\nu \overline{u}\mathrm{\Phi }\overline{u}\overline{u}\mathrm{\Phi }\right)𝑑x=_\mathrm{\Omega }u^0\mathrm{\Phi }_{|t=0}𝑑x$$

where we have noticed that on $`\mathrm{SS}`$, we have

$$\mathrm{\Pi }(\overline{u}\overline{u})=(\overline{u}\overline{u}).$$

Now recalling that $`\overline{u}`$ only depends on the horizontal variable, we deduce the expected equation on $`\overline{u}`$, up to the boundary terms. To get the boundary terms, we simply recall that the limit $`\overline{u}`$ is in $`L_{loc}^2(𝐑^+,H^1(\mathrm{\Omega }_h))`$ hence cannot have a jump on the boundary of $`\mathrm{SS}`$. Then we notice that $`\mathrm{SS}𝒪`$, simply because if $`x\mathrm{SS}`$, then $`b(x)=0`$. So the result follows directly: the boundary condition is a homogeneous Dirichlet boundary condition.

Theorem 3 is proved.


text:
# Where does a cohesive granular heap break?

## I introduction

The economic impact of particle processing is enormous. Better methods for the design and synthesis of unit operations involving divided solids have been identified as a critical need, especially for pharmaceutical, agrochemicals and specialty chemicals. One of the important features in industrial processes is to better understand the transition between a static granular medium and the avalanche process. Avalanches in non-cohesive granular media have been extensively investigated Jaeger et al. (1988); Rajchenbach (1990). A common characteristic of these studies is that granular motion occurs in a relative thin boundary layer (around ten grains) at the surface Jaeger and Nagel (1992) independent of the size of the sample. On the other hand, recent experiments have explored the relatively new subject of cohesive granular media such as“humid granular” media. The presence of capillary bridges between the grains generate adhesive forces which strongly affect the stability of the samples Albert et al. (1997, 1999); Hornbaker et al. (1997); Bocquet et al. (1998); Ursini et al. (2002); Fraysse et al. (1999); Valverde et al. (2000); Quintanilla et al. (2001). These capillary bridges either originate in small amounts of added fluid, and this situations corresponds to “lightly wet granular media”or either are created by a condensable vapor in the atmosphere Bocquet et al. (1998); Ursini et al. (2002); Fraysse et al. (1999). In the latter case the vapor is in chemical equilibrium with the interstitial liquid bridges. These experiments concern the so-called “moist granular media”.

As shown by all these experiments, the behavior of humid granular media strongly departs from the dry materials and many new features appear. In a series of experiments on “lightly wet” granular media, Tegzes et al. Tegzes et al. (1999), measured the angle of repose of a granular media (by the draining crater method) as a function of the liquid content and of the size of the system. They observed an increase of the angle of repose with the liquid content. More precisely, they obtained three regimes: i) at low liquid content, the “granular regime”, the avalanche occurs at the surface and the angle of repose does not depend on the size of the system; ii) at intermediate liquid content, the “correlated regime”, the avalanche takes place more deeply in the heap, the angle of repose depends on the size of the system; iii) at high liquid content, a kind of “plastic flow” is observed. This various regimes are not accounted for by previous theoretical analysis.

In numerical and analytical modelizations, cohesion is usually taken into account as a constant adhesion force. As an example, Forsyth et al. Forsyth et al. (2001) and Valverde et al. Valverde et al. (2000) used a constant van der Waals force to interpret their experimental results, Olivi-Tran et al. Olivi-Tran et al. (2002) used a constant capillary force or a constant force due to a solid bridge in particle dynamics (PD) simulations, Albert et al. Albert et al. (1997) used a constant capillary force in a geometrical model and Nase et al. Nase et al. (2001), Mikami et al. Mikami et al. (1998) used a constant capillary force in PD simulations.

In this paper we show that the assumption of a constant cohesion is not justified in many practical cases. We then show that the dependence of the cohesion on normal stress, usually omitted in the literature, is a key point in determining the stability of a cohesive granular heap. This effect leads to a dependence of the stability on the system size, as found experimentally. To this end, we shall use a continuum analysis Nedderman (1992) to study the stability of cohesive sandpiles taking into account a non constant cohesive stress in the pile. The paper is divided in two parts. We briefly first show that in a many situations, the adhesion force between the grains depends on the normal stress. Then, in the second part, we show how this relationship affects the stability of a heap and the localization of the slip plane upon failure.

## II Contact forces between grains

The adhesion force between two ideal spheres due to solid-solid interactions has been calculated in the case of an elastic contact a long time ago Johnson et al. (1971); Derjaguin et al. (1983); Maugis (1992):

$$F_{SS}=f\pi \gamma R$$
(1)
where $`R`$ is the sphere radius, $`\gamma `$ is the solid/intersticial medium surface tension ($`\gamma =\gamma _{SG}`$ if the spheres are in a gas atmosphere, and $`\gamma =\gamma _{SL}`$ if they are immersed in a liquid), and $`f`$ is a numerical factor between 1.5 and 2. The lower value of $`f=1.5`$ corresponds to the Johnson-Kendall-Roberts calculation (JKR) valid when the attractive forces are strong enough to deform the spheres surfaces. The higher value $`f=2`$ corresponds to the Derjaguin-Muller-Toporov (DMT) calculation, which is valid in the opposite limit of weak attraction and rigid solids. The main hypothesis in obtaining Eq. 1 is that surfaces are ideal. An obvious feature of this result is that adhesion does not depend on the pressure applied on the surfaces before to separate them Israelachvili (1992). However measured values of the adhesion force generally departs from the predicted one Restagno et al. (2002). As we briefly discuss in the next section, the effect of surface roughness leads to qualitatively different results, involving a dependence of the adhesion force on the normal load.

### II.1 Pressure depending adhesion between two grains

Restagno et al. Restagno et al. (2002) have recently measured the adhesion force between two weakly rough surfaces (a sphere and a plane) using a surface forces apparatus. The “pull-off” force $`F_{\text{adh}}`$ to separate two surfaces has been found to depend on the a loading normal force, $`F`$, with a scaling: $`F_{\text{adh}}F^{1/3}`$.

The key-point to interpret this result is the roughness of the surfaces as pointed out 50 years ago by MacFarlane and Tabor Bowden and Tabor (1950): the surfaces are in contact only at the tip of their asperities, and the real contact area $`A_\text{r}`$ is much smaller than the apparent contact area $`A_{\text{app}}`$ (Fig. 1. To obtain the adhesion force they used Eq. 1 with an effective surface tension $`\gamma _{\text{eff}}=\gamma A_\text{r}/A_{\text{app}}`$. On one hand, different ways can be used to calculate $`A_\text{r}`$. The first way is to consider that the pressure in the molecular contact area is so high than it can reach the yield stress of the material Bowden and Tabor (1950). An other way has been proposed by Greenwood Greenwood (1967) considering the statistical distribution of heights of the surface asperities. This two models lead to a simple proportionality between the real contact area $`A_r`$ and the normal load:

$$A_\text{r}F$$
(2)
On the other hand, the calculation of the apparent contact is the elastic contact area calculated by Hertz Landau and Lifchitz (1982) : $`A_{\text{app}}\left(FR/E^{}\right)^{2/3}`$, with $`E^{}=E/[2(1\nu ^2)]`$ the reduced Young modulus of the material, defined in terms of the Young modulus $`E`$ and the Poisson modulus $`\nu `$ of the material. Gathering results, one obtains eventually $`\gamma _{\text{eff}}F^{1/3}`$, so that Eq. 1 is compatible with the observed experimental scaling.

In the previous experiments, the cohesion between two rough surfaces due to the direct interaction between the surfaces has been considered. However, adhesion might also originate in the presence of a cohesive ”binder” in between the surfaces. In many practical situations, the presence of liquid bridges connecting the grains leads to strong adhesive forces between the grains Albert et al. (1997, 1999); Hornbaker et al. (1997); Bocquet et al. (1998); Fraysse et al. (1999); Ursini et al. (2002). When the material is placed in an atmosphere of a condensable vapor Bocquet et al. (1998); Ursini et al. (2002), the size of the liquid bridge is fixed via chemical equilibrium by the value of external humidity. The size of such bridges is very small, around a few nanometers, and one expects that the liquid bridges will condense only in the regions close to real contact: One expects therefore that $`F_{\text{coh}}A_\text{r}`$, the area of real contact area (see e.g. Bocquet et al. (1998) for a more complete discussion). Since as discussed above in Eq. 2, $`A_\text{r}`$ is proportional to the external load, one gets eventually: $`F_{\text{adh}}F`$. This picture is however valid when the size $`r`$ of the bridges is lower than a typical scale of roughness, $`\mathrm{}_R`$: In the opposite case, $`r>\mathrm{}_R`$, roughness is not pertinent any more and the capillary force reduces to a contant value $`2\pi \gamma R`$, with $`R`$ the radius of the beads, and $`\gamma `$ the liquid vapor surface tension.

To sum up, the adhesion force between two surfaces can be written:

$$F_{\text{adh}}F^n$$
(3)
with various expected exponents $`n`$. For humid granular materials, a transition from from $`n=1`$ to $`n=0`$ is expected, depending on the relative size of the capillary bridges (i.e.  on the added amount of liquid) compared to roughness : For small amount of liquid, the value $`n=1`$ is expected, while large amounts of liquid correspond to $`n=0`$. On the other hand, for dry or liquid saturated materials in which cohesion is expected to originate in direct (e.g. van der Waals) interactions, a value $`n=1/3`$ is predicted.

We now derive the local cohesion stress from the adhesion force.

### II.2 from adhesion force to local cohesion stress

The problem of the derivation of the interparticle forces $`F`$ from the bulk stress $`\sigma `$ is an old problem. Quintanilla et al.Quintanilla et al. (2001) have shown, comparing direct adhesion between particules using an AFM and cohesion stress measured by determining the tensile stress of a material, that adhesion increases with the normal load in the material and that a correlation between the microscopic contact force and a macroscopic stress measurement can be performed with a good precision.

For a system of hard monodisperse spherical particles with a random isotropic packing, a continuum-like picture leads to a simple linear relationship between force $`F`$ and stress $`\sigma `$ NOTE  :

$$\sigma =\frac{\varphi k}{4\pi R^2}F$$
(4)
where $`R`$ is the radius of the particles, $`k`$ is the coordination number, which is defined as the average number of contacts per particle, and $`\varphi `$ is the volume fraction. Such a relationship has been widely used in the literature for powders or wet granular media Mason et al. (1999); Quintanilla et al. (2001). From Eq. 4, we can deduce the relationship between the adhesion force $`F_{\text{adh}}`$ and the adhesive stress $`c_{\text{adh}}`$ in the granular medium:

$$c_{\text{adh}}=\frac{\varphi k}{4\pi R^2}F_{\text{adh}}\frac{F_{\text{adh}}}{R^2}$$
(5)
With the same arguments, we can deduce the relationship between the normal load and the normal stress:

$$\sigma =\frac{\varphi k}{4\pi R^2}F\frac{F}{R^2}$$
(6)
Assuming a general relationship of the form of Eq. 3, one gets therefore an adhesive stress which depends on the normal stress as shown in Fig. 2:

$$c_{\text{adh}}=c_0\left(\frac{\sigma }{\sigma _0}\right)^n$$
(7)
where $`n`$, $`\sigma _0`$ and $`c_0`$ characterize the properties of the material. This relationship is our general starting point to study the stability of a cohesive heap.

## III Failure of a cohesive heap

The most simple approach in the continuum limit to study the failure of a granular material is the plastic theory of Mohr-Coulomb. The validity of the continuum approach has been widely discussed. Nevertheless continuum approaches allow to predict the Green’s function of a granular layer Reydellet and Clement (2001) and is widely used in soil mechanics Nedderman (1992); Mason et al. (1999). We consider the stability of a cohesive granular heap, characterized by an adhesive stress depending on the normal stress, as given in Eq. 7 in which we suppose the Mohr-Coulomb criterion valid. Our aim is to locate in such a material the ”slip plane”, where failure occurs. The basis of our analysis is the Coulomb criterion: a granular material is stable if for each surface inside the material the following inequality is obeyed:

$$\tau \mu \sigma $$
(8)
where $`\tau `$ is the shear stress and $`\sigma `$ is the normal stress. For a cohesive material, we assume that this condition can be generalized by adding in Eq. 8 the adhesive stress $`c_{\text{adh}}`$, as defined in the previous section, to the normal stress $`\sigma `$ Bocquet et al. (1998).

We consider the geometry depicted in Fig. 3: A granular heap, with height $`H`$, makes an angle $`\theta `$ with the horizontal. The heap is supposed to be invariant in the direction perpendicular to the figure. We follow and generalize the approach given in Ref. Nedderman (1992) to the case of a cohesive material with the law of cohesion given in Eq. 7. In the simplified description of Ref. Nedderman (1992), the slip surface is assumed to be planar, making an angle $`\alpha `$ against the horizontal.

The force balance in the direction perpendicular and parallel to the slip plane leads to the following conditions:

$$F=P\mathrm{sin}\alpha N=P\mathrm{cos}\alpha +F_{\text{adh}}$$
(9)
where $`P`$ is the weight of the granular material above the slip plane; $`c=F/S_{\text{slip}}`$ and $`\sigma +c_{\text{adh}}=N/S_{\text{slip}}`$ are the normal and tangential stress along the possible slip plane with area $`S_{\text{slip}}`$. In these variables, the stability criterion, Eq. 8 writes $`F<\mu N`$. A simple geometric calculation leads to the relationships Nedderman (1992):

$`P=`$ $`{\displaystyle \frac{1}{2}}\rho gLH^2\left({\displaystyle \frac{1}{\mathrm{tan}\alpha }}{\displaystyle \frac{1}{\mathrm{tan}\theta }}\right)`$ (10)
$`S_{\text{slip}}=`$ $`{\displaystyle \frac{HL}{\mathrm{sin}\alpha }}`$
with $`L`$ the length of the heap in the invariant direction (perpendicular to the figure in Fig. 3). Gathering these results, one gets the following stability criterion:

$$\frac{1}{2}\rho gH^2\left(\frac{1}{\mathrm{tan}\alpha }\frac{1}{\mathrm{tan}\theta }\right)(\mathrm{sin}\alpha \mu \mathrm{cos}\alpha )\mu \frac{c_{\text{adh}}H}{\mathrm{sin}\alpha }$$
(11)
Introducing the angle $`\theta _0`$ defined as $`\mathrm{tan}\theta _0=\mu `$, this relation can be adequately rewritten as:

$$\mathrm{sin}(\theta \alpha )\mathrm{sin}(\alpha \theta _0)2\mathrm{sin}\theta _0\mathrm{sin}\theta \frac{c_{\text{adh}}}{\rho gH}$$
(12)
Now, one has to introduce the dependence of the cohesive stress $`c_{\text{adh}}`$ as a function of the normal stress, as discussed in the previous section. We shall use the general expression given by Eq. 7. Using $`\sigma =P\mathrm{cos}\alpha /S_{\text{slip}}`$, one gets after some algebra:

$$c_{\text{adh}}=c_0\left(\frac{\rho gH}{2\sigma _0}\right)^n\left(\frac{\mathrm{sin}(\theta \alpha )\mathrm{cos}\alpha }{\mathrm{sin}\theta }\right)^n$$
(13)
Introducing the cohesion parameter $`\zeta _{\text{coh}}`$ defined as:

$$\zeta _{\text{coh}}=\frac{c_0}{\rho gH}\left(\frac{\rho gH}{2\sigma _0}\right)^n,$$
(14)
and the function $`f[\alpha ]`$ defined as:

$$f[\alpha ]=\frac{(\mathrm{sin}(\theta \alpha ))^{1n}}{(\mathrm{cos}\alpha )^n}\mathrm{sin}(\alpha \theta _0)$$
(15)
Eq. 12 rewrites:

$$f[\alpha ]\zeta _{\text{coh}}\mathrm{sin}\theta _0(\mathrm{sin}\theta )^{n1}$$
(16)
The heap is therefore stable at an angle $`\theta `$ if this inequality is verified for all possible $`\alpha `$ (with $`0\alpha \theta `$). In the opposite case when this inequality is not verified for some values of $`\alpha `$, the slip plane corresponds to the value of $`\alpha `$ for which this inequality is “first” violated. Note that for a non-cohesive material, $`\zeta _{\text{coh}}=0`$, the maximum angle of stability is simply $`\theta _0`$.

### III.1 $`n1`$ case

First if $`\theta \theta _0`$, the function $`f[\alpha ]`$ is always negative and the stability condition is trivially verified: For $`\theta <\theta _0`$ the heap is always stable.

The opposite case $`\theta >\theta _0`$ is more complex. A typical plot for this function is then given in figure 4.

This function does exhibit a maximum for a value $`\alpha _\text{m}`$, verifying $`f/\alpha (\alpha _\text{m})=0`$, and depending of the heap angle $`\theta `$. According to criterion 16, the heap is stable if $`f[\alpha _\text{m}]\zeta _{\text{coh}}\mathrm{sin}\theta _0(\mathrm{sin}\theta )^{n1}`$. The maximum angle of stability $`\theta _\text{m}`$ thus is solution of the implicit equation:

$$f[\alpha _\text{m}]=\zeta _{\text{coh}}\mathrm{sin}\theta _0(\mathrm{sin}\theta _\text{m})^{n1}.$$
(17)
A straightforward calculation shows that $`\alpha _\text{m}`$ obeys the following relationship:

$$\mathrm{sin}(2\alpha _\text{m}\theta _0\theta )\mathrm{cos}\alpha _\text{m}=n\mathrm{cos}(2\alpha _\text{m}\theta )\mathrm{sin}(\alpha _\text{m}\theta _0)$$
(18)
It is interesting to consider the two limiting cases $`n=1`$ and $`n=0`$.

* for $`n=1`$, the solution of Eq. 18 is: $`\alpha _\text{m}=\theta `$; and $`\theta _\text{m}`$ is found to obey:

$$\frac{\mathrm{sin}(\theta _\text{m}\theta _0)}{\mathrm{cos}\theta _\text{m}\mathrm{sin}\theta _0}=\zeta _{\text{coh}}$$
(19)
This relationship can be rewritten in a more explicit form as :

$$\mathrm{tan}(\theta _\text{m})=\mu (1+\zeta _{\text{coh}})=\mu \left(1+\frac{c_0}{2\sigma _0}\right)$$
(20)
with $`\mu =\mathrm{tan}\theta _0`$. This case presents two important characteristics: (i) failure occurs here at the surface of the heap, as emphasized by the relationship $`\alpha _\text{m}=\theta `$; (ii) whatever the heap height $`H`$, there is always an effect of cohesion on the maximum angle of stability. This regime has been called the “granular regime” according to Tegzes experiments in which a cohesion effect on the avalanche angle is found without dependancy on the heap size. In Tegzes experiments, this regime is found at low liquid content. In section II, we showed that the value $`n=1`$ is indeed expected when cohesion results of small capillary bridges between surfaces asperities.
* for $`n=0`$ (“classic adhesion regime” i.e. constant adhesion force), the solution of Eq. 18 is: $`\alpha _\text{m}=(\theta _0+\theta )/2`$, and $`\theta _\text{m}`$ verifies the equation:

$$\frac{1\mathrm{cos}(\theta _\text{m}\theta _0)}{2\mathrm{sin}\theta \mathrm{sin}\theta _0}=\zeta _{\text{coh}}=\frac{c_0}{\rho gH}$$
(21)
This is a classical result, as obtained e.g. in Ref. Nedderman (1992). It is important to note that in this case, (i) the heap fails deep inside the material, since $`\alpha _\text{m}<\theta _\text{m}`$; (ii) the maximum angle of stability $`\theta _\text{m}`$ depends on the height of the heap $`H`$. Such a dependence is in fact expected when one realizes that for $`n=0`$, a ”capillary length scale” can be defined on dimensional grounds : $`\mathrm{}_{\text{cap}}=c_0/\rho g`$ (see eg Eq. 14). This capillary legnth gives the size of a macroscopic element of the granular medium whose weight equals the adhesion force which act on it. For instance it maybe the maximum size of a powder aggregat which can remain stuck under a horizontal surface. Eq. 21 shows that the maximum angle of stability of the heap reduces to $`\theta _0`$ when the size of the heap is very large compared to the “capillary length”.
* For $`n`$ in between these two limiting cases, one has to solve numerically Eqs. 18 and 17. $`\alpha _\text{m}`$ lies in between $`\theta `$ and $`(\theta +\theta _0)/2`$. We show in Fig. 5 a typical result for the dependence of $`\theta _\text{m}`$ and $`\alpha _\text{m}`$ as a function of $`\zeta _{\text{coh}}`$ for $`n=1/3`$.

The dependence of $`\theta _\text{m}`$ in the limit of small cohesion can be computed analytically. For $`\zeta _{\text{coh}}=0`$, one has $`\alpha _\text{m}=\theta _\text{m}=\theta _0`$ and one may expand the angles around this values for small $`\zeta _{\text{coh}}`$: $`\alpha _\text{m}=\theta _0+\delta \alpha `$ and $`\theta _\text{m}=\theta _0+\delta \theta `$. First linearizing Eq. 18), one gets the relationship $`\delta \alpha =\delta \theta /(2n)`$. Introducing this condition into Eq. 17), one gets eventually:

$$\theta _\text{m}\theta _0=\gamma \zeta _{\text{coh}}^{\frac{1}{2n}}$$
(22)
with $`\gamma =(2n)/(1n)^{(1n)/(2n)}\mathrm{sin}\theta _0(\mathrm{cos}\theta _0)^{n/2n}`$. Inserting the $`H`$ dependence of $`\zeta _{\text{coh}}`$, as defined in Eq. 14, one gets therefore

$$\theta _\text{m}\theta _0H^{\frac{n1}{2n}}$$
(23)
This shows in particular that the dependence of the maximum angle on the heap height is a direct measure of the power $`n`$, characterizing the normal stress dependence of the cohesion stress. This regime is called the “size dependent avalanche” regime.

### III.2 $`n>1`$ case

In this case, it is easy to show that whatever the cohesion, the maximum angle of stability is $`\theta _0`$. This regime is called the “apparent non-cohesive” regime. First as in the previous case, when $`\theta <\theta _0`$, the condition of stability, Eq. 16), is trivially obeyed since the function $`f[\alpha ]`$ is negative. Now, for $`\theta >\theta _0`$, $`f[\alpha ]`$ goes to infinity when $`\alpha \theta `$. The heap is therefore always unstable, whatever the cohesion. As a result, $`\theta =\theta _0`$ is the maximum angle of stability, independently of the cohesion $`\zeta _{\text{coh}}`$

## IV discussion

In this paper we have studied the effect of an adhesive stress on the localization of the slip plane in an avalanche process. More precisely, we have shown that in a lot of practical case, the adhesive stress in a material depends on the normal stress. We have calculated, using a Mohr-Coulomb analysis, the internal angle of slip $`\alpha _\text{m}`$.

A few conclusions can be drawn from these results:

* First, the location of the slip plane, here defined through $`\alpha _\text{m}`$, does also strongly depend on this functional dependence: While for a constant cohesion (i.e. independent of the normal stress) the heap slips deep inside, only surface slip is expected when cohesion is linearly related to normal stress.
* Second the maximum static angle, $`\theta _\text{m}`$, does depend on cohesion via the dimensionless parameter, $`\zeta _{\text{coh}}`$, defined as $`\zeta _{\text{coh}}=\frac{c_0}{\rho gH}\left(\frac{\rho gH}{2\sigma _0}\right)^n`$. An important remark is that for $`n1`$, the maximum static angle depends on the height of the heap, $`H`$. On the other hand, for the specific value $`n=1`$ \- which is expected in some physical situations (see discussion in section II.1)- this dependence disappears and $`\theta _\text{m}`$ is an intrinsic property of the material, independently of the geometry.
* Eventually, a change of regime in the cohesion, e.g. a change from $`n=1`$ to $`n=0`$ as discussed in section II.1, will not only modify the dependence of $`\theta _\text{m}`$ on cohesion, but more dramatically, it will change the localization of the slip plane: For example, while for $`n=1`$ slip occurs at the surface of the the heap ($`\alpha _\text{m}=\theta _\text{m}`$), it will fracture in the interior of the heap for $`n=0`$ ($`\alpha _\text{m}<\theta _\text{m}`$). In other words, any change of slip behaviour reflects a transition of cohesion regime.

This features provide a framework to explain the results of Tegzes et al.. In view of the previous results, it seems natural to interprete the different regimes observed in the experiment, as an indication of a change in the cohesion regime. At low liquid content, the avalanche takes place at the surface of the heap and the angle increases linearly with the adhesion. With our approach, this would correspond to the $`n=1`$ regime, which is indeed expected at low liquid content (where the capillary bridges do not fill the full interstitial space between the grains). At higher liquid content, a change to a regime with $`n<1`$ is expected. This might explain the various observed experimental features, such as the size dependence.

As a conclusion, we hope that this work will motivate further experimental investigation on the stability of cohesive granular materials. The present results suggest that a careful determination of the stability properties and of the failure plane localization yields information on the cohesion forces between grains. In particular the heap height dependence of the maximum angle of stability should provide a direct measure of the cohesion properties. Another interesting geometry is the cylindrical bunker (”Janssen’s problem”), in which cohesion effects, as discussed here, should play a particularly important role.

###### Acknowledgements.

We thank the Région Rhône-Alpes for her financial support (Programme Emergence 02 1892601).