M Bonnivard
email: mathieu.bonnivard@u-paris.fr

I Pažanin
email: pazanin@math.hr

F J Suárez-Grau

A generalized Reynolds equation for micropolar flows past a ribbed surface with nonzero boundary conditions

Keywords: thin-film flow; micropolar fluid, rough boundary, homogenization, unfolding method

Inspired by the lubrication framework, in this paper we consider a micropolar fluid flow through a rough thin domain, whose thickness is considered as the small parameter ε while the roughness at the bottom is defined by a periodical function with period of order ε and amplitude ε δ , with δ > > 1. Assuming nonzero boundary conditions on the rough bottom and by means of a version of the unfolding method, we identify a critical case δ = 3 2 -1 2 and obtain three macroscopic models coupling the effects of the rough bottom and the nonzero boundary conditions. In every case we provide the corresponding micropolar Reynolds equation. We apply these results to carry out a numerical study of a model of squeeze-film bearing lubricated with a micropolar fluid. Our simulations reveal the impact of the roughness coupled with the nonzero boundary conditions on the performance of the bearing, and suggest that the introduction of a rough geometry may contribute to enhancing the mechanical properties of the device.

Introduction

Microfluidics is a multidisciplinary field intersecting engineering, physics, chemistry, microtechnology and biotechnology, with practical applications to the design of systems in which such small volumes of fluids will be used. Microfluidic area emerged in the beginning of the 1980s and is used in the development of inkjet printheads, DNA chips, lab-on-a-chip technology, micro-propulsion, and micro-thermal technologies.

Microfluidics deals with the manipulation of lubricants that are geometrically constrained to a small (typically sub-millimetre) scale, and with the experimental and theoretical study of their mechanical behaviour. This behaviour can differ from "macrofluidic" behaviour since, at the microscale, factors such as surface tension, energy dissipation, and fluidic resistance start to dominate the system. In particular, when a lubricant is in contact with a solid, at a small scale, surfacic effects may become preponderant. As a result, in order to reduce the energy dissipation of microscaled fluid-solid systems, one needs to understand and quantify very precisely the behaviour of the fluid near a solid wall.

From an experimental point of view, an efficient method to reduce the friction consists in using a certain type of rough boundaries, that are called riblets. These riblets are characterized by fast oscillations in the transversal direction, with a low amplitude, and by their constancy in the direction of the flow; they are essentially one-dimensional perturbations of the boundary of the solid. The aim of the use of riblets is to prevent vortices to appear in the neighborhood of the solid wall, and thus to reduce the momentum transfer from the vortices to the solid boundary. By using homogenization techniques, the influence of riblets on the slip behaviour of viscous fluids has been studied recently. In [START_REF] Bucur | On the asymptotic limit of flows past a ribbed boundary[END_REF], starting with perfect slip condition at a highly ribbed surface, it is showed that when the oscillating parameter goes to zero, no-slip condition appears in the transversal direction while perfect slip still holds in the direction of the flow. This means that riblets tend to prevent the fluid from slipping laterally, whereas the motion in the direction of the flow is allowed with no constraint. In the same spirit, it was proved in [START_REF] Casado-Díaz | Asymptotic behavior of the Navier-Stokes system in a thin domain with Navier condition on a slightly rough boundary[END_REF] that surfaces with low amplitude riblets give rise to a friction parameter in the transversal direction and no roughness effects in the direction of the flow.

The mathematical models for describing the motion of the lubricant in a device with small volume usually result from the simplification of the geometry of the lubricant film, i.e. its thickness. Using the film thickness as a small parameter ε, a simple asymptotic approximation can be easily derived providing a well-known Reynolds equation for the pressure of the fluid. Formal derivation goes back to the 19th century and the celebrated work of Reynolds [START_REF] Reynolds | On the theory of lubrication and its applications to Mr. Beauchamp Tower's experiments, including an experimental determination of the viscosity of olive oil[END_REF]. The justification of this approximation, namely the proof that it can be obtained as the limit of the Stokes system (as thickness tends to zero) is provided in [START_REF] Bayada | The transition between the Stokes equations and the Reynolds equation: a mathematical proof[END_REF] for a Newtonian flow between two plain surfaces. Different Reynolds equations for Newtonian fluids including roughness effects have been obtained for example in [START_REF] Bayada | New models in the theory of the hydrodynamic lubrication of rough surfaces[END_REF][START_REF] Benhaboucha | Asymptotic behaviour of pressure and stresses in a thin film flow with a rough boundary[END_REF][START_REF] Bonnivard | Homogenization of a Large Eddy Simulation Model for Turbulent Fluid Motion Near a Rough Wall[END_REF][START_REF] Bresch | Roughness-induced effect at main order on the Reynolds approximation[END_REF][START_REF] Casado-Díaz | A viscous fluid in a thin domain satisfying the slip condition on a slightly rough boundary[END_REF][START_REF] Casado-Díaz | Asymptotic behavior of the Navier-Stokes system in a thin domain with Navier condition on a slightly rough boundary[END_REF][START_REF] Chupin | Rigorous derivation of the thin film approximation with roughness-induced correctors[END_REF][START_REF] Pažanin | Effects of rough boundary on the heat transfer in a thin-film flow[END_REF][START_REF] Pažanin | Homogenization of the Darcy-Lapwood-Brinkman Flow in a Thin Domain with Highly Oscillating Boundaries[END_REF].

Nevertheless, most of the modern lubricants are no longer Newtonian fluids, since the use of additives in lubricants has become a common practice in order to improve their performance. Therefore, several microcontinuum theories [START_REF] Eringen | Simple Microfluids[END_REF] have been proposed to account for the effects of additives. Eringen micropolar fluid theory [START_REF] Eringen | Theory of micropolar fluids[END_REF] ignores the deformation of the microelements and allows for the particle micromotion to take place. From a mathematical point of view, a micropolar Reynolds equation was obtained in [START_REF] Bayada | On micropolar fluids in the theory of lubrication. Rigorous derivation of an analogue of the Reynolds equation[END_REF] for a micropolar flow in a thin film with a plain bottom assuming zero boundary conditions for microrotation. Other related results on the lubrication with a micropolar fluid with zero boundary condition can be found in [START_REF] Marušić-Paloka | An effective model for the lubrication with micropolar fluid[END_REF][START_REF] Prakash | Lubrication theory for micropolar fluids and its application to a journal bearing[END_REF], and some others references including roughness effects in [START_REF] Bonnivard | Effects of rough boundary and nonzero boundary conditions on the lubrication process with micropolar fluid[END_REF][START_REF] Boukrouche | Asymptotic analysis of a micropolar fluid flow in a thin domain with a free and rough boundary[END_REF][START_REF] Pažanin | Analysis of the thin film flow in a rough domain filled with micropolar fluid[END_REF][START_REF] Suárez-Grau | Analysis of the roughness regimes for micropolar fluids via homogenization[END_REF].

In the previously mentioned references, a zero boundary condition for the microrotation is assumed, implying that the fluid elements cannot rotate on the fluid-solid interface. If s is the horizontal velocity of the boundary, these conditions are written as follows:

u = s (u velocity), (1.1)
w = 0 (w microrotation).

(1.2)

However, more general boundary conditions for the microrotation were introduced to take into account the rotation of the microelements on the solid boundary. In the case where the boundary is flat, these conditions read

α 2 (∇ × u) × n = w × n , w • n = 0 , (1.3)
where n is a normal unit vector to the boundary. Conditions (1.3) were effectively proved to be in good accordance with experiments, see [START_REF] Bessonov | A new generalization of the Reynolds equation for a micropolar fluid and its application to bearing theory[END_REF][START_REF] Bessonov | Boundary viscosity conception in hydrodynamical theory of lubrication[END_REF][START_REF] Migun | On hydrodynamic boundary conditions for microstructural fluids[END_REF][START_REF] Prokhorenko | Experimental studies of polar indicator liquids used in capillary penetrant testing[END_REF]. The coefficient α describes the interaction between the given fluid and solid; it characterizes microrotation retardation on the solid surfaces.

In [START_REF] Bessonov | A new generalization of the Reynolds equation for a micropolar fluid and its application to bearing theory[END_REF], a generalized micropolar Reynolds equation is derived by using conditions (1.1), (1.3), and the relevance of the new parameter α regarding the performance of lubricated devices for both load and friction, is established by numerical computations. Nevertheless, it was mathematically proved in [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF] that it is not possible to consider the boundary condition (1.3) and simultaneously retain the no-slip condition (1.1) for the velocity. This would be like considering simultaneously, at the same boundary, a Neumann and a Dirichlet boundary condition. In order to obtain a well-posed variational formulation of the micropolar system, it is straightforward to confirm (see e.g. [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF]) that a velocity condition compatible with (1.3) needs to be introduced. This condition allows a slippage in the tangential direction and retains a non-penetration condition in the normal direction n (δ 0 is a real parameter)

(u -s) × n = δ 0 (∇ × w) × n , u • n = 0 .
(1.4)

It is worth stressing that in most lubrication studies, it is assumed that the speed of the lubricant at the surface equals that of the solid surface. However, it has been found that wall slip occurs, not only in non-Newtonian flows [START_REF] Bair | Shear strength measurements of lubricants at high pressure[END_REF][START_REF] Hatzikiriakos | Wall slip of molten high density polyethylene I -sliding plate rheometer studies[END_REF][START_REF] Jacobson | Non-Newtonian fluid model incorporated into elastohydrodynamic lubrication of rectangular contacts[END_REF][START_REF] Meng | A rheological model for lithinium lubricating grease[END_REF][START_REF] Spikes | The hanlf-wetted bearing. Part 1: extended Reynolds equation[END_REF][START_REF] Wilson | Viscoplastic behavior of silicone oil in a metalforming inlet zone[END_REF], but also in hydrodynamic lubrication or elastohydrodynamic lubrication [START_REF] Bonaccurso | Surface roughness and hydrodynamic boundary slip of a Newtonian fluid in a completely wetting system[END_REF][START_REF] Craig | Shear-dependent boundary slip in a aqueous Newtonian liquid[END_REF][START_REF] Hervet | Friction and slip of a simple liquid at a solid surface[END_REF][START_REF] Kaneta | Observation of wall slip in elastohydrodynamic lubrication[END_REF]. It seems that such phenomenon is linked to physical or chemical interactions of the solid surfaces with the lubricant. Several boundary conditions have been considered in those works to model the observed slippage. Most of them include limited yield stress or retain slippage value proportional to the shear stress. In that context, condition (1.4) appears as a new interpretation of the slippage observed in lubrication with micropolar fluids, expressed in terms of the microrotation field w.

In [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF], by using the nonzero boundary conditions (1.3)-(1.4) described above, in a 2-dimensional thin domain without roughness (see also [START_REF] Pažanin | Asymptotic analysis of the lubrication problem with nonstandard boundary conditions for microrotation[END_REF] for the 3D flow), Bayada et al. derive rigourously a generalized version of the Reynolds equation taking such boundary conditions into account. They perform their study in the critical case where one the non-Newtonian characteristic parameters of the micropolar fluid has specific (small) order of magnitude. The authors provide a comparison with the model in [START_REF] Bessonov | A new generalization of the Reynolds equation for a micropolar fluid and its application to bearing theory[END_REF] that uses the no-slip condition (1.1) for the velocity field, and observe that the introduction of slippage may enhance the performance of a bearing (that is, increase the load and reduce the friction coefficient) if the coupling number of the micropolar fluid and the nondimensional coefficient describing its slippage on the wall, are above a certain value.

Observe that in previous studies, the nonzero boundary condition has been considered on a plain bottom. In this paper, we impose this condition on a surface covered by riblets with low amplitude, and use asymptotic analysis to derive a micropolar Reynolds equation coupling the effects of the nonzero boundary conditions and the riblets. Since we are interested in the effect of the roughness, we adopt a simple geometric setting where the top boundary is plane, given by εh with h > 0. At the bottom we consider a surface covered by periodically distributed riblets with low amplitude, associated with a small parameter ε, where ε δ is the amplitude and ε is the period, where δ > > 1. This type of rough surface has been treated in [START_REF] Bonnivard | On the influence of wavy riblets on the slip behaviour of viscous fluids[END_REF][START_REF] Casado-Díaz | Asymptotic behavior of the Navier-Stokes system in a thin domain with Navier condition on a slightly rough boundary[END_REF][START_REF] Suárez-Grau | Effective boundary condition for a quasi-newtonian fluid at a slightly rough boundary starting from a Navier condition[END_REF][START_REF] Suárez-Grau | Asymptotic behavior of a non-Newtonian flow in a thin domain with Navier law on a rough boundary[END_REF] for fluid flows with Navier slip boundary conditions.

First, we identify a range of values of the coupling parameter N 2 , namely N 2 ≤ 1/2, under which there is existence and uniqueness of solution (Theorem 1). Later, by means of homogenization and reduction of dimension techniques, we identify the critical regime, i.e. δ = 3 2 -1 2 , in 3 which the nonzero boundary conditions make appear two friction parameters reflecting the riblets effect on both the effective velocity and micropolar fields (Theorem 2). Finally, we also obtain a precise description of the corresponding Reynolds equation which implicitly contains the effective nonzero boundary conditions describing the roughness effects (see (4.14) for more details). Moreover, we give the corresponding Reynolds equations corresponding to the sub-critical and super-critical regimes. This constitutes a generalization of the results of [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF] to domains with rough bottom (Theorem 3).

The paper is organized as follows. In Section 2 and 3, we formulate the problem and introduce some notation. In Section 4, we state our main results providing the homogenized model and the generalized Reynolds equation, which are proved in Section 5. The details of certain explicit computations or asymptotic developments are postponed to the Appendix. Finally, in Section 6 we conduct numerical simulations based on the generalized micropolar Reynolds equation obtained for a particular lubrication device: a squeze-film bearing.

Position of the problem

In the following, x ∈ R 3 is decomposed as x = (x , x 3 ) where x = (x 1 , x 2 ) ∈ R 2 and x 3 ∈ R.

We take e 1 , e 2 and e 3 to be the vectors of the canonical basis in R 3 , and e 1 , e 2 to be the vectors of the canonical basis in R 2 . The domain under consideration has the following form

Ω ε = (x , x 3 ) ∈ R 2 × R : x = (x 1 , x 2 ) ∈ L ω, -Ψ ε (x ) < x 3 < h c .
Here L is the characteristic length of the domain, ω ⊂ R 2 is an open subset with smooth boundary, c is the characteristic distance between the plates, h > 0 is an adimensional constant, ε is the ratio ε = c L and Ψ ε is defined by

1 L Ψ ε (x ) = λε δ Ψ 1 Lε x • e 1 (2.1)
see Figure 1, where λ > 0 is an amplitude parameter and δ, > 0 satisfy

1 < < δ. (2.2)
In definition (2.1), Ψ ∈ W 2,∞ # (R) is a R-valued function with period 1 (we use the index # to mean periodicity of period 1), that models the roughness profile on the lower surface, and that is normalized in the sense that

1 0 |∂ z 1 Ψ(z 1 )| 2 dz 1 = 1 . (2.3) Let Γ 0 ε , Γ 1
ε and Γ ε denote the lower, upper and lateral boundaries on Ω ε , namely

Γ 0 ε = (x , x 3 ) ∈ R 2 × R : x ∈ L ω, x 3 = -Ψ ε (x ) , Γ 1 ε = (x , x 3 ) ∈ R 2 × R : x ∈ L ω, x 3 = εhL , Γ ε = ∂Ω ε -(Γ 0 ε ∪ Γ 1 ε ). 4
The exterior normal n ε to Γ 0 ε is defined by

∀x ∈ Lω n ε (x , -Ψ ε (x )) = 1 [1 + ∂ x 1 Ψ ε (x ) 2 ] 1/2 -∂ x 1 Ψ ε (x ), 0, -1 .
(2.4)

For any vector field ξ defined on Γ 0 ε , we note [ξ] tan its tangential part, i.e. is the vector field defined on Γ

0 ε by [ξ] tan = ξ -(ξ • n ε )n ε . x 1 x 2 x 1
Figure 1: Representation of the geometry of riblets (upper) and of their cross section (lower). The riblets are periodic with period ε in the x1 direction, constant in the x2 direction and oscillate with an amplitude of order ε δ in the x3 direction.

The equations and boundary conditions

The micropolar fluid flow is described by the following equations expressing the balance of momentum, mass and angular momentum:

-(ν + ν r )∆u ε + ∇p ε = 2ν r (∇ × w ε ) , (2.5) div u ε = 0 , (2.6) -c r ∆w ε + 4ν r w ε = 2ν r (∇ × u ε ) . (2.7)
In the above system, velocity u ε , pressure p ε and microrotation w ε are unknown. ν is the Newtonian viscosity, while ν r and c r are microrotation viscosities resulting from the asymmetry of the stress tensor. All viscosity coefficients are assumed to be positive constants.

Let V ε be the velocity of the upper plate, and g ε be velocity of the fluid on the lateral boundaries of the domain. As discussed in the Introduction, the following boundary conditions are imposed

u ε = (0, 0, -V ε ), w ε = 0 on Γ 1 ε , (2.8
)

u ε = g ε , w ε = 0 on Γ ε , (2.9
)

u ε • n ε = 0, w ε • n ε = 0 on Γ 0 ε , (2.10) α 2 [Du ε n ε ] tan = w ε × n ε on Γ 0 ε , (2.11)
[Dw ε n ε ] tan = 2ν r c r β u ε × n ε on Γ 0 ε . (2.12)
We remark that along the paper, Du denotes the gradient of a vectorial function u = (u i ) 1≤i≤3 , defined by (Du) i,j = ∂ j u i , and should not be confused with the symmetric part of the gradient. Notice that the usual (Dirichlet) boundary conditions (2.8)-(2.9) for the velocity and microrotation are prescribed on Γ 1 ε ∪ Γ ε . However, on the lower part Γ 0 ε (corresponding to the rough boundary), new type of boundary conditions (2.11)-(2.12) are imposed, together with the non-penetration conditions (2.10). Finally, coefficient β ∈ R + in (2.12) is a friction coefficient that controls the slippage of the fluid at the wall.

Let us stress that conditions (2.11)-(2.12) are adaptations of conditions (1.3)-(1.4) to the present case of an oscillating boundary Γ 0 ε . Since system (2.5)-(2.7) couples a Stokes equation on u ε with an elliptic system on w ε , in the present context of slip boundary conditions, the normal Hence, conditions (2.11)-(2.12) and (1.3)-(1.4) are equivalent in the case of a flat boundary, so the tangential conditions (2.11)-(2.12) can be seen as a generalization of (1.3)- (1.4) to the case of a non flat boundary.

In [START_REF] Bessonov | Boundary viscosity conception in hydrodynamical theory of lubrication[END_REF], it was proposed to define the parameter α appearing in (2.11) as a microrotation retardation at the boundary and to connect it with the different viscosity coefficients. It has been shown experimentally [START_REF] Jacobson | At the boundary between lubrication and wear, Hutchings IM[END_REF][START_REF] Léger | Friction and flow with slip at fluid-solid interfaceds, Interfacial properties on the submicron scal[END_REF] that there are chemical interactions between a solid surface and the nearest fluid layer. This had to be taken into account, especially for a non-Newtonian fluid and a very thin-film thickness. This can be done by introducing a viscosity ν b near the surface which is different from ν and ν r . In [START_REF] Bessonov | Boundary viscosity conception in hydrodynamical theory of lubrication[END_REF], it was proposed to define α by means of this boundary viscosity

ν b by α = ν + ν r -ν b ν r . (2.13)
Following [START_REF] Bessonov | Boundary viscosity conception in hydrodynamical theory of lubrication[END_REF], it is possible to give physical limits to ν b , inducing limits on α:

0 ≤ ν b ≤ ν + ν r ⇒ 0 ≤ α ≤ ν + ν r ν r . ( 2

.14)

The condition α = 0 is equivalent to strong adhesion of the fluid particles to the boundary surface so that they do not rotate relative to the boundary, i.e. w = 0. Thus, from now on, we consider α > 0 so that the stress tensor and the micro-rotation are coupled on the boundary.

It has been observed (see e.g. [START_REF] Bayada | On micropolar fluids in the theory of lubrication. Rigorous derivation of an analogue of the Reynolds equation[END_REF][START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF][START_REF] Pažanin | On the micropolar flow in a circular pipe: the effects of the viscosity coefficients[END_REF]) that the magnitude of the viscosity coefficients appearing in (2.5)-(2.7) may influence the effective flow. Thus, it is reasonable to work with the system written in a non-dimensional form. In view of that, we introduce the characteristic velocity V 0 of the fluid, and define:

x = x L , x 3 = x 3 L , Ψ ε = Ψ ε L , u ε = u ε V 0 , p ε = L V 0 (ν + ν r ) p ε , w ε = L V 0 w ε , g ε = g ε V 0 , V ε = V ε V 0 , N 2 = ν r ν + ν r , R M = c r ν + ν r 1 L 2 .
(2.15) Dimensionless (non-Newtonian) parameter N 2 characterizes the coupling between the equations for the velocity and the microrotation, and is of order O(1) with respect to small parameter ε. Notice that assumption (2.14) yields

1 α ≥ N 2 .
(2.16)

The second dimensionless parameter denoted by R M is related to the characteristic length of the microrotation effects and will be compared with ε. We also assume that friction parameter β is of order O(1).

In view of the above changes of variables, the fluid domain becomes

Ω ε = (x , x 3 ) ∈ R 2 × R : x ∈ ω, -Ψ ε x < x 3 < ε h ,
where according to (2.1), Ψ ε is given by

Ψ ε (x ) = λε δ Ψ 1 ε x • e 1 (2.17)
and the lower, upper and lateral boundaries are now described by

Γ 0 ε = (x , x 3 ) ∈ R 2 × R : x ∈ ω, x 3 = -Ψ ε x , Γ 1 ε = (x , x 3 ) ∈ R 2 × R : x ∈ ω, x 3 = εh , Γ ε = ∂Ω ε -(Γ 0 ε ∪ Γ 1 ε ) .
The exterior normal n ε to Γ 0 ε is now defined by

∀x ∈ ω n ε x , -Ψ ε x = 1 1 + λ 2 ε 2(δ-) ∂ x 1 Ψ 1 ε x • e 1 2 1/2 -λε δ-∂ x 1 Ψ 1 ε x • e 1 , 0, -1 .
The tangential part of a vector field ξ defined on Γ 0 ε is accordingly given by [ξ] tan = ξ -(ξ • n ε )n ε . The flow equations (2.5)-(2.7) now have the following form

-∆u ε + ∇p ε = 2N 2 (∇ × w ε ) in Ω ε , (2.18) div u ε = 0 in Ω ε , (2.19
)

-R M ∆w ε + 4N 2 w ε = 2N 2 (∇ × u ε ) in Ω ε , (2.20)
with boundary conditions

u ε = -V ε e 3 , w ε = 0 on Γ 1 ε , (2.21)
u ε = g ε , w ε = 0 on Γ ε , (2.22)
u ε • n ε = 0, w ε • n ε = 0 on Γ 0 ε , (2.23) α 2 [Du ε n ε ] tan = w ε × n ε on Γ 0 ε , (2.24) R M [Dw ε n ε ] tan = 2N 2 β u ε × n ε on Γ 0 ε . (2.25)
The divergence-free condition (2.19) imposes the following compatibility condition on the boundary data:

Γ ε g ε • n ε dσ = V ε |ω| , (2.26)
where σ stands for the Hausdorff measure of dimension 2, and |ω| is the area of ω.

In the present paper the aim is to derive the macroscopic law describing the effective flow in Ω ε by using rigorous asymptotic analysis with respect to the small parameter ε. In particular, we shall focus on detecting the roughness-induced effects together with the effects of nonzero boundary conditions.

Let us start by defining the notion of weak solution to system (2.18)-(2.25).

Weak formulation of problem (2.18)-(2.25). Let us introduce the functional spaces V ε and V 0 ε defined by

V ε = {ϕ ∈ H 1 (Ω ε ) 3 , ϕ |Γ 1 ε ∪Γ ε = 0, ϕ • n ε = 0 on Γ 0 ε }, V 0 ε = {ϕ ∈ V ε , div ϕ = 0 in Ω ε },
Γ 0 ε [Du ε n ε ] tan • ϕ dσ - Ωε p ε div ϕ dx -2N 2 Ωε (∇ × w ε ) • ϕ dx = 0 .
Hence, boundary condition (2.24) yields

Ωε Du ε : Dϕ dx - Ωε p ε div ϕ dx -2N 2 Ωε (∇ × w ε ) • ϕ dx - 2 α Γ 0 ε (w ε × n ε ) • ϕ dσ = 0 .
Using the integration by part formula

Ωε (∇ × w ε ) • ϕ dx = Ωε (∇ × ϕ) • w ε dx - Γ 0 ε (w ε × n ε ) • ϕ dσ , (2.27)
the previous equality can be rewritten as

Ωε Du ε : Dϕ dx - Ωε p ε div ϕ dx -2N 2 Ωε w ε • (∇ × ϕ) dx -2 1 α -N 2 Γ 0 ε (w ε × n ε ) • ϕ dσ = 0 . (2.28)
Multiplying Equation (2.20) by another test function ψ ∈ V ε , integrating by parts and using boundary condition (2.25), we obtain

R M Ωε Dw ε : Dψ dx -2N 2 β Γ 0 ε (u ε × n ε ) • ψ dσ + 4N 2 Ωε w ε • ψ dx -2N 2 Ωε (∇ × u ε ) • ψ dx = 0 . ( 2
p ε div ϕ dx + R M Ωε Dw ε : Dψ dx -2N 2 Ωε (∇ × u ε ) • ψ dx -2N 2 Ωε w ε • (∇ × ϕ) dx + 4N 2 Ωε w ε • ψ dx -2 1 α -N 2 Γ 0 ε (w ε × n ε ) • ϕ dσ (2.30) -2N 2 β Γ 0 ε (u ε × n ε ) • ψ dσ = 0.
This leads to the following definition.

Definition 1. We say that

(u ε , w ε , p ε ) ∈ H 1 (Ω ε ) 3 × H 1 (Ω ε ) 3 × L 2 0 (Ω ε ) is a weak solution to system (2.18)-(2.25) if (u ε , w ε ) satisfy boundary conditions (2.21)-(2.23), div u ε = 0 in Ω and relation (2.30) holds for any (ϕ, ψ) ∈ V ε × V ε .

Notation

The unitary cube of R 2 will be denoted by Z = (-1 2 , 1 2 ) 2 , and we set Q = Z × (0, +∞). For any M > 0, we define Q M = Z × (0, M ). We introduce the space L 2 # (Z ), which is defined by the functions u in L 2 loc (R 2 ) and Z -periodic. The space L 2 # ( Q) is defined by the functions u in L 2 loc (R 2 × (0, +∞)) and

Q | u| 2 dz < +∞, u(z + k , z 3 ) = u(z), ∀k ∈ Z 2 , a.e. z ∈ R 2 × (0, +∞).
We define L 2 0 (O), with O a bounded and measurable subset of R N , by the functions of L 2 (O) with zero integral.

For every θ = (θ 1 , θ 2 ), we define

[θ ] ⊥ = (-θ 2 , θ 1 ), rot x 3 θ = ∂ x 3 [θ] ⊥ , Rot x θ = ∂ x 1 θ 2 -∂ x 2 θ 1 .
We define the sets

Ω - ε = (ω × (-∞, 0)) ∩ Ω ε , Ω + ε = (ω × (0, +∞)) ∩ Ω ε .
Given k ∈ Z 2 and τ > 0, we define

C k τ = τ Z + τ k , Q k τ = (C k τ × R) ∩ Θ ε ,
where Θ ε = {x ∈ R 2 × R : .Ψ ε (x ) < x 3 < ε}. We consider the function κ : R 2 → Z 2 given by

κ(x ) = k ⇔ x ∈ C k 1 .
We observe that κ is well defined, except for a set of zero measure in R 2 . In addition, for any τ > 0, it holds

κ x τ = k ⇔ x ∈ C k τ .
We denote C ε (x ), for a.e. x ∈ R 2 , by the square

C k ε such that x ∈ C k ε . Given ρ > 0, we take ω ρ = {x ∈ ω : dist(x, ∂ω) > ρ}, (3.1)
I ρ,ε = {k ∈ Z 2 : ω ρ ∩ C k ε = ∅}. By V we define the space of functions ϕ : R 2 × (0, +∞) → R such that ϕ ∈ H 1 # ( Q M )
, for every M > 0, and ∇ ϕ ∈ L 2 # ( Q) 3 . We remark that V is a Hilbert space by considering • V given by

ϕ 2 V = ∇ ϕ 2 L 2 ( Q) 3 + ϕ 2 L 2 (Z ×{0}) .
We observe that when we use O ε , we refer to a generic real sequence which is devoted to tend to zero when ε → 0. Moreover, O ε is allowed to change change from line to line. By C, we denote a generic positive constant, which does not depend on ε and it can also change from line to line.

Main results

As discussed before, different asymptotic behaviours of the flow may be deduced depending on the order of magnitude of the viscosity coefficients. Indeed, if we compare the characteristic number R M defined by (2.15) and appearing in the equation (2.20) with small parameter ε, three different asymptotic situations can be formally identified (see e.g. [START_REF] Bayada | About thin film micropolar asymptotic equations[END_REF][START_REF] Pažanin | On the micropolar flow in a circular pipe: the effects of the viscosity coefficients[END_REF][START_REF] Suárez-Grau | Analysis of the roughness regimes for micropolar fluids via homogenization[END_REF]). The most interesting one is, of course, the one leading to a strong coupling at main order, namely the regime

R M = ε 2 R c , R c = O(1). (4.1)
Hence, we will perform our analysis assuming the above scalings of R M and R c with respect to ε.

Concerning the other parameters, we recall that N 2 , α and β are of order O(1).

Besides, in the case of a squeeze film model, we also assume that the (vertical) velocity of the upper plate V ε is of order ε as ε tends to zero. Hence, we consider the asymptotic regime

V ε = εS, (4.2)
where S is a positive constant.

In order to study the asymptotic behaviour of the solution to system (2.18)-(2.25), we also need to assume a certain regularity on the boundary data g ε , and uniform estimates of relevant norms. A very general way of stating those properties is the following: there exists a sequence of lift functions J ε ∈ H 1 (Ω ε ) 3 satisfying div J ε = 0 in Ω ε , the boundary conditions

J ε = -V ε e 3 on Γ 1 ε , J ε = g ε on Γ ε , J ε • n ε = 0 on Γ 0 ε , (4.3)
and the estimates

∀ε > 0 J ε L 2 (Ωε) 3 ≤ Cε 1 2 , DJ ε L 2 (Ωε) 3×3 ≤ Cε -1 2 , J ε L 2 (Γ 0 ε ) 3 ≤ C, (4.4)
where C > 0 is a universal constant.

Remark 1. One typical construction of a boundary data g ε and the associate lift function J ε is the following, see [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF]. Consider a regular vector field J ∈ H 1 (Ω) 3 , satisfying div J = 0 in Ω, J = -Se 3 on ω × {h}, J = 0 on ω × {0}.

Extending J = (J , J 3 ) by zero on ω × (-∞, 0), we can define

J ε ∈ H 1 (Ω ε ) 3 by J ε (x , x 3 ) = (J (x , x 3 ε ), εJ 3 (x , x 3 ε )) ∀(x , x 3 ) ∈ Ω ε ,
and g ε := J ε |Γ ε in the sense of traces. By the change of variable (x , x 3 ) = (y , εy 3 ), there holds

Ωε |DJ ε | 2 dx dx 3 = ε Ω |D y J | 2 + 1 ε 2 |∂ y 3 J | 2 + ε 2 |∇ y J 3 | 2 + |∂ y 3 J 3 | 2 dy dy 3 , Ωε |J ε | 2 dx dx 3 = ε Ω |J | 2 + ε 2 |J 3 | 2 dy dy 3 ,
so that J ε satisfies all the required properties (4.3)-(4.4).

Since such vector field J is not unique, the lift function J ε and the boundary data g ε are quite arbitrary. In fact, they do not play a significant role in the asymptotic analysis of the problem, provided that conditions (4.3)-(4.4) are satisfied.

Let us start with an existence and uniqueness result for the solution of problems (2.18)-(2.25), whose proof is given in the Section 5.

Theorem 1. Assume that the coupling parameter N 2 satisfies the condition

N 2 ≤ 1 2 , (4.5)
and define the nonnegative parameter γ by

γ = 1 α -N 2 -N 2 β . (4.6)
Assume that the asymptotic regimes (4.1) and (4.2) hold. Then, for any β such that

γ 2 < R c (1 -2N 2 ) h 2 , ( 4.7)
there exists ε 0 > 0 such that for any 0 < ε < ε 0 , there exists a unique weak solution Remark 2. In the case of a flat boundary, Bayada et al. obtained in [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF] less restrictive conditions, namely N 2 < 1 and

(u ε , w ε , p ε ) in H 1 (Ω ε ) 3 × H 1 (Ω ε ) 3 × L 2 0 (Ω ε ) to system (2.
γ 2 < Rc(1-N 2 ) h 2
. However, we stress that the restriction of parameter N (4.5) that it is necessary to guarantee existence and uniqueness of the weak solution, is in fact in agreement with tribology models, where different considerations lead to the same assumption N 2 ≤ 1/2 (see [START_REF] Sinha | Micropolar squeeze films between rough rectangular plates[END_REF][START_REF] Sinha | The three-dimensional Reynolds equation for micropolar-fluid-lubricated bearings[END_REF]).

Rescaling

We also want to describe the asymptotic behaviour of the sequence (u ε , w ε , p ε ) of solution of the micropolar system (2.18)-(2.20) supplemented with boundary conditions (2.24)-(2.25), as ε tends to 0. We start by introducing a change of variables classically used in asymptotic analysis of flows in thin domains: the dilatation

y = x , y 3 = x 3 ε , (4.8)
which changes Ω ε to the set Ω ε of height of order h, defined as follows:

Ω ε = (y , y 3 ) ∈ R 2 × R : y ∈ ω, -Ψ ε (y ) < y 3 < h , (4.9)
where

Ψ ε (y ) = 1 ε Ψ ε (y ) = ε δ-1 Ψ 1 ε y • e 1 .
The lower, upper and lateral boundaries of the rescaled domain Ω ε are now defined by

Γ0 ε = (y , y 3 ) ∈ R 2 × R : y ∈ ω, y 3 = -Ψε (y ) , Γ1 ε = (y , y 3 ) ∈ R 2 × R : y ∈ ω, y 3 = h , Γ ε = ∂ Ωε -( Γ0 ε ∪ Γ 1 ε ).
Accordingly, we define the functions u ε ,

w ε ∈ H 1 ( Ω ε ) 3 and p ε ∈ L 2 0 ( Ω ε ) by u ε (y) = u ε (y , εy 3 ), w ε (y) = w ε (y , εy 3 ), p ε (y) = p ε (y , εy 3 ), a.e. y ∈ Ω ε . (4.10)
Since δ > 1, it is clear that the sequence of domains Ω ε converges (for instance, in the sense of Hausdorff complementary topology) to the limit domain Ω defined by

Ω = ω × (0, h) ⊂ R 2 × R .
We denote by Γ := ω × {0} the lower boundary of Ω.

The next step of the analysis is to identify the effective system satisfied by these rescaled functions.

Effective system

In this subsection, we give the result concerning the asymptotic behaviour of the rescaled functions u ε , w ε , p ε . Depending on the relation between the amplitude parameter δ and the period parameter , we obtain three different regimes, that we call critical, sub-critical and supercritical. Here, we state the result in the critical case; the other cases will be discussed in Remark 3. The proof of the corresponding results is given in Section 5.

Theorem 2. Assume that the asymptotic regimes (4.1) and (4.2) and conditions (4.5) and (4.7) hold. Assume that δ, satisfy the relation δ

= 3 2 -1 2 (critical case). Let (u ε , w ε , p ε ) be a sequence of weak solutions of (2.18)-(2.25). Then, there exist u , w ∈ H 1 (0, h; L 2 (ω)) 2 and p ∈ H 1 (ω) ∩ L 2 0 (ω) such that the rescaled functions u ε , w ε , p ε satisfy u ε ( u , 0) in H 1 (0, h; L 2 (ω)) 3 , ε w ε ( w , 0) in H 1 (0, h; L 2 (ω)) 3 , ε 2 p ε → p in L 2 (Ω) . (4.11)
The triplet ( u , w , p) is the unique solution of the following problem

           -∂ 2 y 3 u + ∇ y p -2N 2 rot y 3 w = 0 in Ω , -R c ∂ 2 y 3 w + 4N 2 w -2N 2 rot y 3 u = 0 in Ω , div y h 0 u (y , y 3 ) dy 3 = S in ω , (4.12)
with the boundary conditions u = 0, w = 0 on ω × {h}, (4.13)

∂ y 3 u = - 2 α [ w ] ⊥ + E λ (u • e 1 ) e 1 on Γ, R c ∂ y 3 w = -2N 2 β[ u ] ⊥ + R c F λ (w • e 1 ) e 1 on Γ . (4.14)
Coefficients E λ , F λ ∈ R appearing in boundary conditions (4.14) are defined by

E λ = Q |D z φ 1,λ | 2 dz, F λ = Q |D z φ 2,λ | 2 dz, (4.15)
where

( φ i,λ , q i,λ ) ∈ V 3 × L 2 ( Q), i = 1, 2, are respectively the solutions of                -∆ z φ 1,λ + ∇ z q 1,λ = 0 in R 2 × R + , div z φ 1,λ = 0 in R 2 × R + , φ 1,λ 3 (z , 0) = λ∂ z 1 Ψ(z • e 1 ) on R 2 × {0}, -∂ z 3 φ 1,λ 1 = 0, -∂ z 3 φ 1,λ 2 = 0 on R 2 × {0} , (4.16
)

and        -∆ z φ 2,λ = 0 in R 2 × R + , φ 2,λ 3 (z , 0) = λ∂ z 1 Ψ(z • e 1 ) on R 2 × {0}, -∂ z 3 φ 2,λ 1 = 0, -∂ z 3 φ 2,λ 2 = 0 on R 2 × {0} .
(4.17)

Remark 3. Theorem 2 can be adapted easily to describe the two other asymptotic regimes:

• In the sub-critical case δ > 3 2 -1 2 , the riblets are so small that there is no effect of roughness, so we obtain the nonzero boundary conditions on Γ,

∂ y 3 u = - 2 α [ w ] ⊥ on Γ, R c ∂ y 3 w = -2N 2 β[ u ] ⊥ on Γ .
Thus, we deduce that the model obtained in [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF][START_REF] Bayada | Wall slip induced by a micropolar fluid[END_REF] even holds for a very slightly rough boundary. • In the super-critical case 1 < δ < 3 2 -1 2 , the effect of the riblets is maximal. Thus, boundary conditions given in (4.14) are replaced by

u • e 1 = w • e 1 = 0 on Γ, ∂ y 3 u • e 2 = ∂ y 3 w • e 2 = 0 on Γ . (4.18)
Thus, we deduce that the roughness is so strong that the fluid adheres to the boundary and fluid elements cannot rotate on the fluid-solid interface in the x 1 -direction.

Generalized micropolar Reynolds equations

In this subsection, we obtain a generalized Reynolds equation associated to the homogenized micropolar system given in Theorem 2 (critical case). For the sake of simplicity, we will consider a 2D domain characteristic of the lubrication assumption. Thus, we consider in Theorem 2 that the flow does not depend on the y 2 -coordinate, and that velocity component u 2 and micropolar component w 1 are zero. Hence, we address the following limit problem posed in Ω = (0, 1) × (0, h):

-∂ 2 y 3 u 1 + ∂ y 1 p + 2N 2 ∂ y 3 w 2 = 0 in Ω, -R c ∂ 2 y 3 w 2 + 4N 2 w 2 -2N 2 ∂ y 3 u 1 = 0 in Ω, (4.19)
completed with the boundary conditions

u 1 = 0, w 2 = 0 on Γ 1 = (0, 1) × {h}, (4.20)
∂ y 3 u 1 = 2 α w 2 + E λ u 1 on Γ = (0, 1) × {0}, R c ∂ y 3 w 2 = -2N 2 β u 1 on Γ , (4.21)
and the incompressibility condition

∂ y 1 h 0 u 1 (y 1 , y 3 ) dy 3 = S in (0, 1) . (4.22)
We give in the Appendix the expression of (ũ 1 , w2 ), solution of system (4.19)-(4.21), in terms of p (see Lemmas 4 and 5). Putting these expressions in (4.22) will lead to the corresponding Reynolds equations that take into account the roughness-induced effects.

Theorem 3. In the critical case δ = 3 2 -1 2 , the pressure p satisfies the following Reynolds equation

1 0 Θ λ ∂ y 1 p(y 1 ) ∂ y 1 θ(y 1 ) dy 1 = 1 0 Sθ(y 1 ) dy 1 , ∀θ ∈ H 1 (0, 1), (4.23)
with Θ λ defined in the case α = 1 by

Θ λ = h 3 3(1-N 2 ) -(1 -η λ ) 3h 3 4(1-N 2 ) -2N 2 k ch(kh)-1 k -η λ hsh(kh) + γα 2 h 2 (1 -2η λ ) -(1 -η λ ) γ α h + 2N 2 k sh(kh) A -2N 2 k sh(kh) k -η λ hch(kh) -(1 -η λ )(1 + ch(kh))h N 2 k B , (4.24
) and in the case α = 1 by

Θ λ = -1 2(1-N 2 ) h 3 3 -µ λ h 3 -(1 -µ λ ) h 2 k(1-N 2 ) 1-ch(kh) sh(kh) -1 1-N 2 h 2 2 -µ λ h 2 + (1 -µ λ ) h k(1-N 2 ) 1-ch(kh) sh(kh) A -2N 2 k sh(kh) k -µ λ hch(kh) -(1 -µ λ ) 2N 2 k h + (1-ch(kh)) 2 k sh(kh) B , (4.25)
where A, A , B and B are defined in Lemmas 4 and 5 in the Appendix.

Remark 4. It is worth mentioning that the effective expressions given in Lemmas 4 and 5, and Theorem 3 are explicitly corrected by the roughness-induced coefficient E λ . Indeed, by putting E λ = 0, which implies η λ = µ λ = 1 (i.e., no roughness introduced), we obtain the same expressions as derived in [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF], which also corresponds to the sub-critical case δ > 3 2 -1 2 .

Using the explicit expressions from Lemma 4, Lemma 5 and Theorem 3, it is possible to develop ũ1 , w2 and p in powers of λ 2 . This will be useful in the numerical computations from Section 6. However, since the corresponding formulas are rather long, we have gathered them in the Appendix for the sake of clarity (see Corollary 2).

Finally, we give the micropolar Reynolds equation corresponding to the super-critical case 1 < δ < 3 2 -1 2 . As in the critical case, its derivation is based on explicit expressions of the velocity and microrotation (see Lemma 6 in the Appendix).

Theorem 4. In the super-critical case 1 < δ < 3 2 -1 2 , the pressure p satisfies the following Reynolds equation

1 0 Θ∂ y 1 p(y 1 ) ∂ y 1 θ(y 1 ) dy 1 = 1 0 Sθ(y 1 ) dy 1 , ∀θ ∈ H 1 ((0, 1)), (4.26)
with Θ defined by

Θ = h 3 12(1 -N 2 ) - 2N 2 k ch(kh) -1 k - h 2 sh(kh) A - 2N 2 k sh(kh) k - h 2 (ch(kh) + 1) B ,
where A and B are defined in Lemma 6.

Proofs of the main results

We start by proving the existence and uniqueness of solution of problem (2.18)-(2.25).

Proof of Theorem 1.

Let J ε ∈ H 1 (Ω ε ) 3 be a sequence of free divergence lift functions satisfying (4.3)-(4.4). Replacing u ε by v ε + J ε in the weak formulation (2.30), we see that (u ε , w ε , p ε ) is a weak solution to system (2.18)-(2.25) if and only if (v ε , w ε , p ε ) ∈ V 0 ε × V ε × L 2 0 (Ω ε ) and satisfies for any (ϕ, ψ) ∈ V ε × V ε Ωε Dv ε : Dϕ dx - Ωε p ε div ϕ dx + R M Ωε Dw ε : Dψ dx -2N 2 Ωε (∇ × v ε ) • ψ dx -2N 2 Ωε w ε • (∇ × ϕ) dx + 4N 2 Ωε w ε • ψ dx -2 1 α -N 2 Γ 0 ε (w ε × n ε ) • ϕ dσ -2N 2 β Γ 0 ε (v ε × n ε ) • ψ dσ = - Ωε DJ ε : Dϕ dx + 2N 2 Ωε (∇ × J ε ) • ψ dx + 2N 2 β Γ 0 ε (J ε × n ε ) • ψ dσ . (5.1)

Equation (5.1) justifies the introduction of the bilinear forms

A ε : (V ε × V ε ) 2 → R and B ε : (V ε × V ε ) × L 2 0 (Ω ε ) → R respectively defined by A ε ((v, w), (ϕ, ψ)) = Ωε Dv : Dϕ dx - Ωε p ε div ϕ dx + R M Ωε Dw : Dψ dx -2N 2 Ωε (∇ × v) • ψ dx -2N 2 Ωε w • (∇ × ϕ) dx + 4N 2 Ωε w • ψ dx -2 1 α -N 2 Γ 0 ε (w × n ε ) • ϕ dσ -2N 2 β Γ 0 ε (v × n ε ) • ψ dσ , (5.2)
and

B ε ((v, w), q) = - Ωε q div v dx , (5.3)
and of the linear form

L ε : V ε × V ε → R defined by L ε (ϕ, ψ) = - Ωε DJ ε : Dϕ dx + 2N 2 Ωε (∇ × J ε ) • ψ dx + 2N 2 β Γ 0 ε (J ε × n ε ) • ψ dσ . (5.4) Hence, (u ε , w ε , p ε ) is a weak solution to system (2.18)-(2.25) if and only if (v ε , w ε , p ε ) ∈ V ε × V ε × L 2 0 (Ω ε
) and satisfies the following mixed formulation

A ε ((v, w), (ϕ, ψ)) + B ε ((ϕ, ψ), p ε ) = L ε (ϕ, ψ) ∀(ϕ, ψ) ∈ V ε × V ε , (5.5)
B ε ((v ε , w ε ), q) = 0 ∀q ∈ L 2 (Ω ε ) . (5.6)
The existence and uniqueness of the solution (v ε , w ε , p ε ) to the mixed formulation (5.5)-(5.6) is established in [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF], in the case where the oscillating boundary Γ 0 ε is replaced by a flat boundary ω × {0}. For the sake of completeness, we recall the main steps of the proof, highlighting the differences that are implied by the oscillations of the lower boundary Γ 0 ε . First, let us state some useful quantitative inequalities.

Trace inequality on Γ 0 ε . Since the lower boundary Γ 0 ε is not flat, one needs to take into account the variations of the normal direction n ε in order to estimate the L 2 -norm of the trace of a function ψ ∈ V ε . To this end, we introduce the quantity τ ε defined by

τ ε := sup x ∈ω 1 + |∇ x Ψ ε (x )| 2 .
We also denote by h ε the height of the domain Ω ε , defined by

h ε := sup x ∈ω εh + Ψ ε (x ) = ε sup x ∈ω h + λε δ-1 Ψ( 1 ε x • e 1 ) .
In particular, lim ε→0 h ε /ε = h. With this notation, there holds the trace inequality

∀ψ ∈ V ε ψ L 2 (Γ 0 ε ) 3 ≤ τ ε h ε Dψ L 2 (Ωε) 3×3 . (5.7)
In fact, this inequality holds true for any vector field

ψ ∈ H 1 (Ω ε ) 3 vanishing on Γ 1 ε . By density, it is enough to prove it for any ψ ∈ H 1 (Ω ε ) 3 ∩ C 1 (Ω ε ) 3 such that ψ |Γ 1 ε = 0.
Integrating on vertical lines, we obtain

Γ 0 ε |ψ| 2 = ω |ψ(x , -Ψ ε (x ))| 2 1 + |∇ x Ψ ε (x )| 2 dx ≤ τ ε ω |ψ(x , -Ψ ε (x ))| 2 dx ≤ τ ε ω εh -Ψε(x ) ∂ x 3 ψ(x , x 3 ) dx 3 2 dx ≤ τ ε ω (εh + sup ω Ψ ε ) εh -Ψε(x ) |∂ x 3 ψ(x , x 3 )| 2 dx 3 dx ≤ τ ε h ε Ωε |Dψ| 2 dx .
This proves (5.7).

Poincaré inequality. In the same fashion, Poincaré inequality in

V ε reads ∀ϕ ∈ V ε ϕ L 2 (Ωε) 3 ≤ h ε Dϕ L 2 (Ωε) 3×3 . (5.8)
Relation between ∇ × ϕ L 2 (Ωε) 3 and Dϕ L 2 (Ωε) 3×3 . Let us recall that for any vector field ϕ ∈ V ε ,

Ωε |div ϕ| 2 + |∇ × ϕ| 2 dx = Ωε |Dϕ| 2 dx + Γ 0 ε ((ϕ • ∇)n ε ) • ϕ dσ , (5.9)
(see, for instance, [START_REF] Boyer | Mathematical tools for the study of the incompressible Navier-Stokes equations and related models[END_REF], formula (IV.23)). In particular, if the lower boundary Γ 0 ε was flat, the identity ∇ × ϕ 2 L 2 (Ωε) 3 = Dϕ 2 L 2 (Ωε) 3×3 would hold for any ϕ ∈ V 0 ε , since the remaining term

Γ 0 ε ((ϕ • ∇)n ε )
• ϕ dσ would vanish. However, in the present geometric configuration, one cannot expect this term to be zero in general. In fact, a classical estimate reads

Γ 0 ε ((ϕ • ∇)n ε ) • ϕ dσ ≤ Lip(n ε ) ϕ 2 L 2 (Γ 0 ε ) 3 ,
where Lip(n ε ) is the Lipschitz constant of the normal vector field n ε , locally extended in a neighborhood of the surface {x 3 = Ψ ε (x )} (in the sense of [START_REF] Boyer | Mathematical tools for the study of the incompressible Navier-Stokes equations and related models[END_REF]Section 3.4]). However, using definition of Ψ ε given in (2.17) and condition (2.2) on parameters δ, , it turns out that in the general case where ∂ 2 11 Ψ ∞ > 0, Lip(n ε ) is of order ε δ-2 , hence diverging since ε δ-goes to zero. As a result, we cannot use identity (5.9) to estimate the L 2 norms of div ϕ and ∇ × ϕ over Ω ε by the L 2 norm of Dϕ, as is done in [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF] in the case of a flat boundary.

Instead, we rely on the following elementary estimate:

∇ × ϕ L 2 (Ωε) 3 ≤ √ 2 Dϕ L 2 (Ωε) 3×3 ∀ϕ ∈ H 1 (Ω ε ) 3 .
(5.10)

The presence of the constant √ 2 in the previous estimate is at the origin of the term 2N 2 in condition (4.7), which was simply N 2 in the case of a flat surface, as established in [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF].

Existence and uniqueness of the solution of the mixed formulation (5.5)-(5.6). Using Cauchy-Schwarz inequality and inequalities (5.7), (5.8) and (5.10), it is easy to see that A ε , B ε and L ε are continuous on their respective domains of definition, for any fixed value of parameter ε. Hence, noticing that by definition of V 0 ε ,

V 0 ε × V ε = (ϕ, ψ) ∈ V ε × V ε , B ε ((ϕ, ψ), q) = 0 for any q ∈ L 2 0 (Ω ε ) ,
and denoting by (•, •) Vε×Vε the norm defined by

(ϕ, ψ) Vε×Vε = Dϕ 2 L 2 (Ωε) 3×3 + Dψ 2 L 2 (Ωε) 3×3 1/2
the existence and uniqueness of the solution (v ε , w ε , p ε ) to the mixed formulation (5.5)-(5.6) result from the following properties (see [34, paragraph 4.1 p. 57]) : (i) coerciveness of A ε : there exists η = η(ε) > 0 such that

∀(ϕ, ψ) ∈ V 0 ε × V ε A ε ((ϕ, ψ), (ϕ, ψ)) ≥ η Dϕ 2 L 2 (Ωε) 3×3 + Dψ 2 L 2 (Ωε) 3×3 , (ii) inf-sup condition: there exists c = c(ε) > 0 such that inf q∈L 2 0 (Ωε) sup (ϕ,ψ)∈Vε×Vε B ε ((ϕ, ψ), q) (ϕ, ψ) Vε×Vε q L 2 0 (Ωε) ≥ c .
The inf-sup condition (ii) can be proved using the exact same arguments as in the proof of Theorem 2.2 in [START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF], that relies on the solvability in H 1 0 (Ω ε ) 3 of equation div ϕ = q, for an arbitrary q ∈ L 2 0 (Ω ε ), with natural estimates.

To establish the coerciveness condition (i), we use Hölder inequality, Poincaré inequality (5.8), the trace inequality (5.7) and estimate (5.10) to obtain the lower estimate

A ε ((ϕ, ψ), (ϕ, ψ)) = Ωε |Dϕ| 2 dx + R M Ωε |Dψ| 2 dx -4N 2 Ωε (∇ × ϕ) • ψ dx + 4N 2 Ωε |ψ| 2 dx -2 1 α -N 2 Γ 0 ε (ψ × n ε ) • ϕ dσ -2N 2 β Γ 0 ε (ϕ × n ε ) • ψ dσ ≥ Dϕ 2 L 2 (Ωε) 3×3 + R M Dψ 2 L 2 (Ωε) 3×3 -4N 2 √ 2 Dϕ L 2 (Ωε) 3×3 ψ L 2 (Ωε) 3 (5.11) + 4N 2 ψ 2 L 2 (Ωε) 3 -2γ τ ε h ε Dϕ L 2 (Ωε) 3×3 Dψ L 2 (Ωε) 3×3 ,
where γ is the defined by (4.6). Now, by condition (4.7), there exists

c 1 > 0 satisfying γh R c < c 1 < 1 -2N 2 γh , (5.12)
and by Young inequality,

Dϕ L 2 (Ωε) 3×3 Dψ L 2 (Ωε) 3×3 ≤ c 1 2ε Dϕ 2 L 2 (Ωε) 3×3 + ε 2c 1 Dψ 2 L 2 (Ωε) 3×3 .
By continuity, there exists a real number c 2 satisfying 0 < c 2 < 1, and such that

c 1 < 1 -2N 2 c 2 γh , (5.13)
and we also have

Dϕ L 2 (Ωε) 3×3 ψ L 2 (Ωε) 3 ≤ √ 2 4c 2 Dϕ 2 L 2 (Ωε) 3×3 + c 2 √ 2 ψ 2 L 2 (Ωε) 3 .
Going back to estimate (5.11), we obtain

A ε ((ϕ, ψ), (ϕ, ψ)) ≥ A ε Dϕ 2 L 2 (Ωε) 3×3 + ε 2 B ε Dψ 2 L 2 (Ωε) 3×3 + 4N 2 (1 -c 2 ) ψ 2 L 2 (Ωε) 3 , (5.14)
where A ε , B ε are defined by

A ε = 1 - 2N 2 c 2 -c 1 γτ ε h ε ε , B ε = R c - γτ ε c 1 h ε ε .
In particular, there holds

lim ε→0 A ε = 1 - 2N 2 c 2 -c 1 γh , lim ε→0 B ε = R c - γ c 1 h .
Using conditions (5.12)-(5.13), we conclude that for ε small enough, A ε is coercive.

A priori estimates and convergences

In this subsection we give a priori estimates and convergence results for the rescaled functions u ε , w ε , p ε . Also, in order to take into account the effects of the rough boundary, we will introduce the unfolding method before proceeding with the proof of the Theorem 2.

Proposition 1. Assume that the asymptotic regimes (4.1) and (4.2) and conditions (4.5) and (4.7) hold. Then, there exists a constant C > 0, independent on ε, such that for any 0 < ε < ε 0 , the solution (u ε , w ε , p ε ) satisfies the estimates

u ε L 2 (Ωε) 3 ≤ Cε 1 2 , Du ε L 2 (Ωε) 3×3 ≤ Cε -1 2 ,
(5.15)

w ε L 2 (Ωε) 3 ≤ Cε -1 2 , Dw ε L 2 (Ωε) 3×3 ≤ Cε -3 2 , (5.16
)

p ε L 2 (Ωε) ≤ Cε -3 2 , ∇p ε H -1 (Ωε) 3 ≤ Cε -1 2 .
(5.17)

Moreover, the rescaled solution ( u ε , w ε , p ε ) defined by (4.10) satisfies the estimates:

u ε L 2 ( Ωε) 3 ≤ C, ∂ y 3 u ε L 2 ( Ωε) 3 ≤ C, D y u ε L 2 ( Ωε) 3×2 ≤ Cε -1 , (5.18)
w ε L 2 ( Ωε) 3 ≤ Cε -1 , ∂ y 3 w ε L 2 ( Ωε) 3 ≤ Cε -1 , D y w ε L 2 ( Ωε) 3×2 ≤ Cε -2 , (5.19)
p ε L 2 ( Ωε) ≤ Cε -2 , ∂ y 3 p ε H -1 ( Ωε) ≤ Cε -1 , ∇ y p ε H -1 ( Ωε) 2 ≤ Cε -2 .
(5.20)

Proof of Proposition 1. According to Theorem 1, there exists ε 0 > 0 such that for any 0 < ε < ε 0 , there exists a unique weak solution (u ε , [START_REF] Casado-Díaz | Asymptotic behavior of the Navier-Stokes system in a thin domain with Navier condition on a slightly rough boundary[END_REF]). Now, we obtain the estimates on the velocity and microrotation and then, we obtain the estimates for the pressure.

w ε , p ε ) ∈ H 1 (Ω ε ) 3 × H 1 (Ω ε ) 3 × L 2 0 (Ω ε ) to system (2.18)-(2.
Estimates on u ε , w ε , u ε and w ε . Using Hölder inequality, inequalities (5.7) and (5.8), and estimate (5.10), we obtain for any (ϕ, ψ) ∈ V ε × V 0 ε :

L ε (ϕ, ψ) ≤ DJ ε L 2 (Ωε) 3×3 Dϕ L 2 (Ωε) 3×3 + 2N 2 ∇ × J ε L 2 (Ωε) 3 h ε Dψ L 2 (Ωε) 3×3 + 2N 2 β J ε L 2 (Γ 0 ε ) 3 √ τ ε h ε Dψ L 2 (Ωε) 3×3 (5.21) ≤ DJ ε L 2 (Ωε) 3×3 Dϕ L 2 (Ωε) 3×3 + 2N 2 h ε C ε Dψ L 2 (Ωε) 3×3 ,
where C ε is defined by

C ε = √ 2 √ h ε DJ ε L 2 (Ωε) 3×3 + β √ τ ε J ε L 2 (Γ 0 ε ) 3 .
Notice that, by properties (4.4) satisfied by J ε , there exists a constant C > 0 such that C ε ≤ C for any ε > 0.

To obtain estimates (5.15) and (5.16), we test against (ϕ, ψ) = (v ε , w ε ) in (5.5) and use inequalities (5.21) and (5.14) to obtain

A ε Dv ε 2 L 2 (Ωε) 3×3 + ε 2 B ε Dw ε 2 L 2 (Ωε) 3×3 ≤ DJ ε L 2 (Ωε) 3×3 Dv ε L 2 (Ωε) 3×3 + 2N 2 h ε C ε Dw ε L 2 (Ωε) 3×3 .

Last inequality can be written equivalently as

A ε Dv ε L 2 (Ωε) 3×3 - 1 2A ε DJ ε L 2 (Ωε) 3×3 2 + ε 2 B ε Dw ε L 2 (Ωε) 3×3 - N 2 ε 2 B ε h ε C ε 2 ≤ 1 4A ε DJ ε 2 L 2 (Ωε) 3×3 + N 4 ε 2 B ε h ε C 2 ε .
Since τ ε , b ε /ε, A ε , B ε and C ε are uniformly bounded, using assumptions (4.4), we deduce than the right hand side of the previous inequality is bounded by C/ε for a certain constant C > 0. This implies the following bounds:

Dv ε L 2 (Ωε) 3×3 ≤ Cε -1/2 , Dw ε L 2 (Ωε) 3×3 ≤ Cε -3/2 .
Hence, using Poincaré inequality (5.8), the relation u ε = v ε + J ε and properties (4.4) satisfied by J ε , we obtain the desired estimates (5.15)-(5.16). Finally, estimates (5.18)-(5.19) are direct consequences of the rescaling (4.10).

Estimates on p ε and p ε . In order to estimate ∇p ε in

H -1 (Ω ε ) 3 , we test against ϕ ∈ H 1 0 (Ω ε ) 3 in (2.28) : ∇p ε , ϕ H -1 (Ωε) 3 ×H 1 0 (Ωε) 3 = - Ωε p ε div ϕ dx = - Ωε Du ε : Dϕ dx + 2N 2 Ωε w ε • (∇ × ϕ) dx .
Using Hölder inequality and estimate (5.10), we deduce

∇p ε , ϕ H -1 (Ωε) 3 ×H 1 0 (Ωε) 3 ≤ Du ε L 2 (Ωε) 3×3 + 2N 2 √ 2 w ε L 2 (Ωε) 3×3 Dϕ L 2 (Ωε) 3×3 .
We conclude from the upper bounds (5.15)-(5.16

) that ∇p ε H -1 (Ωε) 3 ≤ Cε -1/2
, where C depends only on N . Finally, to estimate p ε in L 2 0 (Ω ε ), we apply the following inequality, whose proof is given in [25, Corollary 4.2]:

p ε L 2 0 (Ωε) ≤ C ε ∇p ε H -1 (Ωε) 3 .
This proves (5.17). Finally, estimates (5.20) are direct consequences of the rescaling (4.10), which concludes the proof of Proposition 1.

As a consequence of the a priori estimates stated in Proposition 1, and the fact that Ω ⊂ Ω ε and | Ω ε \ Ω| → 0, we have the following convergences for rescaled solutions u ε , w ε restricted to the limit domain Ω.

Lemma 1. Assume that the asymptotic regimes (4.1) and (4.2) and conditions (4.5) and (4.7) hold. Then, for a subsequence of ε, still denoted by ε, there exist u , w ∈ H 1 (0, h; L 2 (ω)) 2 with u (x , h) = w (x , h) = 0 for a.e. x ∈ ω, and

div y h 0 u (y , y 3 ) dy 3 = S in H -1 (ω), (5.22)
such that

u ε|Ω ( u , 0) in H 1 (0, h; L 2 (ω)) 3 , (5.23) ε w ε|Ω ( w , 0) in H 1 (0, h; L 2 (ω)) 3 .
(5.24)

Proof. The space H 1 (0, h; L 2 (ω)) is a Hilbert space for the norm

v H 1 (0,h;L 2 (ω)) = v 2 L 2 (Ω) + ∂ y 3 v 2 L 2 (Ω) 1/2 .
By estimates (5.18)-(5.19), u ε|Ω and ε w ε|Ω are bounded in H 1 (0, h; L 2 (ω)) 3 , so there exist u and w such that for a subsequence of ε, still denoted by ε, we have u ε|Ω u and ε w ε|Ω w in H 1 (0, h; L 2 (ω)) 3 .

(5.25)

By continuity of the trace operator from H 1 (0, h; L 2 (ω)) into L 2 (ω×{h}), the conditions u ε (x , h) = -εSe 3 and ε w ε (x , h) = 0 for a.e. x ∈ ω pass to the limit, yielding u(x , h) = w(x , h) = 0 for a.e.

x ∈ ω. Now, we prove that u 3 = 0. Since u ε is divergence free, using definition (4.10), the rescaled function u ε satisfies

div y u ε + 1 ε ∂ y 3 u ε,3 = 0 a.e. in Ω ε , (5.26
)

so for any φ ∈ C ∞ c (Ω), 0 = Ω div y u ε φ dy + Ω 1 ε ∂ y 3 u ε,3 φ dy = - Ω u ε • ∇ y φ dy + Ω 1 ε ∂ y 3 u ε,3 φ dy .
Hence, Ω ∂ y 3 u ε,3 φ dx = ε Ω u ε • ∇ y φ dx and using (5.25), we deduce that Ω ∂ y 3 u 3 φ dx = 0. As a result, u 3 does not depend on y 3 , and since it vanishes on y 3 = h, it is identically null. Next, we prove the divergence equation (5.22). Using condition (2.19), integration by parts, boundary conditions (2.21), (2.23) and the change of variables (4.8), we have for any

φ ∈ C ∞ c (ω) 0 = Ωε (div u ε )φ(x ) dx = - Ωε u ε • ∇ x φ dx + Γ 1 ε u ε,3 φ dσ = -ε Ωε u ε • ∇ y φ(y ) dy -εS ω φ(y ) dy .
Noticing that, by the bound (5.18) and Hölder inequality, there holds lim ε→0 Ωε\Ω | u ε | 2 dy = 0, we deduce that -

Ω u ε • ∇ y φ(y ) dy = S ω φ(y ) dy + O ε .
Using the weak convergence (5.25), we can pass to the limit in the previous equality and obtain

S ω φ(y ) dy = - Ω u • ∇ y φ(y ) dy = ω div y h 0 u (y , y 3 ) dy 3 φ(y ) dy ,
which proves (5.22). Finally, it remains to prove that w 3 = 0. To do this, for any ψ ∈ C ∞ c (Ω), we consider ψ ε = εψ(x , x 3 /ε)e 3 as test function in the variational formulation (2.29). Applying the change of variables (4.8) and extending the integrals to Ω, we get

R c Ω ε∂ y 3 w ε,3 ∂ y 3 ψ 3 dy + 4N 2 Ω ε w ε,3 ψ dy = 2N 2 Ω εRot x u ε ψ 3 dy + O ε .
Integrating by parts the right-hand side, we get

R c Ω ε∂ y 3 w ε,3 ∂ y 3 ψ 3 dy + 4N 2 Ω ε w ε,3 ψ dy = 2N 2 Ω ε[ u ε ] ⊥ ∇ y ψ 3 dy + O ε .
Using convergences (5.25), when ε tends to zero, we get A similar argument works for the first term of (5.28) works replacing ϕ(y ) by λε δ-∂ 1 Ψ( 1 ε y • e 1 )ϕ(y ), which goes to 0 in L ∞ (ω).

Then, from the above, passing to the limit in (5.28), we get ω w 3 (y , 0)ϕ(y )dy = 0, which is equivalent to w 3 (y , 0) = 0 for a.e. y ∈ ω.

Finally, from (5.27) and taking into account that w 3 (y , h) = w 3 (y , 0) = 0 for a.e. y ∈ ω, it is easily deduced that w 3 = 0, which ends the proof.

In order to give the convergence of the rescaled pressure p ε , let us give a more accurate estimate for pressure p ε . For this, we need to recall a decomposition result for p ε whose proof can be found in [START_REF] Casado-Díaz | Asymptotic behavior of the Navier-Stokes system in a thin domain with Navier condition on a slightly rough boundary[END_REF]Corollary 4.2].

Proposition 2. The following decomposition for

p ε ∈ L 2 0 (Ω ε ) holds p ε = p 0 ε + p 1 ε , (5.29)
where p 0 ε ∈ H 1 (ω), which is independent of x 3 , and p 1 ε ∈ L 2 (Ω ε ). Moreover, the following estimates hold

p 0 ε H 1 (ω) ≤ Cε -3 2 ∇p ε H -1 (Ωε) 3 , p 1 ε L 2 (Ωε) ≤ C ∇p ε H -1 (Ωε) 3 .
(5.30)

From this result, we are able to give the convergence result for p ε . We denote by p 1 ε the rescaled function associated with p 1 ε , defined by p 1 ε (y) = p 1 ε (y , εy 3 ) for a.e. y ∈ Ω ε . Corollary 1. Previous result implies the existence of p ∈ H 1 (ω) and p 1 ∈ L 2 (Ω) satisfying

ε 2 p 0 ε p in H 1 (ω), ε p 1 ε|Ω p 1 in L 2 (Ω), (5.31)
and moreover

ε 2 p ε|Ω → p in L 2 (Ω). (5.32)
Proof. From (5.30) and (5.17) 2 , we get

p 0 ε H 1 (ω) ≤ Cε -2 , p 1 ε L 2 (Ωε) ≤ Cε -1 2 , (5.33)
and after rescaling p 1 ε , last inequality becomes

p 1 ε L 2 ( Ωε) ≤ Cε -1 .
Previous estimates and the fact that Ω ⊂ Ω ε and | Ω ε \Ω| → 0 imply (5.31). The strong convergence (5.32) for the complete pressure p ε is a direct consequence of (5.31) and the decomposition (5.29).

Unfolding method

In order to capture the behaviour of u ε , w ε and p 1 ε (introduced in Proposition 2) near the rough boundary Γ ε , we need to introduce a new change of variables, which is adapted from the unfolding method (see [START_REF] Arbogast | Derivation of the double porosity model of single phase flow via homogenization theory[END_REF][START_REF] Casado-Díaz | Asymptotic behavior of the Navier-Stokes system in a thin domain with Navier condition on a slightly rough boundary[END_REF][START_REF] Cioranescu | Periodic unfolding and homogenization[END_REF]). To do this, for u ε ,

w ε ∈ H 1 (Ω ε ) 3 satisfying boundary conditions (2.21)- (2.23), p 1 ε ∈ L 2 (Ω ε )
, and ρ > 0, we set u ε , w ε and p 1 ε by

u ε (x , z) = u ε ε κ x ε + ε z , ε z 3 , (5.34)
w ε (x , z) = w ε ε κ x ε + ε z , ε z 3 , (5.35
)

p 1 ε (x , z) = p 1 ε ε κ x ε + ε z , ε z 3 , (5.36)
for a.e. (x , z) ∈ ω ρ × Z ε , where ω ρ is defined by (3.1) and

Z ε = z ∈ Z × R : -ε δ-Ψ(z • e 1 ) < z 3 < ε 1-h .
Remark 5. For every k ∈ I ρ,ε , the functions u ε , w ε and p 1 ε restricted to C k ε × Z ε are independent of x . However, as functions depending on z, they are obtained from their original functions by means of

z = x -ε k ε , z 3 = x 3 ε , (5.37)
that converts Q k ε in Z ε . Thanks to the estimates satisfied by u ε , w ε and p 1 ε given in (5.15), (5.17) and (5.33) 2 , respectively, we have the following compactness results. Lemma 2. Consider two sequences u ε , w ε ∈ V ε satisfying (5.15) 2 and (5.16) 2 , respectively. Define u ε , w ε ∈ H 1 ( Ω ε ) 3 by (4.10), so that (5.23) and (5.24) hold. Set δ ≤ 3 2 -1 2 . Then,

(i) If δ < 3 2 -1 2 , it holds u 1 (x , 0)∂ y 1 Ψ(z • e 1 ) = 0, a.e. (x , z ) ∈ ω × Z , (5.38
)

w 1 (x , 0)∂ y 1 Ψ(z • e 1 ) = 0, a.e. (x , z ) ∈ ω × Z .
(5.39)

(ii) If δ = 3 2 -1 2 , there exist u, w ∈ L 2 (ω, V 3 ), where V is the space of functions ϕ : R 2 × (0, +∞) → R such that ϕ ∈ H 1 # ( Q M ), for every M > 0, and ∇ ϕ ∈ L 2 # ( Q) 3 , satisfying u 3 (x , z , 0) = -λ∂ y 1 Ψ(z • e 1 ) u 1 (x , 0), a.e. (x , z ) ∈ ω × Z , (5.40)
w 3 (x , z , 0) = -λ∂ y 1 Ψ(z • e 1 ) w 1 (x , 0), a.e. (x , z ) ∈ ω × Z , (5.41)
and such that, for any ρ, M > 0, the sequences u ε and w ε , respectively given by (5.34) and (5.35), satisfy

ε 1- 2 D z u ε D z u in L 2 (ω ρ × Q M ) 3×3 , ε 3- 2 D z w ε D z w in L 2 (ω ρ × Q M ) 3×3 . (5.42)
Moreover, if one assumes div u ε = 0 in Ω ε , then u satisfies

div z u = 0 in ω × Q. (5.43)
Proof. This result is a direct consequence of [START_REF] Casado-Díaz | Asymptotic behavior of the Navier-Stokes system in a thin domain with Navier condition on a slightly rough boundary[END_REF]Lemma 5.4], applied to the sequences ε 2 u ε and ε 3 w ε . Lemma 3. We consider p 1 ε ∈ L 2 (Ω ε ) such that (5.33) 2 holds. Then, there exists p 1 ∈ L 2 (ω × Q) satisfying, up to a subsequence, the convergence

ε 1+ 2 p 1 ε p 1 in L 2 ω ρ × Q M ∀ρ, M > 0.
(5.44)

Proof. It is a direct consequence of [25, Lemma 5.5], applied to the sequence ε 2 p 1 ε .

We are now in position to prove Theorem 2 in the critical case δ = 3 2 -1 2 , which is the most relevant from the mechanical point of view since it describes the coupling effects between the riblets and nonzero boundary conditions.

Proof of Theorem 2. Let us consider δ = 3 2 -1 2 with > 1. First of all, Lemma 1 and Corollary 1 implies the existence of u , w ∈ H 1 (0, h; L 2 (ω)) 2 such that (4.13), and p ∈ H 1 (ω) so that convergences of u ε , w ε and p ε given in (4.11) hold. Also, we have that the divergence condition (4.12) 3 holds. From Corollary 1, the sequences p 0 ε and p1 ε satisfy convergences given in (5.31).

We recall the variational formulation given by (2.28) and (2.29). For ϕ, ψ ∈ V ε , (u ε , w ε , p ε ) satisfies

Ωε Du ε : Dϕ dx - Ωε p ε div ϕ dx -2N 2 Ωε w ε • (∇ × ϕ) dx -2 1 α -N 2 Γ 0 ε (w ε × n ε ) • ϕ dσ = 0 , (5.45)
ε 2 R c Ωε Dw ε : Dψ dx -2N 2 β Γ 0 ε (u ε × n ε ) • ψ dσ + 4N 2 Ωε w ε • ψ dx -2N 2 Ωε (∇ × u ε ) • ψ dx = 0 .
(5.46)

Now, we want to pass to the limit in the above variational formulations. To do this, we will use appropriate test functions ϕ ε , ψ ε . We divide the proof in four steps.

Step 1. Definition of the test functions. Lemma 2 gives the existence of u, w ∈ L 2 (ω; V 3 ) satisfying (5.40), (5.41) and (5.43). Thanks to this, we consider the following test functions. For any ϕ, ψ ∈ C 1 c (ω × (-h, h)) 3 , with

ϕ 3 = ψ 3 = 0, ϕ, ψ ∈ C 1 c (ω; C 1 # ( Q)) 3 satisfying              D z ϕ(x , z) = 0 a.e. in {z 3 > M } for some M > 0, ϕ (y , y 3 ) = ϕ (y , 0) when y 3 ≤ 0, ϕ(x , z , z 3 ) = ϕ(x , z , 0) when z 3 ≤ 0, λ ∂ z 1 Ψ(z • e 1 ) ϕ 1 (y , 0) + ϕ 3 (y , z , 0) = 0, (5.47)                D z ψ(x , z) = 0 a.e. in {z 3 > M } for some M > 0, ψ (y , y 3 ) = ψ (y , 0) when y 3 ≤ 0, ψ(x , z , z 3 ) = ψ(x , z , 0) when z 3 ≤ 0, λ ∂ z 1 Ψ(z • e 1 ) ψ 1 (y , 0) + ψ 3 (y , z , 0) = 0, (5.48) and a function ζ ∈ C ∞ (R) such that ζ(s) =      1 when s < 1 3 , 0 when s > 2 3 , (5.49) we set ϕ ε , ψ ε ∈ H 1 (Ω ε ) 3 as follows        ϕ ε (x) = ε ϕ x , x 3 ε + ε 1+ 2 ϕ x , x ε ζ x 3 ε , ϕ ε,3 = ε 1+ 2 ϕ 3 x , x ε ζ x 3 ε -ε ϕ 1 x , x ε λ ∂ x 1 Ψ 1 ε x • e 1 ζ x 3 ε ,        ψ ε (x) = ψ x , x 3 ε + ε -1 2 ψ x , x ε ζ x 3 ε , ψ ε,3 = ε -1 2 ψ 3 x , x ε ζ x 3 ε -ε -1 ψ 1 x , x ε λ ∂ x 1 Ψ 1 ε x • e 1 ζ x 3 ε .
Since ϕ (x), ψ (x), ϕ (x , z) and ψ(x , z) are zero when x is out of a compact subset of ω, (5.47) and (5.48), then ϕ ε , ψ ε are such that

ϕ ε = ψ ε = 0 on ∂Ω ε \ Γ 0 ε , ϕ ε • n ε = ψ ε • n ε = 0 on Γ 0 ε .
So, we are able to consider ϕ ε and ψ ε , respectively, as test functions in (5.45) and (5.46). The difficulty now is to obtain the limit of every terms of (5.45) and (5.46). For this, we observe that from the conditions D z ϕ = D z ψ = 0 a.e. in {z 3 > M } and (5.49), it follows

ϕ ε (x) = ε ϕ x , x 3 ε , 0 + g ε in Ω ε , (5.50)
Dϕ ε (x) = 2 i=1 ∂ y 3 ϕ i x , x 3 ε e i ⊗ e 3 + ε 1- 2 D z ϕ x , x ε + h ε (x) in Ω ε , (5.51)
ψ ε (x) = ψ x , x 3 ε , 0 + gε in Ω ε , (5.52
)

Dψ ε (x) = ε -1 2 i=1 ∂ y 3 ψ i x , x 3 ε e i ⊗ e 3 + ε -1+ 2 D z ψ x , x r ε + hε (x) in Ω ε , (5.53)
where

g ε , gε ∈ C 0 (Ω ε ) 3 , h ε , hε ∈ C 0 (Ω ε ) 3×3 (thanks to > 1) are such that ε -3 Ωε |g ε | 2 dx ≤ C ε +1 + ε 3( -1) = O ε , (5.54)
ε -2 Γ 0 ε |g ε | 2 dσ ≤ Cε -1 = O ε , (5.55)
ε -1 Ωε |h ε | 2 dx ≤ Cε 3 ε -2 + ε -4 + 1 ε = O ε , (5.56)
ε -1 Ωε |g ε | 2 dx ≤ C ε -1 + ε 3( -1) = O ε , (5.57)
Γ 0 ε |g ε | 2 dσ ≤ Cε -1 = O ε , (5.58) ε Ωε | hε | 2 dx ≤ Cε 3 ε -2 + ε -4 + 1 ε = O ε . (5.59)
We remark that functions g ε , gε , h ε and hε and previous estimates are devoted to identify terms of the variational formulation that are negligible in the asymptotic analysis.

Step 2. Passing to the limit in variational formulation (5.45). We can pass to the limit in every term of (5.45).

• 1st term of (5.45): From (5.15), (5.51) and (5.56), we deduce

Ωε Du ε : Dϕ ε dx = Ω + ε ∂ x 3 u ε (x) • ∂ y 3 ϕ x , x 3 ε dx +ε 1- 2 Ω + ε Du ε (x) : D z ϕ x , x ε dx + O ε .
(5.60)

Observe that main order terms are defined in Ω + ε = ω × (0, εh). To obtain this, here we have used that, in

Ω - ε = ω × (-Ψ ε (x ), 0), ε 1- 2 Ω - ε Du ε (x) : D z ϕ x , x ε dx = O ε ,
Last estimate results from Cauchy-Schwarz inequality, the estimate of Du ε given in (5.15) and the estimate

Ω - ε |ε 1- 2 D z ϕ| 2 dx ≤ Cε 1-|Ω - ε | = Cε 1-+δ = Cε 1+ 2 .
Also, we have used that

Ωε Du ε : h ε dx ≤ Du ε L 2 (Ωε) 3×3 h ε L 2 (Ωε) 3×3 ≤ Cε -1 2 h ε L 2 (Ωε) 3×3 = O ε .
Now, using the dilatation (4.8) and (4.11) 1 , we get

Ω + ε ∂ x 3 u ε (x) • ∂ y 3 ϕ x , x 3 ε dx = Ω ∂ y 3 u ε (y) • ∂ y 3 ϕ (y) dy = Ω ∂ y 3 u (y) • ∂ y 3 ϕ (y) dy + O ε .
Next, from the unfolding (5.37), the hypothesis of the support of D z ϕ and (5.42) 1 , we deduce

ε 1- 2 Ω + ε Du ε (x) : D z ϕ x , x ε dx = ω× Q M D z ε 1- 2 u ε (x , z) : D z ϕ(x , z) dx dz + O ε = ω× Q D z u(x , z) : D z ϕ(x , z) dx dz + O ε ,
(5.61) where Q M = Z × (0, M ) and Q = Z × (0, +∞). For more details of the unfolding change, we refer to [START_REF] Casado-Díaz | Asymptotic behavior of the Navier-Stokes system in a thin domain with Navier condition on a slightly rough boundary[END_REF]. Then, we have that (5.60) is given by

Ωε Du ε : Dϕ dx = Ω ∂ y 3 u (y) • ∂ y 3 ϕ (y) dy + ω× Q D z u(x , z) : D z ϕ(x , z) dx dz + O ε . (5.62)
• 2nd term of (5.45). Using the decomposition (5.29), (5.50) and (5.54), we have that

- Ωε p ε div ϕ ε dx = Ωε ∇ x p 0 ε (x ) • ϕ ε (x) dx - Ωε p 1 ε (x) div ϕ ε (x) dx .
Applying the change of variables (4.8) and (5.31) 1 to the first integral, we get

Ωε ∇ x p 0 ε (x ) • ϕ ε (x) dx = ε Ω + ε ∇ x p 0 ε (x ) • ϕ x , x 3 ε dx + O ε = Ω ε 2 ∇ y p 0 ε (y ) • ϕ (y) dy + O ε = Ω ∇ y p(y ) • ϕ (y) dy + O ε .
Here, we have used that

Ωε ∇ x p 0 ε (x )g ε (x)dx ≤ ∇ x p 0 ε (x ) L 2 (Ωε) g ε L 2 (Ωε) ≤ Cε -3 2 g ε L 2 (Ωε) = O ε .
For the second integral, using the change of variables (5.37) and (5.44) (see [START_REF] Casado-Díaz | Asymptotic behavior of the Navier-Stokes system in a thin domain with Navier condition on a slightly rough boundary[END_REF] for more details), we obtain

Ωε p 1 ε (x)div ϕ ε (x) dx = ε 1- 2 Ω + ε p 1 ε (x) div z ϕ x , x ε dx + O ε = ω× Q M ε 1+ 2 p 1 ε (x , z) div z ϕ(x , z) dx dz + O ε = ω× Q p 1 (x , z) div z ϕ(x , z) dx dz + O ε .
Then, we get

- Ωε p ε div ϕ ε dx = Ω ∇ y p(y ) • ϕ (y) dy - ω× Q p 1 (x , z) div z ϕ(x , z) dx dz + O ε .
(5.63)

• 3rd term of (5.45). Using (5.51), (5.56), the change of variables (4.8) and (5.24), we have

-2N 2 Ωε w ε • (∇ × ϕ ε ) dx = -2N 2 Ω + ε w ε (x) • rot y 3 ϕ x , x 3 ε dx + O ε = -2N 2 Ω ε w ε (y) • rot y 3 ϕ (y) dy + O ε = -2N 2 Ω w (y) • rot y 3 ϕ (y) dy + O ε .
(5.64)

Among others, here we have used that

Ωε w ε • (∇ × g ε )(x) dx ≤ w ε L 2 (Ωε) 3 h ε L 2 (Ωε) 3×3 ≤ Cε -1 2 h ε L 2 (Ωε) 3×3 = O ε .
Integrating by parts by using the formula (2.27), we then get

-2N 2 Ωε w ε • (∇ × ϕ ε ) dx = -2N 2 Ω rot y 3 w (y) • ϕ (y) dy -2N 2 Γ [ w ] ⊥ (y) • ϕ (y) dσ + O ε . (5.65)
• 4th term of (5.45). From w ε = 0 on Γ 1 ε and estimate (5.16), we obtain

Γ 0 ε |w ε | 2 dσ ≤ Cε Ωε |Dw ε | 2 dx ≤ Cε -2 .
Then, from (5.50), (5.55) and (5.47) 2 , we deduce

-2 1 α -N 2 Γ 0 ε (w ε × n ε ) • ϕ ε dσ = -2 1 α -N 2 ε ω w ε x , -λε 3 -1 2 Ψ x 1 ε × n ε • ϕ(x , 0) 1 + λ 2 ε -1 2 2 ∂ x 1 Ψ 1 ε x • e 1 2 dx + O ε = -2 1 α -N 2 ε ω w ε x , -λε 3 -1 2 Ψ x ε × n • ϕ(x , 0) dx + O ε ,
where n = (0, 0, -1). Here we have used that

Γ 0 ε (w ε × n ε ) • ϕ ε dσ ≤ w ε L 2 (Γ 0 ε ) 3 g ε L 2 (Γ 0 ε ) 3 ≤ Cε -1 g ε L 2 (Γ 0 ε ) 3 = O ε .
By means of integration in the variable x 3 , we get

ω εw ε x , -λε 3 -1 2 Ψ x 1 ε -εw ε (x , 0) 2 dx ≤ Cε 3 +3 2 Ωε |Dw ε | 2 dx ≤ Cε 3 -1 2 .
(5.66)

Then, from w ε (x , 0) = w ε (x , 0) and (5.24), we obtain

-2 1 α -N 2 Γ 0 ε (w ε × n ε ) • ϕ ε dσ = -2 1 α -N 2 ω (ε w ε (y , 0) × n) • ϕ(y , 0) dy + O ε = -2 1 α -N 2 ω [ w ] ⊥ (y , 0) • ϕ (y , 0) dy + O ε .
(5.67) By considering (5.62), (5.63), (5.65) and (5.67), then we get that u , w , p, u and p 1 satisfy

Ω ∂ y 3 u (y) • ∂ y 3 ϕ (y) dy + ω× Q D z u(x , z) : D z ϕ(x , z) dx dz + Ω ∇ y p(y ) • ϕ (y) dy - ω× Q p 1 (x , z) div z ϕ(x , z) dx dz -2N 2 Ω rot y 3 w (y) • ϕ (y) dy - 2 α ω [ w ] ⊥ (y , 0) • ϕ (y , 0) dy = 0 (5.68) for any ϕ ∈ C 1 c (ω × (-h, h)) 2 , ϕ ∈ C 1 c (ω; C 1 # ( Q)) 3 satisfying (5.

47), and then, by arguments of density, for any

ϕ ∈ H 1 (0, h; L 2 (ω)) 2 and any ϕ ∈ L 2 (ω; V) 3 satisfying ϕ(x , h) = 0 a.e. x ∈ ω, λ∂ z 1 Ψ(z • e 1 ) ϕ 1 (x , 0) + ϕ 3 (x , z , 0) = 0 a.e. (x , z ) ∈ ω × Z .
By considering ϕ = 0 in (5.68), we get that ( u,

p 1 ) ∈ V 3 × L 2 # ( Q) solves              -∆ z u + ∇ z p 1 = 0 in R 2 × R + , div z u = 0 in R 2 × R + , u 3 (x , z , 0) = -λ∂ z 1 Ψ(z • e 1 ) u 1 (x , 0) on R 2 × {0}, ∂ z 3 u = 0 on R 2 × {0}, (5.69)
a.e. x in ω. Setting ( φ 1,λ , q 1,λ ) by (4.16), we derive

D z u(x , z) = -u 1 (x , 0)D z φ 1,λ (z), a.e. in R 2 × R + , p 1 (x , z) = -u 1 (x , 0) q 1,λ (z), a.e. in R 2 × R + .
(5.70)

Next, by considering ϕ ∈ H 1 (0, h; L 2 (ω)) 2 , ϕ (x , h) = 0, a.e. x ∈ ω, as test function in (5.68), setting ϕ by ϕ(x , z) = -ϕ 1 (x , 0) φ 1,λ (z), and taking into account (5.70), we obtain

Ω ∂ y 3 u (y) • ∂ y 3 ϕ (y) dy + Ω ∇ y p(y ) • ϕ (y) dy -2N 2 Ω rot y 3 w (y) • ϕ (y) dy +E λ ω u 1 (y , 0) ϕ 1 (y , 0) dy - 2 α ω [ w ] ⊥ (y , 0) • ϕ (y , 0) dy = 0, (5.71)
where E λ ∈ R is given by (4.15) 1 .

Step 3. Passing to the limit in the variational formulation (5.46). This step is similar to the previos step, so we will only give some details.

• 1st term of (5.46). Analogously to the first step of the previous variational formulation, we deduce

ε 2 R c Ωε Dw ε : Dψ ε dx = εR c Ω + ε ∂ x 3 w ε (x) • ∂ y 3 ψ x , x 3 ε dx +ε 3- 2 R c Ω + ε Dw ε (x) : D z ψ x , x ε dx + O ε , (5.72)
and by using the changes of variables, we can pass to the limit in every terms by obtaining

ε 2 R c Ωε Dw ε : Dψ ε dx = R c Ω ∂ y 3 w (y) • ∂ y 3 ϕ (y) dy + R c ω× Q D z w(x , z) : D z ψ(x , z) dx dz + O ε .
(5.73) • 2nd term of (5.46). From u ε = 0 on Γ 1 ε and the estimate (5.15), we get

Γ 0 ε |u ε | 2 dσ ≤ Cε Ωε |Du ε | 2 dx ≤ C.
Using this, and proceeding analogously to the development of the fourth term in (5.45), we get

-2N 2 β Γ 0 ε (u ε × n ε ) • ψ ε dσ = -2N 2 β ω [ u (y , 0)] ⊥ • ψ (y , 0) dy + O ε .
(5.74)

• 3rd term of (5.46). Applying (5.57), the change of variables (4.8) and convergence (5.24), we have 4N 2

Ωε

w ε • ψ ε dx = 4N 2 Ω w ε (x) • ψ x , x 3 ε dx + O ε = 4N 2 Ω ε w ε (y) • ψ (y) dy + O ε = 4N 2 Ω w (y) • ψ (y) dy + O ε .
(5.75)

• 4th term of (5.46). Similarly to (5.64) and taking into account convergence (5.23), we have

-2N 2 Ωε (∇ × u ε ) • ψ ε dx = -2N 2 Ω rot y 3 u (y) • ψ (y) dy + O ε .
(5.76) Finally, from (5.73) to (5.76), we have obtained

R c Ω ∂ y 3 w (y) • ∂ y 3 ψ (y) dy + R c ω× Q D z w(x , z) : D z ψ(x , z) dx dz + 4N 2 Ω w (y) • ψ (y) dy -2N 2 Ω rot y 3 u (y) • ψ (y) dy -2N 2 β ω [ u (y , 0)] ⊥ • ψ (y , 0) dy = 0, (5.77) for any ψ ∈ C 1 c (ω × (-1, 1)) 2 , ψ ∈ C 1 c (ω; C 1 # ( Q)) 3 satisfying (5.

47), and by argument of density, for any

ψ ∈ H 1 (0, 1; L 2 (ω)) 2 , ψ ∈ L 2 (ω; V) 3 satisfying ψ(x , 1) = 0, a.e. x ∈ ω, λ∂ 1 Ψ(z 1 ) ψ 1 (x , 0) + ψ 3 (x , z , 0) = 0, a.e. (x , z ) ∈ ω × Z .
Similarly tot he previous step, we eliminate w from (5.77). To do this, we consider ψ = 0 in (5.77), which gives that

w ∈ V 3 × L 2 # ( Q) 3 solves        -∆ z w = 0 in R 2 × R + , w 3 (x , z , 0) = -λ∂ 1 Ψ(z 1 ) w 1 (x , 0) on R 2 × {0}, ∂ z 3 w = 0 on R 2 × {0}, (5.78)
a.e. x in ω. Setting φ 2,λ by (4.17), we derive

D z w(x , z) = -w 1 (x , 0)D z φ 2,λ (z), a.e. in R 2 × R + , (5.79)
Next, by considering ψ ∈ H 1 (0, h; L 2 (ω)) 2 , ψ (x , h) = 0, a.e. x ∈ ω as test function in (5.77), defining ψ by ψ(x , z) = -ψ 1 (x , 0) φ 2,λ (z), and from (5.79), we deduce

R c Ω ∂ y 3 w (y) • ∂ y 3 ψ (y) dy + 4N 2 Ω w (y) • ψ (y) dy -2N 2 Ω rot y 3 u (y) • ψ (y) dy +R c F λ ω w 1 (y , 0) ψ 1 (y , 0) dy -2N 2 β ω [ u (y , 0)] ⊥ • ψ (y , 0) dy = 0 , (5.80)
where F λ is given by (4.15) 2 .

Step 4. Conclusion. Since ϕ and ψ are arbitrary, we derive from (5.71) and (5.80) that ( u , w , p) satisfies the system (4.12) 1,2 with boundary conditions (4.14). To ensure that the whole sequences u ε , w ε , ε p ε converge, it remains to prove the existence and uniqueness of weak solution of the effective system (4.12)- (4.14).

Note that we can always reduce the non vanishing divergence problem (4.12)-(4.14) to a free divergence problem by considering the lift function J ∈ H 1 (Ω) 3 such that

       div y J = 0 in Ω, J • n = 0 on Γ , J = -Se 3 on Γ 1 , J = 0 on Γ,
and hence using the change of unknowns U = u -J . Therefore, it is sufficient to study the existence and uniqueness of the weak solution of the problem

           -∂ 2 y 3 u + ∇ y p -2N 2 rot y 3 w = ∂ 2 y 3 J in Ω, -R c ∂ 2 y 3 w + 4N 2 w -2N 2 rot y 3 u = 2N 2 rot y 3 J in Ω, div y h 0 u (y , y 3 ) dy 3 = 0 in ω.
with the boundary conditions u = 0, w = 0 on ω × {h},

∂ y 3 u = - 2 α [ w ] ⊥ + E λ (u • e 1 ) e 1 on Γ, R c ∂ y 3 w = -2N 2 β[ u ] ⊥ + R c F λ (w • e 1 ) e 1 on Γ .
The existence and uniqueness of this problem follows the lines of the proof of Theorem 1 with a flat bottom, taking into account suitable spaces and the new boundary terms E λ u • e 1 and R c F λ w • e 1 and the new source terms ∂ 2 y 3 J and 2N 2 rot y 3 J in the variational formulation.

The proofs of Lemmas 4, 5 and 6 and Corollary 2 are given in the Appendix. We finish this section by providing the proofs of Theorems 3 and 4.

Proof of Theorem 3. Using (4.22), we obtain

- 1 0 h 0 u 1 (y 1 , y 3 ) dy 3 ∂ y 1 θ(y 1 ) dy 1 = 1 0 Sθ(y 1 ) dy 1 , ∀ θ ∈ H 1 ((0, 1)).
From Lemmas 4 and 5, by averaging (6.12) and (6.14) we obtain

h 0 u 1 (y 1 , y 3 ) dy 3 = -Θ λ ∂ y 1 p(y 1 ).
Proof of Theorem 4. Using (4.22), we obtain

- 1 0 h 0 u 1 (y 1 , y 3 ) dy 3 ∂ y 1 θ(y 1 ) dy 1 = 1 0 Sθ(y 1 ) dy 1 , ∀ θ ∈ H 1 ((0, 1)).
From Lemma 6, by averaging (6.20) 1 , we obtain h 0 u 1 (y 1 , y 3 ) dy 3 = -Θ∂ y 1 p(y 1 ).

Application to squeeze-film bearing

As a an application of the results presented in Section 4, we consider in this section a squeezefilm bearing composed of two parallel plates separated by a micropolar fluid film. The lower surface is at rest and composed of a rough material, while the smooth upper surface is under normal squeeze motion. Hence, the distance between the plates is a decreasing function of time, which is expected to go to zero as the fluid is squeezed out of the gap. We assume that the motion of the upper plate is slow, so that the inertial effects can be neglected and the behaviour of the bearing can be captured by a quasi-static model.

Derivation of the model

Let T 0 > 0 be the characteristic time of the motion. Following the notation from Section 2, we denote by L the horizontal dimension of the bearing and by hc(T 0 t) the distance between the plates at time t, where t stands for the dimensionless time variable and h > 0 is an adimensional constant. We introduce the small parameter ε(t) = c(T 0 t)/L and assume that, at each time t, the fluid flow is described by the system (2.5)-(2.12) with ε = ε(t), where all the constant ν, ν r , c a , c d , α, β are fixed. However, the velocity V of the upper plate is related to the load W applied on the bearing, which is assumed to be independant on t, through the implicit relation

W = Γ ε(t) p ε(t) . (6.1)
We assume that the rough bottom is composed of periodically distributed riblets, described by a given function Ψ of the form

Ψ(x 1 ) = LΛΨ 1 LM x 1 (6.2)
where Λ, M > 0 correspond respectively to the amplitude and period of the ribbed surface, divided by L, and where the function Ψ is 1-periodic and regular, and satisfies (2.3).

Using the notation introduced in Section 2, for a given time t, the fluid domain Ω ε(t) is given by Ω

ε(t) = (x , x 3 ) ∈ (Lω) × R , -Ψ ε(t) (x 1 ) < x 3 < h L ε(t) ,
where according to definition (2.1), the function Ψ ε(t) takes the form

Ψ ε(t) (x 1 ) = Lλ(t)ε(t) δ(t) Ψ 1 Lε(t) (t) x 1 . (6.3)
Since the geometry of the lower plate is, in fact, independent on time, Ψ ε(t) coincides with the fixed profile Ψ given by (6.2). Hence, parameters λ(t), δ(t) satisfy

λ(t) ε(t) δ(t) = Λ and ε(t) (t) = M . (6.4)
Considering the critical regime δ(t) = 3 2 (t)-1 2 , the parameters (t), δ(t), λ(t) can thus be expressed as functions of ε(t), M, L, and in particular the parameter λ(t) is given by λ(t) = Λ ε(t)

1 2 -3 2 ln M ln ε(t) , which simplifies to λ(t) = Λ M 3/2 ε(t) 1/2 . (6.5)
Analogously to (2.15), we introduce the dimensionless pressure p ε(t) defined by p ε(t) = L V 0 (ν+νr) p ε(t) , and accordingly, the adimensional load W = W LV 0 (ν+νr) . Defining also the rescaled pressure p ε(t) , the constraint (6.1) can be rephrased as ω p ε(t) (y , h) dy = W .

(6.6)

Since the system is independent on the x 2 -direction, we take ω = (0, 1) and apply Theorems 2 and 3 to approximate p ε(t) by p λ(t),S(t) /ε(t) 2 , where p λ(t),S(t) is the solution of Reynolds equation (4.26), with λ = λ(t) and S = S(t), S(t) being fixed consistently with relation (6.6). Indeed, for a given time t, if one assumes that ε(t) has been computed, all roughness parameters λ(t), (t), δ(t) are then prescribed by relations (6.4) and by the critical relation δ(t) = 3 2 (t) -1 2 . Setting λ = λ(t), = (t), δ = δ(t), the state of the system at time t is thus described by (2.18)-(2.25), where ε takes the value ε(t). All other parameters being fixed (equal to their values at instant t), we can define a sequence of problems (2.18)-(2.25), where ε ∈ (0, ε(t)) is the only parameter that is let to zero. Since δ and are bound by relation δ = 3 2 -1 2 , and ε(t) is small, we can apply the asymptotic analysis result from Theorems 2 and 3 to replace the solution of problem (2.18)-(2.25), with ε = ε(t), by the solution of the effective problem (4.12)-(4.14).

Let us point that this strategy aims at approximating the configuration of the system at each instant t, from a computational point of view, by the limit provided by Theorems 2 and 3. The study of the physical limit of the system when ε(t) goes to zero (which implies that λ(t) also goes to zero by relation (6.5)) is beyond the scope of the paper. In particular, the value of the physical parameter ε(t), albeit small, remains greater than the positive quantity ε 0 /2 throughout the simulations.

As is usual in the lubrication field, we impose Dirichlet boundary conditions for the pressure on x 1 ∈ {0, 1}, instead of the Neumann boundary conditions implicitly contained in the weak formulation (4.26). We obtain the relation

1 0 p λ(t),S(t) (y 1 ) dy 1 = ε(t) 2 W . (6.7)
Since Reynolds equation (4.26) is linear with respect to S, there holds p λ(t),S(t) = S(t) p λ(t),1 (where p λ(t),1 satisfies (4.26) with λ = λ(t) and S = 1) so that S(t) is given by

S(t) = ε(t) 2 W 1 0 p λ(t),1 (y 1 ) dy 1 . (6.8)
Using relations (4.2) and (2.15) 5 , the normal velocity of the upper plate at time t, which is given by -Lh T 0 ε (t), can thus be expressed by -Lh T 0 ε (t) = V 0 ε(t) S(t), yielding the differential equation

ε (t) = - T 0 V 0 W Lh ε(t) 3 1 0 p λ(t),1 (y 1 ) dy 1 , (6.9)
where λ(t) depends on ε(t) through relation (6.5). Since p λ(t),1 satisfies

-Θ λ(t) ∂ 2 y 1 p λ(t),1 = 1 in (0, 1) , p λ(t),1 (0) = p λ(t),1 (1) = 0 ,
is can be expressed by the explicit formula p λ(t),1 (y 1 ) = y 1 (1-y 1 )

2Θ λ(t) so the ODE (6.9) can be rewritten

ε (t) = - 12T 0 V 0 W Lh Θ λ(t) ε(t) 3 . (6.10)
Equations of motion. Denoting respectively by κ and χ the dimensionless constants

κ = Λ 2 M 3 , χ = 12T 0 V 0 W Lh
the model can be summarized by the system of equations

ε (t) = -χ Θ λ(t) ε(t) 3 λ(t) = κε(t) (6.11)
which needs to be completed with the definition of an initial state ε 0 > 0 such that ε(0) = ε 0 .

To facilitate the comparison with the numerical results presented in [START_REF] Bessonov | A new generalization of the Reynolds equation for a micropolar fluid and its application to bearing theory[END_REF][START_REF] Bayada | New models in micropolar fluid and their application to lubrication[END_REF][START_REF] Bonnivard | Effects of rough boundary and nonzero boundary conditions on the lubrication process with micropolar fluid[END_REF], we introduce the dimensionless parameter

δ = R c 2N 2 β
and define the relative viscosity coefficient νb by νb = ν b /ν, where ν b is the boundary viscosity and ν is the classical viscosity of the fluid. By definition of α (see (2.13)), this coefficient can be expressed as

νb = 1 -αN 2 1 -N 2 .
Hence, the solution t → ε(t) to the system (6.11) depends on the following set of parameters:

-parameters N, R c characterizing the physical properties of the micropolar fluid, -νb , δ characterizing the interaction between the micropolar fluid and the upper wall, -E, κ related to the geometry of the rough pattern, -h, χ associated with characteristic dimensions of the model, -and the initial datum ε 0 . Since we are mostly interested in the combined effects of micropolarity and roughness, we will simply put h = 1 and χ = 1 in the sequel.

In the aim of estimating the influence of the parameters on the performance of the squeeze-film bearing, for each set of parameters, we solve numerically the associated system (6.11) (with χ = 1) by a second-order Runge-Kutta method, and compute the "half-life time" T half , i.e. the first instant t such that ε(t) < ε 0 /2, which means that the width of the bearing has been divided by two. In our analysis, the configurations giving rise to the highest values of T half will be considered the most efficient from a mechanical perspective.

Numerical results

Determination of the initial width ε 0

In order to ensure the stability of our numerical method, we choose to determine a small initial value ε 0 that guarantees that, for each tested set of parameters, the computed value of ε(t) decreases during the simulation, which means by the ODE ε (t) = -Θ λ(t) ε(t) 3 that the function Θ λ(t) should remain positive. To this aim, we take advantage of the asymptotic development of function Θ λ as λ goes to zero, given in Corollary 2 : Table 1: Second derivative of λ → Θ λ at λ = 0, computed using the exact formula of Θ λ or the asymptotic development Θ0 -CjEλ 2 Θ1 , for the set of parameters given in Fig. 2.

Θ λ = Θ 0 -C j Eλ 2 Θ 1 + O(λ 4 ) , 36 0 0.1 0.2 0.3 0.4 0.5 -2 -1.5 -1 -0.5 0 0.5 λ Case α = 1 Θ λ Θ 0 -C α Eλ 2 Θ 1 0 0.1 0.2 0.3 0.4 0.5 -0.1 0 0.1 0.2 λ Case α = 1 Θ λ Θ 0 -C α Eλ 2 Θ 1
νb ∂ 2 Θ λ ∂λ 2 |λ=0 -2C j EΘ 1 0.1 -17.7344 -17.7346 1 -2.97108 -2.97109 with j = α if α = 1, j = N if α = 1.
As a validation of the above formula, we have plotted in Fig. 2 the exact quantity Θ λ (given by (4.24)) and its approximation Θ 0 -C j Eλ 2 Θ 1 as functions of λ, using 2 different sets of parameters corresponding respectively to a case where α = 1 and α = 1. We have also computed numerically the second derivative

∂ 2 Θ λ ∂λ 2 |λ=0
and, in each case,

given the expected value -2C j EΘ 1 (see Table 1).

For each set of parameters, we use the approximation Θ λ ≈ Θ 0 -C α Eλ 2 Θ 1 to estimate an initial value ε 0 such that the associated value of λ, given by λ 0 = (κε 0 ) 1/2 (see (6.5)), satisfies Θ λ 0 > 0. Keeping in mind that ε 0 should be small so that the Reynolds equation gives a good approximation of the pressure, we introduce a threshold ε max and set

ε 0 = min Θ 0 κC j E Θ 1 , ε max .
We check a posteriori that the solution ε(t) to the ODE (6.9) is indeed a decreasing function of time.

6.2.2. Influence of parameters N, R c , νb , δ, E Since the model under study depends on many parameters, in order to perform comparisons, we have chosen to unify the presentation of the numerical results by plotting the half-life time T half as a function of N ∈ [0, 0.7] (which ensures that condition N 2 ≤ 1/2 is fulfilled), after normalization by its value for N = 0, for different values of R c ∈ {0.025, 0.05, 0.1, 0.2}, using various sets of parameters νb , δ, E.

Influence of νb . We have plotted in Fig. 3 and Fig. 4 the results obtained with E = 0 and E = 10 respectively, considering different values of parameter νb ∈ {0.05, 0.1, 0.2, 0.4} and for a fixed δ = 1. It appears that modifying νb does not have much of an impact of the computed value of T half , at least qualititatively. Consequently, we have decided to fix this parameter and impose νb = 0.1 in the rest of the simulations.

Influence of N, R c . On the opposite, T half is very sensitive to the couple of parameters (N, R c ). It appears that increasing N from the initial value N = 0 leads, at first, to a slow increasing of T half , up to a certain value of N . Then, the dependence of T half on N is strongly affected by the value of R c . For instance, in the case E = 10 (Fig. 4), the maximal value R c = 0.2 leads to a gradual increasing of the slope of the curve, up to N = 0.7, whereas the minimal value R c = 0.025 corresponds to a reduction of T half as N increases, ending up with a division by a factor 2 with respect to the initial value for N = 0.

Influence of E. In order to estimate the possible values of E encountered in practical applications, we consider the example of three riblet profiles that are often used in the engineering literature, namely, the V -shape, the U -shape and the blade riblets. We have plotted in Fig. 7 the corresponding Ψ functions, normalized by (2.3), and in Fig. 8, a possible rough geometry defined as the graph of the function -Ψ ε , defined by (6.3) and (6.4). To compute the energy E associated with each of these riblet profiles, we have solved system (4.16) (in case λ = 1) by a finite element method using FreeFem++ software [START_REF] Hecht | New Development in FreeFem++[END_REF], implementing a Taylor-Hood approximation for the velocity-pressure pair, i.e., P 2 elements for the velocity field and P 1 elements for the pressure. The results are summarized in Table 2.

Even though the computed values of E associated with the previous examples of riblets are of order 10 to 100, it turns out that the results of our simulations are very stable when E exceeds 10, so we have chosen to represent in Fig. 5 the range of values of E that produces the most significant changes in the behaviour of the model, which is E ∈ [0, 10] for δ = 1 and νb = 0.1. In the case E = 0, i.e. in the absence of roughness, raising N from N = 0 results at first in a slight increase of T half . Then, this behaviour reverses and for larger values of N , we observe a significant reduction of T half , for all the tested values of R c . Raising the value of E produces a visible change in the model, whose behaviour gets more and more dependent on parameter R c . As E increases, the function N → T half becomes monotonic for the highest values of R c (R c = 0.1 and R c = 0.2). This means that for certain micropolar fluids, the roughness of the upper plate may contribute to enhance the performance of the squeeze-film bearing, in the sense that T half becomes larger.

Influence of δ. We conclude this numerical study by investigating the impact of the slip length δ on the behaviour of the bearing. We have used values from δ = 0.7 to δ = 10, a range for which a significant change occurs in the simulations (see Fig. 6). We observe that small values of δ favor an enhancement of the performance of the bearing, for R c ∈ {0.1, 0.2}. On the opposite, large values of δ such as δ = 10 will typically reduce the performance of the micropolar fluid lubrication, with respect to lubrication with a Newtonian fluid, with the exception of the largest value or parameter R c = 0.2, where a slight increase of T half occurs.

Interpretation in terms of pressure using Corollary 2

As presented in Subsection 6.1, we have assumed in this study that λ and ε satisfy a relation of the form λ = (κε) 1/2 , with κ = Λ 2 /M 3 (see (6.5)), where ε is small. Consequently, λ itself is a small parameter, so one can apply the asymptotic developments obtained in Corollary 2 to shed a Riblet profile V -shape U -shape Blade Value of E 12.12 62.85 93.24

Table 2: Value of E associated with three riblet profiles given in Fig. 7: the V -shape, U -shape and blade riblets. light on the observed behaviour of the model. Indeed, for a small value of λ, equation (6.8) shows that the vertical velocity of the upper plate is inversely proportional to 1 0 p λ,1 (y 1 ) dy 1 , which can be approximated by

1 0 p λ,1 (y 1 ) dy 1 ≈ 1 0 p 0 (y 1 ) dy 1 + C α Eλ 2 p 1 (y 1 ) dy 1 ≈ 1 0 (1 + C α E Θ 1 Θ 0 λ 2 )p 0 (y 1 )dy 1 ,
where p 0 , p 1 are the respective solutions of (6.17) and (6.18) with S = 1. In particular, at first order in λ, there holds (in case α = 1)

1 0 p λ,1 (y 1 ) dy 1 - 1 0 p 0 (y 1 ) dy 1 1 0 p 0 (y 1 ) dy 1 ≈ C α E Θ 1 Θ 0 λ 2 .
As a result, the behaviour of the system appears to be driven by the factor C α E Θ 1 Θ 0 , in the sense that increasing this quantity should result in reducing the velocity of the upper plate subject to a given load, thus increasing T half and enhancing the performance of the bearing. However, as shown in Figs. 9 and 10, the behaviour of T half as a function of N and R c seems more appropriately described by the ratio Θ 1 Θ 0 than by the quantity C α E Θ 1 Θ 0 itself. This suggests that computing the quotient Θ 1 Θ 0 using the explicit formulas provided by Corollary 2 may be a valuable tool to predict and compare the relative performance of squeeze-film bearings, depending on the physical parameters associated with the system.

Conclusion

We may conclude from these simulations that, at least for a moderate range of N values and for small values of slip length δ, the model predicts an enhancement of the performance of the bearing lubricated by a micropolar fluid, with respect to a bearing lubricated by a Newtonian fluid. This feature is generally affected by the value of parameter R c , in the sense that large values of R c favor the enhancement of the performance of the squeeze-film bearing. Also, the introduction of a rough geometry for the lower plate, characterized by roughness parameter E, usually results in enlarging the range of N values for which lubrication with a micropolar fluid gets more efficient than lubrication with a Newtonian fluid.

ũ1 (y 1 , y 3 ) = 2N 2 k sh(ky 3 ) -η λ sh(kh) + γ α (y 3 -η λ h) -(1 -η λ ) γ α + 2N 2 kh sh(kh) A + 2N 2 k ch(ky 3 ) -η λ ch(kh) + (1 -η λ ) 2N 2 k -1 + y 3 h (1 -ch(kh)) B + 1 2(1-N 2 ) y 2 3 -h 2 + (1 -η λ )(y 3 h + h 2 ) ∂ y 1 p(y 1 ) , (6.12) w 2 (y 1 , y 3 ) = ch(ky 3 ) + γα 2 -(1 -η λ ) γα 2 + N 2 kh sh(kh) A + sh(ky 3 ) + (1 -η λ )(1 -ch(kh)) N 2 kh B + 1 2(1-N 2 ) y 3 + (1 -η λ ) h 2 ∂ y 1 p(y 1 ) ,
where

k = 2N 1 -N 2 R c , γ α 2 = 1 -αN 2 α -1 , η λ = 1 + αh α -1 E λ -1 , (6.13)
A = L 2(1-N 2 ) h 2 (1 + η λ ) 4N 4 η λ (1 -ch(kh)) + Rc β k 2 -ksh(kh) + (1 -η λ )(1 -ch(kh)) N 2 h Rc β -2N 2 h 2 η λ , B = L 2(1-N 2 ) N 2 h(1 + η λ )η λ γ α kh + 2N 2 sh(kh) +k Rc β -2N 2 h 2 η λ ch(kh) + γα 2 -(1 -η λ ) γα 2 + N 2 kh sh(kh) , L = -γα 2 + ch(kh) -(1 -η λ ) γα 2 + N 2 kh sh(kh) 4N 4 η λ (1 -ch(kh)) + Rc β k 2 +2N 2 η λ sh(kh) + (1 -η λ )(1 -ch(kh)) N 2 kh γ α kh + 2N
) = 2N 2 k ch(ky 3 ) -µ λ ch(kh) -(1 -µ λ ) 2N 2 k 1 -(1 -ch(kh)) sh(ky 3 ) sh(kh) B + 1 2(1-N 2 ) y 2 3 -µ λ h 2 -(1 -µ λ ) h 2 1-N 2 sh(ky 3 ) sh(kh) + 1 1-N 2 (y 3 -µ λ h) -(1 -µ λ ) h 1-N 2 sh(ky 3 ) sh(kh) A ∂ y 1 p(y 1 ) , w2 (y 1 , y 3 ) = sh(ky 3 ) + (1 -µ λ )(1 -ch(kh)) ch(ky 3 ) sh(kh) B + y 3 2(1-N 2 ) -(1 -µ λ ) kh 2 2N 2 (1-N 2 ) ch(ky 3 ) sh(kh) + 1 2(1-N 2 ) -(1 -µ λ ) kh 2N 2 (1-N 2 )
ch(ky 3 ) sh(kh) A ∂ y 1 p(y 1 ) , (6.14) where

µ λ = 1 - N 2 1 -N 2 sh(kh) k E λ -1 , (6.15)
A = L 4N 4 µ λ (1 -ch(kh)) + Rc β k 2 h -(1 -µ λ )coth(kh) h 2 k 2N 2 -k sh(kh) + (1 -µ λ )coth(kh)(1 -ch(kh)) Rc β -2N 2 h 2 µ λ B =k L 2(1-N 2 ) 2N 2 h 2 µ λ + Rc β -(1 -µ λ )coth(kh) hk N 2 Rc β , L = -1 -(1 -µ λ )coth(kh) hk N 2 4N 4 µ λ (1 -ch(kh)) + Rc β k 2 +4N 2 hkµ λ sh(kh) + (1 -µ λ )coth(kh)(1 -ch(kh)) -1
.

Let us define

C α = αh α -1 , C N = N 2 1 -N 2 sh(kh) k , E = Q |D z φ 1,1 | 2 dz.
as consequence of the previous results, we get the following developments in powers of λ, which will be useful for the numerical part.

Corollary 2. We obtain the following developments:

ũ1 (y 1 , y 3 ) = v 0 (y 3 )∂ y 1 p 0 (y 1 ) + C j Eλ 2 v 1 (y 3 )∂ y 1 p 0 (y 1 ) + v 0 (y 3 )∂ y 1 p 1 (y 1 ) + O(λ 4 ), w2 (y 1 , y 3 ) = 0 (y 3 )∂ y 1 p 0 (y 1 ) + C j Eλ 2 1 (y 3 )∂ y 1 p 0 (y 1 ) + 0 (y 3 )∂ y 1 p 1 (y 1 ) + O(λ 4 ), p(y 1 ) = p 0 (y 1 ) + C j Eλ 2 p 1 (y 1 ) + O(λ 4 ), Θ λ = Θ 0 -C j Eλ 2 Θ 1 + O(λ 4 ), (6.16)
where j = α if α = 1 and j = N if α = 1. Here, p 0 satisfies the following equation

1 0 Θ 0 ∂ y 1 p 0 (y 1 )∂ y 1 θ(y 1 ) dy 1 = 1 0
Sθ(y 1 ) dy 1 , ∀θ ∈ H 1 (0, 1), (6.17)

and p 1 the following one

1 0 Θ 0 ∂ y 1 p 1 (y 1 )∂ y 1 θ(y 1 ) dy 1 = 1 0 Θ 1 ∂ y 1 p 0 (y 1 )∂ y 1 θ(y 1 ) dy 1 , ∀θ ∈ H 1 (0, 1). ( 6

.18)

In particular, p 1 is given by p 1 (y 1 ) = Θ 1 Θ 0 p 0 (y 1 ).

For α = 1, v i , i , Θ i , i = 0, 1, are defined as follows

v 0 (y 3 ) = 1 2(1-N 2 ) -2N 2 k (sh(ky 3 ) -sh(kh)) + γ α (y 3 -h) L 0 A 0 -2N 2 k (ch(ky 3 ) -ch(kh)) L 0 B 0 + y 2 3 -h 2 , v 1 (y 3 ) = 1 2(1-N 2 ) -2N 2 kh sh(kh) + γ α (h -1)L 0 A 0 -2N 2 kh (1 -ch(kh)) y 3 -h L 0 B 0 -2N 2 k (sh(ky 3 ) -sh(kh)) + γ α (y 3 -h) L 0 (A 1 + L 0 L 1 A 0 ) -2N 2 k ch(ky 3 ) -ch(kh) L 0 (B 1 + L 0 L 1 B 0 ) + y 3 h + h 2 , 0 (y 3 ) = 1 2(1-N 2 ) -ch(ky 3 ) + γ α L 0 A 0 -sh(ky 3 )L 0 B 0 + y 3 , 1 (y 3 ) = 1 2(1-N 2 ) -ch(ky 3 ) + γ α L 0 (A 1 + L 0 L 1 A 0 ) + γα 2 + N 2 kh sh(kh) L 0 A 0 -sh(ky 3 )L 0 (B 1 + L 0 L 1 B 0 ) -(1 -ch(kh)) N 2 kh L 0 B 0 + h 2 Θ 0 = h 3 3(1-N 2 ) -1 2(1-N 2 ) 2N 2 k ch(kh)-1 k -hsh(kh) -γα 2 h 2 L 0 A 0 + 2N 2 k sh(kh) k
hch(kh) L 0 B 0 , For α = 1, v i , i , Θ i , i = 0, 1, are defined as follows v 0 (y 3 ) = -N 2 (1-N 2 ) ch(ky 3 ) -ch(kh) L 0 B 0 + 1 2(1-N 2 ) y 2 3 -h 2 -1 1-N 2 (y 3 -h)L 0 A 0 , v 1 (y 3 ) = -N 2 1-N 2 ch(ky 3 ) -ch(kh) L 0 (B 1 -L 0 L 1 B 0 ) -N 2 1-N 2 (1 -ch(kh)) sh(ky 3 ) sh(kh) L 0 B 0 + h 2 1-N 2 -1 2 + sh(ky 3 ) sh(kh) -1 1-N 2 (y 3 -h)L 0 (A 1 -L 0 L 1 A 0 ) -L 0 A 0 -1 + sh(ky 3 ) sh(kh) , 0 (y 3 ) = -k 2(1-N 2 ) sh(ky 3 )L 0 B 0 + y 3 2(1-N 2 ) -

Θ 1 = 3h 3 4(1-N 2 ) -1 2(1-N 2 ) 2N 2 k ch(kh)-1 k -hsh(kh) -γα 2 h 2 L 0 (A 1 + L 0 L 1 A 0 ) + -2N 2 k sh(kh)(1 -h) -γ α h(1 -h) L 0 A 0 + 2N 2
1 2(1-N 2 ) L 0 A 0 , 1 (y 3 ) = -k 2(1-N 2 ) sh(ky 3 )L 0 (B 1 -L 0 L 1 B 0 ) + k 2(1-N 2 ) (1 -ch(kh)) ch(ky 3 ) sh(kh) L 0 B 0 kh 2 2N 2 (1-N 2 ) ch(ky 3 ) sh(kh) - 1 2(1-N 2 ) L 0 (A 1 -L 0 L 1 A 0 ) - kh 2N 2 (1-N 2 )
ch(ky 3 ) sh(kh) L 0 A 0 ,

Θ 0 = h 3 3(1-N 2 ) -h 2 2(1-N 2 ) L 0 A 0 + N 2 1-N 2 sh(kh) k -h ch(kh) L 0 B 0 , Θ 1 = -h 2 1-N 2 h 2 + 1-ch(kh)
k sh(kh)

+ h 2 2(1-N 2 ) L 0 A 1 -L 0 L 1 A 0 + h 1-N 2 h + 1-ch(kh) k sh(kh) L 0 A 0 -N 2 1-N 2 sh(kh) k -h ch(kh) L 0 B 1 -L 0 L 1 B 0 -N 2 1-N 2 -h ch(kh) + h + (1-ch(kh)) 2
k sh(kh)

L 0 B 0 .

where A i , B i , L i , i = 0, 1 are given by with γ α defined in (6.13). Using the condition u 1 = 0 on Γ 1 in (6.25), (6.27) and taking into account η λ defined in (6.13), we obtain K 2 (y 1 ) = -η λ 2N 2 k sh(kh) + hγ α Ã(y 1 )

A 0 = h 4N 4 (1 -ch(kh)) + Rc β k 2 -k sh(kh) Rc β -2N 2 h 2 , A 1 = h 4N 4
-η λ 2N 2 k ch(kh) + E λ αh α -1 2N 2 k B(y 1 ) -η λ h 2 2(1 -N 2 )
∂ y 1 p(y 1 ).

Plugging this expression in (6.27) yields

K 1 (y 1 ) = (1 -N 2 ) γ α -E λ α α -1 η λ 2N 2 k sh(kh) + hγ α Ã(y 1 ) +E λ α(1 -N 2 ) α -1 2N 2 k 1 -η λ ch(kh) + E λ αh α -1 B(y 1 ) -E λ αh 2 2(α -1)
η λ ∂ y 1 p(y 1 ).

Taking into account that

η λ = 1 + αh α -1 E λ -1 ⇔ η -1 λ -1 = αh α -1 E λ ⇔ 1 -η λ = αh α -1 E λ η λ ,
we rewrite K 1 and K 2 as follows From condition w 2 = 0 on Γ 1 and (4.21) 2 , we obtain, using (6.28) and (6.29), the following system

K 1 (y 1 ) = (1 -N 2 ) γ α - 1 -η λ h
Q Ã B =      - h 4(1 -N 2 ) (1 + η λ ) Rc β -2N 2 h 2 η λ 2(1 -N 2 )      ∂ y 1 p(y 1 )
where Q is the matrix defined by

Q =     γα 2 + ch(kh) -1-η λ 2h hγ α + 2N 2 k sh(kh) sh(kh) + 1-η λ h (1 -ch(kh)) N 2 k 2N 2 η λ γ α h + 2N 2 k sh(kh) -4N 4 k η λ (1 -ch(kh)) -Rc β k     .
(6.30)

The, the solution of this system is given by Ã(y 1 ) = A∂ y 1 p(y 1 ), B(y 1 ) = B∂ y 1 p(y 1 ),

54

where A and B are solution of

Q A B = 1 2(1 -N 2 )     - h 2 (1 + η λ ) R c β -2N 2 h 2 η λ     ,
Computing A, B, then u 1 and w 2 are obtained by (6.25) and (6.26) as functions of p and of known data.

Proof of Lemma 5. The beginning of the proof is as in Lemma 4. Using (6.24), (6.25) and (6.26) in the boundary condition (4.21) 1 , for α = 1, we obtain .31) with k given in (6.13). Using the condition u 1 = 0 on Γ 1 in (6.25), taking (6.31) into account, we obtain

Ã(y 1 ) = - E λ 2(1 -N 2 ) K 2 (y 1 ) - N 2 1 -N 2 E λ k B(y 1 ) . ( 6
K 2 (y 1 ) = -µ λ h 1 -N 2 K 1 (y 1 ) + 2N 2 k µ λ (1 -ch(kh)) -1 B(y 1 ) -µ λ h 2 2(1 -N 2 )
∂ y 1 p(y 1 ) , (

with µ λ defined in (6.15). Using the definition of µ λ , we rewrite A as follows From the conditions w 2 = 0 on Γ 1 and (4.21) 2 , using (6.32) and (6.33) the following system is obtained

Q K 1 B =      1 2(1 -N 2 ) -h + (1 -µ λ )coth(kh) h 2 k 2N 2 Rc β -2N 2 h 2 µ λ 2(1 -N 2 )      ∂ y 1 p(y 1 ) ,
where Q is the matrix defined by

Q =   1 2(1-N 2 ) 1 -(1 -µ λ )coth(kh) kh N 2 sh(kh) + (1 -µ λ )coth(kh)(1 -ch(kh)) 2N 2 h 1-N 2 µ λ -4N 4 k µ λ (1 -ch(kh)) -Rc β k   .
The solution of this system is given by K 1 (y 1 ) = A ∂ y 1 p(y 1 ), B(y 1 ) = B ∂ y 1 p(y 1 ) ,

55

where A and B are solution of

Q A B =      1 2(1 -N 2 ) -h + (1 -µ λ )coth(kh) h 2 k 2N 2 Rc β -2N 2 h 2 µ λ 2(1 -N 2 )      .
Computing A , B , then u 1 and w 2 are obtained as functions of p and of known data.

Proof of Corollary 2. We first remark that the roughness parameter E λ given in Theorem 2 satisfies

E λ = Q |D z φ 1,λ | dz = λ 2 Q |D z φ 1,1 | dz = λ 2 E.
We explain the case α = 1 (for the case α = 1 proceed similarly). Using power series of λ 2 and omitting terms of order O(λ 4 ), there holds

η λ = (1 + C α Eλ 2 ) -1 ∼ 1 -C α Eλ 2 , with C α = αh α -1 , η 2 λ = (1 + C α Eλ 2 ) -2 ∼ 1 -2C α Eλ 2 .
Using the development of η λ in terms of λ 2 in A, B and L given in Lemma 4, we get

A = -1 2(1-N 2 ) L 0 A 0 + C α Eλ 2 [A 1 + L 0 L 1 A 0 ] , B = -1 2(1-N 2 ) L 0 B 0 + C α Eλ 2 [B 1 + L 0 L 1 B 0 ] , L = -L 0 1 + L 0 L 1 C α Eλ 2 ,
with A i , B i and L i , i = 0, 1 given by (6.19). Next, developing the pressure as p(y 1 ) = p 0 (y 1 ) + C α Eλ 2 p 1 (y 1 ) + O(λ 4 ) and using previous development of A, B and L in the expressions of ũ1 and w2 given in (6.12), we get the expressions (6.16).

Finally, using again the development of A, B and L given above, the development of p and the development of ũ1 in the Reynolds equation (4.26), we deduce that p 0 and p 1 satisfy (6.17) and (6.18), respectively. Combining the equations on p 0 and p 1 , one easily sees that Θ 0 Θ 1 p 1 and p 0 satisfy the same equation, hence they are equal by uniqueness of the solution to (6.17).

Proof of Lemma 6. The beginning of the proof is as in Lemma 4. In this case, we consider the boundary conditions given in (4.18). Thus, using (6.25) and boundary conditions ũ1 (y 1 , 0) = 0 and ũ1 (y 1 , h) = 0, respectively, we have with k given in (6.13). From the boundary conditions w2 (y 1 , h) = 0 and ∂ y 3 w2 (y 1 , 0) = 0, the following system is obtained

K 2 (y 1 ) = -
Q Ã B =     -h 4(1 -N 2 ) - 1 2(1 -N 2 )k     ∂ y 1 p(y 1 ),

  18)-(2.25) (in the sense of Definition 1).

∂ y 3 w 3 ∂ 0 - 0 - 2 =

 3002 y 3 ψ dy + 4N 2 Ω w 3 ψ dy = 0.(5.27)Next, we prove that w 3 (y , 0) = 0 for a.e. x ∈ ω. The condition w ε • n ε = 0 on Γ 0 ε can be rewritten as followsε wε,1 (y , -Ψε (y ))λε δ-∂ 1 Ψ( 1 ε y • e 1 ) + ε wε,3 (y , -Ψε (y )) = 0 a.e. y ∈ ω.Multiplying this equality by ψ ∈ C ∞ c (ω) and integrating on ω, we getω ε wε,1 (y , -Ψε (y ))λe δ-∂ 1 Ψ( 1 ε y • e 1 )ψ(y )dy + ω ε wε,3 (y , -Ψε (y ))ψ(y )dy = 0. (5.28) We can write the second term of (5.28) as follows ω ε wε,3 (y , -Ψε (y ))ψ(y )dy = ω ε wε,3 (y , 0)ψ(y )dyω Ψε(y ) ε∂ y 3 wε,3 (y , s)ds ψ(y )dy . Then, since ε wε | Ω is bounded in H 1 (0, h; L 2 (ω)) 3 , by continuity of the trace operator from H 1 (0, h; L 2 (ω)) 3 into L 2 (ω × {0}) 3 and convergence (5.25), we have that ω ε wε,3 (y , 0)dy = ω ε w3 (y , 0)dy + O ε . Moreover, from the Cauchy-Schwarz inequality, estimate (5.19) 2 and | Ω - ε | → 0 (recall that Ψ ε L ∞ ≤ Cε δ-1 ), then ω Ψε(y ) ε∂ y 3 (y , s)dsϕ(y )dy ≤ Ωε |ε∂ y 3 wε,3 | 2 O ε , and so, we get that ω ε wε,3 (y , -Ψε (y ))dy = ω ε w3 (y , 0)dy + O ε .

Figure 2 :

 2 Figure 2: Example of function Θ λ and its asymptotic development Θ0 -CjEλ 2 Θ1 plotted against λ, for the set of parameters N = 0.3, Rc = 0.1, δ = 1, E = 10, and with νb = 0.1 (case α = 1, left) and νb = 1 (case α = 1, right).

4 R 2 Figure 3 :

 423 Figure 3: T half plotted against N , for δ = 1, E = 0 and different values of Rc ∈ {0.025, 0.05, 0.1, 0.2}, with νb ∈ {0.05, 0.1, 0.2, 0.4}.

4 R 2 Figure 4 : 10 R

 42410 Figure 4: T half plotted against N , for δ = 1, E = 10 and different values of Rc ∈ {0.025, 0.05, 0.1, 0.2}, with νb ∈ {0.05, 0.1, 0.2, 0.4}.

Figure 5 : 10 R

 510 Figure 5: T half plotted against N , for νb = 0.1, Rc ∈ {0.025, 0.05, 0.1, 0.2} and δ = 1. From left to right and top to bottom: E = 0, E = 1, E = 3, E = 5, E = 7 and E = 10.

Figure 6 :Figure 7 :

 67 Figure 6: T half plotted against N , for νb = 0.1, Rc ∈ {0.025, 0.05, 0.1, 0.2} and E = 5. From left to right and top to bottom: δ = 0.7, δ = 0.8, δ = 1,δ = 1.2, δ = 2 and δ = 10.

1 - 8 • 10 - 2 - 6 • 10 - 2 - 4 • 10 - 2 - 2 Figure 8 :

 181026102410228 Figure 8: Examples of rough geometries associated with functions Ψε(x1) defined by (6.3) and (6.4), with L = 2, M = 0.5 and Λ = M 3/2 . Each riblet profile is described by the corresponding function Ψ plotted in Fig. 7.

  R c = 0.05 R c = 0.1 R c = 0.2

Figure 9 : 2 Figure 10 :

 9210 Figure 9: CαΘ1/Θ0 (top, left) and Θ1/Θ0 (top, right) plotted against N , normalized by the value at N = 0, for νb = 0.1, δ = 1, E = 5 and different values of Rc ∈ {0.025, 0.05, 0.1, 0.2}. Bottom: T half plotted against N for the same sets of parameters.

Lemma 4 .

 4 For α = 1, the solutions of system (4.19) with boundary conditions (4.20)-(4.21) are:

- 2 2- 1 ,L 1 = 4N 4 ( 1 -

 21141 hch(kh) L 0 (B 1 + L 0 L 1 B 0 ) + 2N 2 k hch(kh) -(1 + ch(kh))h N 2 k L 0 B 0 ,where A i , B i , L i , i = 0, 1 are given byA 0 = h 4N 4 (1 -ch(kh)) + Rc β k 2 -ksh(kh) Rc β -2N 2 h 2 , A 1 = -h 4N 4 (1 -ch(kh)) + Rc β k 2N 2 h 2 ksh(kh) -(1 -ch(kh)) N 2 h Rc β , B 0 = 2N 2 h γ α kh + 2N 2 sh(kh) + k Rc β -2N 2 h 2 ch(kh) + γα 2 , B 1 = -γα 2 k 2N 2 h 2 + Rc β -sh(kh) 4N 4 h + Rc β N 2 h + 2N 2 h 2 kch(kh), L 0 = γα 2 + ch(kh) 4N 4 (1 -ch(kh)) + Rc β k 2 + 2N 2 sh(kh) γ α kh + 2N 2 sh(kh) ch(kh)) γ α + ch(kh) + N 2 kh sh(kh) + Rc β k 2 γα 2 + N 2 kh sh(kh) +2N 2 sh(kh) -N 2kh (1 -ch(kh)) γ α kh + 2N 2 sh(kh) .

( 1 - 4 ( 1 - 1 ,L 1 = 2 Lemma 6 . 2 kch(ky 3 ) 1 2( 1 -N 2 ) y 2 3 - 1 2( 1 -N 2 ) y 3 -h 2 ∂

 1411126231131122 ch(kh)) + 2N 2 hksh(kh) + Rc β k coth(kh) 1 -ch(kh) + k 2 h 2 2N 2 , B 0 = 2N 2 h 2 + Rc β , B 1 = 2N 2 h 2 + coth(kh) ch(kh)) + Rc β k 2 + 4N 2 hksh(kh) -4N 4 (1 -ch(kh)) + coth(kh) k 3 h N Rc β + 4khN 2 sh(kh). In the super-critical case 1 < δ < 3 2 -12 , the solutions of system (4.19) with boundary conditionsu 1 = w 2 = 0 on Γ 1 , u 1 = ∂ y 3 w 2 = 0 on Γ , are ũ1 (y 1 , y 3 ) = 2N 2 k sh(ky 3 ) -y 3 h sh(kh) A + 2N -y 3 h (ch(kh) -1) -1 B + y 3 h ∂ y 1 p(y 1 ) , w2 (y 1 , y 3 ) = ch(ky 3 ) -N 2 kh sh(kh) A+ sh(ky 3 ) -N 2 kh (ch(kh) -1) B + y 1 p(y 1 ) ,(6.20)

( 1 - 2 ∂ y 1 p(y 1 K 2 1 ) -η λ h 2 2 ( 1 -

 12112121 N 2 )(1 -ch(kh))(1 -η λ ) B(y 1 ) -(1 -η λ ) h (y 1 ) = -η λ 2N 2k sh(kh) + hγ α Ã(y 1 ) -2N 2 k (η λ (ch(kh) -1) + 1) B(y N 2 ) ∂ y 1 p(y 1 ).(6.29)

Ã(y 1 )

 1 = -(1 -µ λ ) 1 2N 2 (1 -N 2 ) kh sh(kh) K 1 (y 1 ) + (1 -µ λ ) 1 -ch(kh) sh(kh) B(y 1 )

2N 2 k

 2 B(y 1 ) , K 1 (y 1 ) = -2N 2 (1 -N 2 ) kh sh(kh) Ã(y 1 ) -2N 2 (1 -N 2 ) kh (ch(kh) -1) B(y 1 ) -h 2 ∂ y 1 p(y 1 ) ,

  endowed with the norm Dϕ L 2 (Ωε) 3 . Assume that (u ε , w ε , p ε ) is a classical solution to system (2.18)-(2.25). Multiplying (2.18) by a test function ϕ ∈ V ε , integrating by parts and taking into account the boundary conditions and the free divergence condition satisfied by ϕ, we obtain

Ωε

Du ε : Dϕ dx -

conditions (2.10) must be completed by tangential conditions on [Du ε n ε ] tan and [Dw ε n ε ] tan , in the aim of obtaining a well-posed problem.To obtain conditions (2.11)-(2.12), we have interpreted the rotational terms (∇ × u) × n, (∇ × w) × n appearing in the initial formulatation of the tangential boundary conditions (1.3)-(1.4), as being respectively equal to [Du n] tan and [Dw n] tan . This is indeed the case for a flat boundary Γ of normal n, since for regular vector fields u, v satisfying u• n = v • n = 0 on Γ, there holds Γ [(∇ × u) × n] • v dσ = Γ [Du n] tan • v dσ .Last equality is obtained by writing(∇ × u) × n = Du n -(Du) T n and using that v • ∇(u • n) = 0 on Γ, which gives Γ [(∇ × u) × n] • v dσ = Γ [Du n] • v -Γ [(Du) T n] • v dσ = Γ [Du n] tan • v dσ -Γ [(v • ∇)u] • n dσ = Γ [Du n] tan • v dσ + Γ [(v • ∇)n] • u dσ = Γ [Du n] tan • v dσbecause (v • ∇)n = 0 since n is a constant vector.

Acknowledgements

The first author has been supported by a public grant as part of the Investissement d'avenir project, reference ANR-11-LABX-0056-LMH, LabEx LMH.

The second author has been supported by the Croatian Science Foundation under the project Multiscale problems in fluid mechanics -MultiFM (IP-2019-04-1140).

The third author has been partially supported by the Ministerio de Economía y Competitividad (Spain), under the project Proyecto Excelencia MTM2014-53309-P.

The authors would also like to thank the referees for their helpful remarks and suggestions.

Appendix

This section completes Section 4.3 by providing the explicit formulas for the solution (ũ 1 , w2 ) to system (4.19), in the critical case (Lemmas 4 and 5) and in the super-critical case (Lemma 6). In the critical case, we also give the first terms of the developments of ũ1 , w2 , p and the coefficient Θ λ appearing in Reynolds equation (4.26) in powers of λ 2 (see Corollary 2). We conclude the Appendix by the proofs of these statements.

where

.

In the rest of the appendix, we give the proofs of Lemmas

where K 1 is an unknown function. Putting (6.21) into (4.19) 2 , we obtain

whose solution can be written as

where c i (y 1 ), i = 1, 2 are unknown functions and k = 2N 1 -N 2 R c .

Putting (6.23) in (6.21), we obtain

(6.24)

Integrating (6.24), we get

where Ã(y 1 ) = c 1 (y 1 ) + c 2 (y 1 ), B(y 1 ) = c 1 (y 1 ) -c 2 (y 1 ), and K 2 (y 1 ) is an unknown function. w 2 in (6.23) can also be written as follows:

w 2 (y 1 , y 3 ) = Ã(y 1 )ch(ky 3 ) + B(y 1 )sh(ky 3 ) +

. (6.26)

Using (6.24) and (6.26) in the boundary condition (4.21) 1 , we obtain

where Q is the matrix defined by

.

The solution of this system is given by Ã(y 1 ) = A ∂ y 1 p(y 1 ), B(y 1 ) = B ∂ y 1 p(y 1 ) ,

where A and B are solution of

Computing A , B , then u 1 and w 2 are obtained as functions of p and of known data.