%~Mouliné par MaN_auto v.0.29.1 2023-11-13 19:21:00
\documentclass[CRMATH,Unicode,XML]{cedram}

\TopicFR{\'Equations aux dérivées partielles, Mécanique}
\TopicEN{Partial differential equations, Mechanics}

\usepackage[nameinlink]{cleveref}
\usepackage{mathrsfs}
\usepackage{mathtools}
\usepackage{empheq}

\crefname{prop}{Proposition}{Propositions}
\crefname{propc}{Proposition}{Propositions}
\crefname{theo}{Theorem}{Theorems}
\crefname{defi}{Definition}{Definitions}
\crefname{rema}{Remark}{Remarks}
\crefname{remac}{Remark}{Remarks}
\crefname{lemm}{Lemma}{Lemmas}
\crefname{coro}{Corollary}{Corollaries}
\crefname{nota}{Notation}{Notations}
\crefname{clai}{Claim}{Claims}
\crefname{exem}{Example}{Example}
\crefname{exam}{Example}{Example}
\crefname{conj}{Conjecture}{Conjectures}
\crefname{section}{Section}{Sections}
\crefname{subsection}{Section}{Sections}
\crefname{appendix}{Appendix}{Appendices}
\crefname{cdrenonce}{Problem}{Problems}
\crefname{figure}{Figure}{Figures}

\newcommand*{\dd}{\mathrm{d}}
\newcommand*{\dx}{\dd x}
\newcommand*{\dt}{\dd t}

\DeclareMathOperator{\Reop}{Re}
\DeclareMathOperator{\We}{We}
\DeclareMathOperator{\Fr}{Fr}

\newenvironment{smallpmatrix}{\left(\begin{smallmatrix}}{\end{smallmatrix}\right)}

\makeatletter
\def\@setafterauthor{%
  \vglue3mm%
%  \hspace*{0pt}%
\begingroup\hsize=12cm\advance\hsize\abstractmarginL\raggedright
\noindent
%\hspace*{\abstractmarginL}\begin{minipage}[t]{10cm}
   \leftskip\abstractmarginL
  \normalfont\Small
  \@afterauthor\par
\endgroup
\vskip2pt plus 3pt minus 1pt
}
\makeatother

\DeclareMathOperator{\sym}{sym}
%\DeclarePairedDelimiter\paren\lparen\rparen
\newcommand*\paren[1]{\left(#1\right)}


\DeclarePairedDelimiter\abs\lvert\rvert
\DeclarePairedDelimiter\norm\lVert\rVert
\DeclarePairedDelimiter\braket\langle\rangle
\DeclarePairedDelimiter\brackets{\lbrack}{\rbrack}
\DeclarePairedDelimiter\braces\{\}


\newcommand{\spdim}{N}
\newcommand{\red}[1] {{\color[rgb]{1,0,0}{{#1}}}}

\newcommand{\half}[1][1]{\dfrac{#1}{2}}
\newcommand{\height}{h}
\newcommand{\velocity}[1][\boldsymbol]{#1{v}}
\newcommand{\stress}[1][\boldsymbol]{#1{\sigma}}
\newcommand{\substrate}{\Sigma_s}
\newcommand{\freesurface}{\Sigma_f}
\newcommand{\body}{\mathcal{B}}
\newcommand{\R}{\mathbb{R}}
\newcommand{\mean}[2][\height]{\braket{#2}_{#1}}
\newcommand{\bigmean}[2][\height]{\braket*{#2}_{#1}}
\newcommand{\var}{\mathrm{var}_{\height}}
\newcommand{\cov}{\mathrm{cov}_{\height}}
\newcommand{\Div}{\mathrm{div}}
\newcommand{\inv}[1]{\dfrac{1}{#1}}
\newcommand{\pdd}[3][]{\ensuremath{\dfrac{\partial\ifx\relax#1\relax\else^{#1}\fi#2}{\partial #3\ifx\relax#1\relax\else^{#1}\fi}}}
\newcommand{\textpdd}[3][]{\ensuremath{\partial\ifx\relax#1\relax\else^{#1}\fi#2 / \partial #3\ifx\relax#1\relax\else^{#1}\fi}}
\newcommand{\pddt}[2][]{\ensuremath{\pdd[#1]{#2}{t}}}
\newcommand{\textpddt}[2][]{\ensuremath{\textpdd[#1]{#2}{t}}}
\newcommand{\advection}[2][]{\ensuremath{\paren{#2 \cdot \nabla\ifx\relax#1\relax\else_{#1}\fi}}}
\newcommand{\unknown}{u}
\newcommand{\flux}[1][\boldsymbol]{#1{j}}

\newcommand*\Dscr{\mathscr{D}}
\newcommand*\ubold{\boldsymbol u}

\NewDocumentCommand{\materialdd}{o o m}{
\IfNoValueTF{#2} {
\IfNoValueTF {#1} {
\ensuremath{\dfrac{D #3}{D t}} } {
\ensuremath{\dfrac{D_{#1} #3}{D t}} } } {
\ensuremath{\dfrac{D_{#1}^{#2} #3}{D t}} } }
\NewDocumentCommand{\textmaterialdd}{o o m}{
\IfNoValueTF{#2} {
\IfNoValueTF{#1} {
\ensuremath{D #3 / D t} } {
\ensuremath{D_{#1} #3 / D t} } } {
\ensuremath{D_{#1}^{#2} #3 / D t} } }

\NewDocumentCommand{\gsdd}{o o m}{
\IfNoValueTF{#2} {
\IfNoValueTF {#1} {
\ensuremath{\dfrac{\Dscr #3}{\Dscr t}} } {
\ensuremath{\dfrac{\Dscr_{#1} #3}{\Dscr t}} } } {
\ensuremath{\dfrac{\Dscr_{#1}^{#2} #3}{\Dscr t}} } }
\NewDocumentCommand{\textgsdd}{o o m}{
\IfNoValueTF{#2} {
\IfNoValueTF{#1} {
\ensuremath{\Dscr #3 / \Dscr t} } {
\ensuremath{\Dscr_{#1} #3 / \Dscr t} } } {
\ensuremath{\Dscr_{#1}^{#2} #3 / \Dscr t} } }


\newcommand\MTkillspecial[1]{
\bgroup
\catcode`\&=9
\let\\\relax
\scantokens{#1}
\egroup }


%\reDeclarePairedDelimiterInnerWrapper\paren{star}{
%\mathopen{#1\vphantom{\MTkillspecial{#2}}\kern-\nulldelimiterspace\right.} #2
%\mathclose{\left.\kern-\nulldelimiterspace\vphantom{\MTkillspecial{#2}}#3}}


\reDeclarePairedDelimiterInnerWrapper\abs{star}{
\mathopen{#1\vphantom{\MTkillspecial{#2}}\kern-\nulldelimiterspace\right.} #2
\mathclose{\left.\kern-\nulldelimiterspace\vphantom{\MTkillspecial{#2}}#3}}


\reDeclarePairedDelimiterInnerWrapper\norm{star}{
\mathopen{#1\vphantom{\MTkillspecial{#2}}\kern-\nulldelimiterspace\right.} #2
\mathclose{\left.\kern-\nulldelimiterspace\vphantom{\MTkillspecial{#2}}#3}}


\reDeclarePairedDelimiterInnerWrapper\braket{star}{
\mathopen{#1\vphantom{\MTkillspecial{#2}}\kern-\nulldelimiterspace\right.} #2
\mathclose{\left.\kern-\nulldelimiterspace\vphantom{\MTkillspecial{#2}}#3}}


\reDeclarePairedDelimiterInnerWrapper\brackets{star}{
\mathopen{#1\vphantom{\MTkillspecial{#2}}\kern-\nulldelimiterspace\right.} #2
\mathclose{\left.\kern-\nulldelimiterspace\vphantom{\MTkillspecial{#2}}#3}}


\reDeclarePairedDelimiterInnerWrapper\braces{star}{
\mathopen{#1\vphantom{\MTkillspecial{#2}}\kern-\nulldelimiterspace\right.} #2
\mathclose{\left.\kern-\nulldelimiterspace\vphantom{\MTkillspecial{#2}}#3}}

\DeclarePairedDelimiterX\innerprod[2]{\langle}{\rangle}{#1,#2}
\DeclarePairedDelimiterX\altinnerprod[2]{\lparen}{\rparen}{#1,#2}

\DeclarePairedDelimiterX{\dbrackets}[1]{\llbracket}{\rrbracket}{
\nhphantom{\lbrack}{#1} \delimsize\lbrack \mathopen{} #1 \mathclose{} \delimsize\rbrack \nhphantom{\rbrack}{#1} }
\DeclarePairedDelimiterX{\dbraces}[1]{\lbrace}{\rbrace}{
\nhphantom{\lbrace}{#1} \delimsize\lbrace \mathopen{} #1 \mathclose{} \delimsize\rbrace \nhphantom{\rbrace}{#1} }

\graphicspath{{./figures/}}

\newcommand*{\mk}{\mkern -1mu}
\newcommand*{\Mk}{\mkern -2mu}
\newcommand*{\mK}{\mkern 1mu}
\newcommand*{\MK}{\mkern 2mu}

\hypersetup{urlcolor=purple, linkcolor=blue, citecolor=red}

\newcommand*{\relabel}{\renewcommand{\labelenumi}{(\theenumi)}}
\newcommand*{\romanenumi}{\renewcommand*{\theenumi}{\roman{enumi}}\relabel}
\newcommand*{\Romanenumi}{\renewcommand*{\theenumi}{\Roman{enumi}}\relabel}
\newcommand*{\alphenumi}{\renewcommand*{\theenumi}{\alph{enumi}}\relabel}
\newcommand*{\Alphenumi}{\renewcommand*{\theenumi}{\Alph{enumi}}\relabel}
\let\oldtilde\tilde
\renewcommand*{\tilde}[1]{\mathchoice{\widetilde{#1}}{\widetilde{#1}}{\oldtilde{#1}}{\oldtilde{#1}}}
\let\oldhat\hat
\renewcommand*{\hat}[1]{\mathchoice{\widehat{#1}}{\widehat{#1}}{\oldhat{#1}}{\oldhat{#1}}}
\let\oldforall\forall
\renewcommand*{\forall}{\mathrel{\oldforall}}


\title{A new framework for shallow approximations of incompressible flows}

\author{\firstname{Nathan} \lastname{Shourick}\IsCorresp}
\address{Lab. Jean Kuntzmann -- CNRS and Université Grenoble-Alpes, F-38041 Grenoble, France}
\address{Univ. Grenoble Alpes, CNRS, UMR 5525, VetAgro Sup, Grenoble INP, TIMC, 38000 Grenoble, France}
\email{nathan.shourick@univ-grenoble-alpes.fr}

\author{\firstname{Ibrahim} \lastname{Cheddadi}}
\address[2]{Univ. Grenoble Alpes, CNRS, UMR 5525, VetAgro Sup, Grenoble INP, TIMC, 38000 Grenoble, France}
\email{ibrahim.cheddadi@univ-grenoble-alpes.fr}

\author{\firstname{Pierre} \lastname{Saramito}}
\address[1]{Lab. Jean Kuntzmann -- CNRS and Université Grenoble-Alpes, F-38041 Grenoble, France}
\email{pierre.saramito@imag.fr}
\subjclass{35Q35}

\begin{abstract}
A new framework for the asymptotic analysis of incompressible flows of complex non-Newtonian materials is presented in this paper. It allows both to avoid redundant mathematical hypotheses and to dramatically reduce the amount of tedious formal calculations. The starting point of the proposed framework is a generic equation, easily adaptable to most problems of continuum mechanics, for which a thin-layer approximation is developed. We then show how to treat the so-called Gordon--Schowalter derivative, a general objective time derivative involved in non-Newtonian viscoelastic fluids. As a proof of concept of our framework, we apply it to the Maxwell viscoelastic model.
\end{abstract}

\begin{document}
\maketitle

\section{Introduction}

The scientific literature addresses a wide range of thin-layer flow problems, suitable for geophysical, biological and industrial applications. The corresponding approximation leads to reduce a tridimensional time-dependent and free-surface flow problem to a bidimensional one with an explicit description of the elevation of the free-surface. The problem is then much simpler, from both theoretical and numerical points of view.

This approach was first introduced in 1871 by Saint-Venant~\cite{Sai-1871} for the so-called shallow-water equations that approximate the Euler's fluid equations. Since then, significant progress has been made in obtaining reduced equations. Especially, asymptotic analysis allows a rigorous formal development, up to an arbitrary order, of quantities with respect to~$\epsilon$, the small aspect ratio of the thin geometry (see \cref{fig-asymptotic}). Among others, it has allowed to obtain thin-layer approximations of Navier--Stokes model~\cite{GerPer-2000,Mar-2006-sv}, of the Bingham viscoplastic fluid model~\cite{LiuMei-1990,BalCra-1999,BerSarSmu-2014,FerNobVil-2010}, and more recently of the Maxwell viscoelastic model~\cite{BouBoy-2013}. See e.g.~\cite{BouBoy-2016} for a recent review, including slowly varying topography effects. Note that the reduced equations obtained after asymptotic analysis strongly depend upon the boundary conditions imposed on the substrate on which the fluid is moving. In the literature, two main cases are considered: a slip-friction together with a no-penetration condition~\cite{Mar-2006-sv,GerPer-2000,FerNobVil-2010,BouBoy-2013,BouBoy-2016} and a no-slip condition~\cite{LiuMei-1990,BalCra-1999}.

Our aim is to propose a general framework for facilitating the formal derivation of the thin-layer approximation of incompressible fluids, either Newtonian or non-Newtonian, with slip-friction and no-penetration boundary condition on a flat substrate. It allows both to avoid redundant mathematical hypotheses and to dramatically reduce most tedious formal calculations. While we focus on fluid mechanics applications, the present framework could be suitable to a wider class of partial differential equations from mathematical physics. The idea to provide some kind of generalization has already been suggested in the literature. In~2004, Bouchut and Westdickenberg~\cite{BouWes-2004} studied general varying topographies and, in~2016, Bouchut and Boyaval~\cite{BouBoy-2016} proposed for the first time a unified framework for deriving thin-layer models for shallow free-surface flows driven by gravity, with slip-friction and no-penetration boundary condition, and under the motion by slices hypothesis. The present framework requires a minimal number of assumptions. In particular, we neither assume that the system is driven by gravity nor that the motion by slices hypothesis is satisfied.

The viscoelastic Maxwell model is studied as a proof of concept for the proposed framework. It writes~\cite{Old-1950}:
\begin{subequations}\label{eq-maxwell}
\begin{empheq}[left={\empheqlbrace}]{align}\label{eq-maxwell-evolution}
\lambda\paren{\pddt{\boldsymbol{\tau}} + \advection{\velocity}{\boldsymbol{\tau}} - \nabla\velocity\cdot\boldsymbol{\tau} - \boldsymbol{\tau}\cdot\nabla\velocity^\top} + \boldsymbol{\tau} &= 2\eta D(\velocity), \\
\label{eq-maxwell-momentum}
\rho\paren{\pddt{\velocity} + \advection{\velocity}{\velocity}} - \Div\,\stress &= \rho\boldsymbol g, \\
\Div\,\velocity &= 0,
\end{empheq}
\end{subequations}
where the total stress tensor is $\stress = \boldsymbol{\tau} - p\boldsymbol I$, and the three unknowns are~$\boldsymbol{\tau}$, the elastic stress tensor, $\velocity$, the velocity field and~$p$, the pressure. Here, $\rho$ is the density of the fluid, assumed to be constant, and~$\boldsymbol g$ is the constant gravity vector. The parameter~$\eta$ is a viscosity and~$\lambda$ is a relaxation time. The notation~\mbox{$D(\velocity)$} stands for the rate of deformation tensor, that is the symmetric part of the velocity gradient, defined by~\mbox{$D(\velocity) \coloneqq \left(\nabla\velocity + \nabla\velocity^\top\right)/2$}. Note that when~\mbox{$\lambda=0$}, the Maxwell model reduces to the Navier--Stokes one. Rather than directly performing the asymptotic analysis of the system~\eqref{eq-maxwell}, we propose in this paper to treat the following general equation:
\begin{equation}
\mathcal{A}(\unknown) - \Div\,\boldsymbol j = 0,\label{eq-general}
\end{equation}
where~$\mathcal{A}$ is a given scalar operator depending on the scalar unknown~$\unknown$ and~$\boldsymbol j$ is the equivalent of a flux. This equation is scalar but it can correspond component by component to vectorial conservation equations such as~\eqref{eq-maxwell-momentum}, and tensorial evolution equations such as~\eqref{eq-maxwell-evolution}, and, in general, to a large class of equations from mathematical physics.


After a presentation of the general mathematical setting in~\cref{sec-math-setting}, the asymptotic analysis of~\eqref{eq-general} will be presented along with sufficient conditions for the motion by slices. This section ends with the asymptotic analysis of the general objective Gordon--Schowalter derivatives of both vector and tensor fields. Finally, \cref{sec-examples} turns to applications to Navier--Stokes and Maxwell systems. To the best of our knowledge, the tools presented here are new and we hope that they will \emph{(i)} help to clarify the mathematical hypothesis required for asymptotic analysis, and \emph{(ii)} make formal calculations easier for future applications in the study of complex materials.


\section{Mathematical setting}\label{sec-setting}
\label{sec-math-setting}

\subsection{Geometry}

We consider here a thin-layer fluid domain~\mbox{$\body\subset\R^3$} moving at velocity~\mbox{$\velocity~\colon~\body\times\R_+ \to \R^3$}, with a free surface characterized by its height~\mbox{$\height~\colon~\body\times\R_+ \to \R_+$}, assumed small compared to the flow length. Let~$\substrate$ be the three-dimensional rigid substrate on which the fluid flows. It is defined as the set of points~$\boldsymbol x$ with~\mbox{$x_3 = 0$}, and~$\Omega\subset\R^2$ its two-dimensional counterpart such that~\mbox{$\substrate \coloneqq \Omega \times \{0\}$}, see \cref{fig-asymptotic}. The free surface is defined as the set $\freesurface(t) \coloneqq \Omega \times \{\height(t)\}$. The thin-layer approximation aims at obtaining a two-dimensional problem defined over~$\Omega$.
\begin{figure}[htbp]
\includegraphics[width=\textwidth]{asymptotic}
\caption{Asymptotic analysis of a free-surface flow of a shallow fluid.}\label{fig-asymptotic}
\end{figure}

\subsection{Evolution equations}\label{subsec-setting-equations}

A large part of the equations derived from mechanics can be written as one of the following equations:
\begin{subequations}\label{eq-general-time}
\begin{align}
\lambda\materialdd{\varphi} + b(\varphi) - \Div\,{\flux} &= f, \\
\lambda\gsdd[a]{\boldsymbol\unknown} + \boldsymbol b(\boldsymbol\unknown) - \Div\,{\boldsymbol\sigma} &= \boldsymbol f, \\
\text{or } \quad
\lambda\gsdd[a]{\boldsymbol\tau} + \boldsymbol\beta(\boldsymbol\tau) - \Div\,{\boldsymbol{\mathsf{J}}} &= \boldsymbol\phi,
\end{align}
where~$\lambda$ is a time relaxation, $b$ (resp.~$\boldsymbol b$ and~$\boldsymbol\beta$) is a generic scalar-valued (resp. vector and symmetric tensor) source term, $\flux$ (resp.~$\boldsymbol\sigma$ and~$\boldsymbol{\mathsf{J}}$) is a vector-valued (resp. tensor-valued and order 3 tensor-valued) flux and~$f$ (resp.~$\boldsymbol f$ and~$\boldsymbol\phi$) is a scalar-valued (resp. vector and symmetric tensor) source term independent of the scalar (resp. vector and symmetric tensor) unknown~$\varphi$ (resp.~$\boldsymbol\unknown$ and~$\boldsymbol\tau$). Three time derivatives have also been introduced. The first one is the Lagrangian derivative~\mbox{$\textmaterialdd{\square} \coloneqq \textpddt{\square} + \advection{\velocity}\square$}. The two others are the vector and the tensor Gordon--Schowalter derivatives, parameterized by~\mbox{$a \in \R$} and defined by
\begin{align*}
\gsdd[a]{\ubold} &\coloneqq \materialdd{\boldsymbol\unknown} - \paren{W(\velocity) + aD(\velocity)} \cdot \boldsymbol\unknown, \\
\gsdd[a]{\boldsymbol\tau} &\coloneqq \materialdd{\boldsymbol\tau} - 2\sym\big(\paren{W(\velocity) + aD(\velocity)} \cdot \boldsymbol\tau\big),
\end{align*}
where~$\sym({\,\cdot\,})$ is the tensor symmetric part operator and~\mbox{$W(\velocity) \coloneqq (\nabla\velocity - \nabla\velocity^\top)/2$} is the vorticity tensor, i.e. the skew-symmetric part of the velocity gradient. Note that these time derivatives are defined for a general velocity field~$\velocity$, not necessarily divergence-free.
\end{subequations}

Note that all the three variants of~\eqref{eq-general-time} reduce component-by-component to the general case~\eqref{eq-general} for a suitable choice of the~$\mathcal A$ that includes both the time derivative, the source and the right-hand-side terms. This approach leads to split the difficulties for the asymptotic analysis of the next \cref{sec-framework}: we first study~\eqref{eq-general} with a general operator~$\mathcal A$ before we deal with the complicated terms involved by Gordon--Schowalter derivatives.

\subsection{Boundary conditions}

Only the boundary conditions satisfied on the free surface and on the substrate are involved by the asymptotic analysis. Indeed, conditions on vertical boundary walls have no effects on the reduced models.

\subsubsection{On the free surface}

For fluid flow systems such as Navier--Stokes or Maxwell~\eqref{eq-maxwell}, the boundary conditions on the free surface~$\freesurface(t)$ write~\cite{Mar-2006-sv,BouBoy-2016}:
\begin{subequations}\label{eq-setting-bc-h}
\begin{align}
\stress \cdot \boldsymbol n - \gamma\kappa\boldsymbol n &= \boldsymbol 0,\label{eq-setting-bc-h-surface-tension}
\\
\text{and } \quad
\pddt{\height} + \velocity[]_1\pdd{\height}{x_1} + \velocity[]_2\pdd{\height}{x_2} &= \velocity[]_3\label{eq-setting-bc-h-free-surface}
\end{align}
\end{subequations}
where~$\boldsymbol n$ is the unit outward normal vector, along the direction~\mbox{$ (-\textpdd{\height}{x_1}, -\textpdd{\height}{x_2}, 1)$}. Eqn.~\eqref{eq-setting-bc-h-surface-tension} expresses the normal-stress continuity condition where~$\gamma$ is the surface tension at the air/fluid interface and~\mbox{$\kappa \coloneqq - \Div\,\boldsymbol n$} is the local mean curvature while~\eqref{eq-setting-bc-h-free-surface} is the usual free surface kinematic equation. Note that~\eqref{eq-setting-bc-h-surface-tension} is a non-homogeneous Neumann boundary condition. By extension, for our general scalar equation~\eqref{eq-general}, we consider the general condition on free surface~$\freesurface(t)$:
\begin{equation}
\flux \cdot \boldsymbol n = q.\label{eq-setting-bc-h-diffusion}
\end{equation}


\subsubsection{On the substrate} On the substrate~$\Sigma_s$, we consider:
\begin{subequations}\label{eq-setting-bc-0}
\begin{align}
\velocity \cdot \boldsymbol n &= 0,\label{eq-setting-bc-0-no-penetration}
\\
\stress_{nt} &= k\velocity_t,\label{eq-setting-bc-0-friction}
\end{align}
\end{subequations}
where~\mbox{$\stress_{nt} \coloneqq \stress \cdot \boldsymbol n - (\boldsymbol n \cdot \stress \cdot \boldsymbol n)\boldsymbol n$} is the tangential part of the normal stress,
\mbox{$\velocity_t \coloneqq \velocity - (\velocity \cdot \boldsymbol n)\boldsymbol n$} is the tangential part of the velocity. Relation~\eqref{eq-setting-bc-0-no-penetration} expresses the \emph{no-penetration} of the flow material across the substrate surface while~\eqref{eq-setting-bc-0-friction} is the tangential
\emph{slip-friction} condition, where~$k$ is the friction coefficient, which is constant for laminar flows.

Condition~\eqref{eq-setting-bc-0} is the usual starting point of most \emph{motion by slices} approximate models. Note that a popular alternative is the \emph{no-slip} condition~\mbox{$\velocity = \boldsymbol 0$}, which is the starting point of the so-called \emph{lubrication} approximate models~\cite{LiuMei-1990,BalCra-1999,BerSarSmu-2014}. For simplicity, the later will not be considered and we assume~\eqref{eq-setting-bc-0} all along this paper.


\section{The framework}\label{sec-framework}

\subsection{Notations}

\subsubsection{Asymptotic analysis}\label{sec-framework-asymptotic}

Let~$H$, $L$ and~$U$ be the characteristic height, length and velocity, respectively. For instance~$L$ is the substrate length. The asymptotic analysis consists in exploiting the small, yet fixed, aspect ratio
\begin{equation*}
\epsilon \coloneqq H / L
\end{equation*}
in order to reduce the three-dimensional problem to a bidimensional one.
For this purpose, we are going to scale the spatial coordinates with this parameter, distinguishing the contribution from planar and horizontal dimensions. We then assume each \emph{dimensionless} unknown admits a Taylor expansion with respect to~$\epsilon$ at second order, which, to fix ideas, writes
\begin{equation}
\varphi = \varphi^{(0)} + \varphi^{(1)}\epsilon + \varphi^{(2)}\epsilon^2 + O(\epsilon^3). \label{eq-asymptotic-expansion}
\end{equation}
By convention, we consider that both the height and the~$x_3$-coordinate are exactly of order one while both the~$x_1$-coordinate and~$x_2$-coordinate do not depend on~$\epsilon$:
\begin{equation}
\boldsymbol x_s = L\boldsymbol x_s^{(0)},\quad x_3 = Lx_3^{(1)}\epsilon,\quad t = \dfrac{L}{U}t^{(0)},\label{eq-asymptotic-expansion-coords}
\end{equation}
where the~$s$ index refers to the planar components of a given variable,
e.g.~\mbox{$\boldsymbol x_s \coloneqq (x_1, x_2)$}. We deduce the following scaling for differential operators:
\begin{equation}
\nabla_s = \inv{L}\nabla_s^{(0)},\quad
\pdd{}{x_3} = \inv{\epsilon L}\pdd{}{x_3^{(1)}},\quad
\pddt{} = \dfrac{U}{L}\pdd{}{t^{(0)}}, \label{eq-asymptotic-expansion-diff}
\end{equation}
where~\mbox{$\nabla_s \coloneqq (\textpdd{}{x_1}, \textpdd{}{x_2})^\top$} is the planar gradient. We also define the planar divergence operator~\mbox{$\Div_s\square\coloneqq \textpdd{\square_{i_1\ldots i_{p-1}1}}{x_1} + \textpdd{\square_{i_1\ldots i_{p-1}2}}{x_2}$}, where~$p$ is the order of the given tensor, and the planar Laplacian~\mbox{$\Delta_s \coloneqq \nabla_s\cdot\nabla_s$}. Similarly, we assume that the dimensionless velocity field $\velocity / U$ and the dimensionless height field $\height / L$ satisfy an expansion of the form~\eqref{eq-asymptotic-expansion} with~\mbox{$\velocity[]_3^{(0)} = 0$} and \mbox{$\height^{(0)} = 0$}. This choice is justified by our wish to have $\epsilon U$ as the characteristic vertical speed and $H = \epsilon L$ as the characteristic height.

Observe that the notation~$\varphi^{(k)}$ makes it possible to immediately recognise the order with respect to $\epsilon$ of the considered term. For instance, the product~\mbox{$\varphi^{(k)}\psi^{(l)}$} will be associated to a term of order~\mbox{$O\left(\epsilon^{k+l}\right)$}, while the quotient~\mbox{$\varphi^{(k)}/\psi^{(l)}$} will be associated to a term of order~\mbox{$O\left(\epsilon^{k-l}\right)$}.

From now until \cref{sec-examples}, we assume that every field, parameter and differential operator is dimensionless, for convenience. In particular, when we refer to the dimensionless height (resp. vertical velocity), we mean the height scaled by $L$ (resp. $U$), in order to keep the property presented in the previous paragraph.

\subsubsection{Depth averaging}

Depth averaging is the final stage of the reduction of equations. It is defined through the following operator:

\begin{defi}[depth average]
Let $\varphi~\colon~\body\times\R_+ \to \R$ a continuous scalar field with respect to $x_3$. Its depth average is defined by the quantity
\begin{equation*}
\mean{\varphi}(\boldsymbol x_s, t) \coloneqq
\begin{dcases}
\inv{\height}\int_0^\height
\varphi((\boldsymbol x_s, x_3), t) \dx_3,\quad \forall (\boldsymbol x_s, t) \in \Omega \times \R_+ & \text{when } \height(\boldsymbol x_s, t) \neq 0,
\\[0.25em]
\varphi((\boldsymbol x_s, 0), t) &\text{otherwise.}
\end{dcases}
\end{equation*}
\end{defi}

This linear operator extends componentwise to any vector or tensor field. We also define the depth variance and covariance:
\begin{subequations}
\begin{defi}[depth variance and covariance]
Let $\varphi$ and $\psi$ be two scalar fields defined in $\body\times\R_+$, continuous with respect to $x_3$.
\begin{itemize}
\item The depth variance of $\varphi$ is defined by the quantity
\begin{equation*}
\var(\varphi) \coloneqq \mean{(\varphi - \mean{\varphi})^2}.
\end{equation*}
\item The depth covariance of $\varphi$ and $\psi$ is defined by the quantity
\begin{equation*}
\cov(\varphi, \psi) \coloneqq \mean{(\varphi - \mean{\varphi})(\psi - \mean{\psi})}.
\end{equation*}
\end{itemize}
In particular, $\cov(\varphi, \varphi) = \var(\varphi)$.
\end{defi}

Those operators satisfy the Cauchy--Schwartz inequality
\begin{equation}
\vert\cov(\varphi, \psi)\vert^2 \leqslant \var(\varphi)\var(\psi),\label{eq-asymptotic-mean-cauchy-schwarz}
\end{equation}
and the König--Huygens theorem-like
\begin{align}
\var(\varphi) &= \mean{\varphi^2} - \mean{\varphi}^2,
\\
\cov(\varphi, \psi) &= \mean{\varphi\psi} - \mean{\varphi}\mean{\psi}.\label{eq-asymptotic-mean-cov-konig-huygens}
\end{align}
\end{subequations}
Note the following useful identities:
\begin{equation}
\mean{x_3} = \height / 2
\quad
\text{ and }
\quad
\var{x_3} = \height^2 / 12.\label{eq-asymptotic-mean-examples}
\end{equation}

\subsection{Analysis of the general scalar problem}\label{sec-asymptotic-diffusion}

Considering equation~\eqref{eq-general} together with the Neumann condition~\eqref{eq-setting-bc-h-diffusion}, let us perform the asymptotic analysis of the following general scalar problem:

\begin{enonce}
{Problem}[general scalar]\label{pb-asymptotic-diffusion}
Find $\unknown$, defined in $\body \times \R_+$, such that for any time $t \in \R_+$
\begin{subequations}\label{eq-asymptotic-diffusion}
\begin{empheq}[left=\empheqbiglbrace]{alignat=2}
\mathcal{A}(\unknown(t)) - \Div\,{\flux(t)} &= \boldsymbol 0\; &\text{ in }& \body(t),\label{eq-asymptotic-diffusion-pde}
\\
\flux(t) \cdot \boldsymbol n &= q\; &\text{ on }& \Sigma_f(t).\label{eq-asymptotic-diffusion-bc}
\end{empheq}
\end{subequations}
\end{enonce}

Recall that this scalar problem extends componentwise to a vector or a tensor one. Assume that~$\unknown$ is a solution of~\eqref{eq-asymptotic-diffusion} and that the following expansions hold:
\begin{subequations}\label{eq-asymptotic-diffusion-expansion}
\begin{align}
\mathcal A &= \mathcal A^{(-1)}\inv{\epsilon} + \mathcal A^{(0)} + \mathcal A^{(1)}\epsilon + O(\epsilon^2), \\
\flux &= \flux^{(-1)}\inv{\epsilon} + \flux^{(0)} + \flux^{(1)}\epsilon + \flux^{(2)}\epsilon^2 + O(\epsilon^3), \\
q &= q^{(-1)}\inv{\epsilon} + q^{(0)} + q^{(1)}\epsilon + q^{(2)}\epsilon^2 + O(\epsilon^3)
\end{align}
\end{subequations}
Then~\eqref{eq-asymptotic-diffusion-pde} leads to
\begin{equation*}
- \pdd{\flux[]_3^{(-1)}}{x_3^{(1)}}\inv{\epsilon^2} + \sum_{k=-1}^{1}\paren{\mathcal A^{(k)} - \Div_s\flux_s^{(k)} - \pdd{\flux[]_3^{(k+1)}}{x_3^{(1)}}}\epsilon^k + O(\epsilon^2) = 0,
\end{equation*}
while the boundary condition~\eqref{eq-asymptotic-diffusion-bc} becomes
\begin{equation*}
\paren{\flux[]_3^{(-1)} - q^{(-1)}}\inv{\epsilon} + \sum_{k=-1}^{1}\paren{\flux[]_3^{(k+1)} - \flux_s^{(k)} \cdot \nabla_s\height^{(1)} - q^{(k+1)}}\epsilon^{k+1} + O(\epsilon^2) = 0\; \text{ on } \freesurface(t).
\end{equation*}
By identification, we get, for any~\mbox{$k \in \{-1,0,1\}$}
\begin{subequations}\label{eq-asymptotic-diffusion-identification-pb}
\begin{alignat}{2}
\pdd{\flux[]_3^{(-1)}}{x_3^{(1)}} &= 0,&&\label{eq-asymptotic-diffusion-identification--1}
\\
\mathcal A^{(k)} - \Div_s\flux_s^{(k)} &= \pdd{\flux[]_3^{(k+1)}}{x_3^{(1)}},&&\label{eq-asymptotic-diffusion-identification}
\\
\flux[]_3^{(-1)} &= q^{(-1)} &&\text{ on } \freesurface(t),\vphantom{\pdd{\flux[]_3^{(-1)}}{x_3^{(1)}}}\label{eq-asymptotic-diffusion-identification-cl-freesurface--1}
\\
\flux[]_3^{(k+1)} &= \nabla_s\height^{(1)} \cdot \flux_s^{(k)} + q^{(k + 1)}\; &&\text{ on } \freesurface(t). \vphantom{\pdd{\flux[]_3^{(-1)}}{}}\label{eq-asymptotic-diffusion-identification-cl-freesurface-k}
\end{alignat}
\end{subequations}
From~\eqref{eq-asymptotic-diffusion-identification--1} and~\eqref{eq-asymptotic-diffusion-identification-cl-freesurface--1}, we also deduce
\mbox{$\flux[]_3^{(-1)} = q^{(-1)}$ in $\Omega \times \R_+$}.

\begin{defi}[hypotheses $\mathcal{H}^{(k)}$]
For any~\mbox{$k \in \{-1, 0, 1\}$}, we say hypothesis~$\mathcal{H}^{(k)}$ holds if and only if~\mbox{$\mathcal A^{(k)} = 0$} and~\mbox{$\flux_s^{(k)} = \boldsymbol 0$}.
\end{defi}

\goodbreak
\begin{lemm}[vertical flux and Neumann BC]\label{lem-asymptotic-diffusion-substrate-flux}
Assume expansions~\eqref{eq-asymptotic-diffusion-expansion} are satisfied. Then,
\begin{subequations}
\begin{itemize}
\item the term of order $-1$ of the vertical flux is given by
\begin{equation}
\flux[]_3^{(-1)} = q^{(-1)}\quad \text{in } \Omega \times \R_+.\label{eq-asymptotic-vertical-stress--1}
\end{equation}
\item if hypothesis~$\mathcal{H}^{(k)}$ also holds, for~\mbox{$k \in \{-1, 0, 1\}$},
\begin{equation}
\flux[]_3^{(k+1)} = q^{(k+1)}\quad \text{in } \Omega \times \R_+.\label{eq-asymptotic-vertical-flux}
\end{equation}
\end{itemize}
\end{subequations}
\end{lemm}

\begin{proof}
The first point has already been proven in the main text. Under hypothesis $\mathcal{H}^{(k)}$, \eqref{eq-asymptotic-diffusion-identification} reduces to
\begin{equation*}
\pdd{\flux[]_3^{(k+1)}}{x_3^{(1)}} = 0 \quad \text{in } \Omega \times \R_+,
\end{equation*}
whence the result, from the boundary condition~\eqref{eq-asymptotic-diffusion-identification-cl-freesurface-k}.
\end{proof}

This lemma is the starting point for obtaining the property of motion by slices (see~\cref{cor-motion-by-slices}).

It remains to average in depth the problem defined by equations~\eqref{eq-asymptotic-diffusion-identification--1}--\eqref{eq-asymptotic-diffusion-identification} together with boundary conditions~\eqref{eq-asymptotic-diffusion-identification-cl-freesurface--1}--\eqref{eq-asymptotic-diffusion-identification-cl-freesurface-k}. This can be performed thanks to~\cref{cor-asymptotic-leibniz-nabla} in the appendix, which gives the

\begin{enonce}{Problem}[reduced general scalar]\label{pb-asymptotic-diffusion-eps}
Find $\mean{\unknown^{(0)}}$ and $\mean{\unknown^{(1)}}$, defined in $\Omega \times \R_+$, such that for any time $t \in \R_+$ and any index $k \in \{-1, 0, 1\}$
\begin{equation}
\height^{(1)}\mean{\mathcal A^{(k)}} - \Div_s\paren{\height^{(1)}\mean{\flux_s^{(k)}}} = q^{(k+1)} - \flux[]_3^{(k+1)}(x_3 = 0).\label{eq-asymptotic-diffusion-eps}
\end{equation}
\end{enonce}

\begin{theo}[asymptotic analysis of the general problem]\label{th-asymptotic-diffusion}
Let $u$ be the solution of \cref{pb-asymptotic-diffusion}, which is assumed to be statisfying expansion~\eqref{eq-asymptotic-diffusion-expansion}, and $\mean{u^{(0)}}$ and $\mean{u^{(1)}}$ be solutions of \cref{pb-asymptotic-diffusion-eps}. Then,
\begin{equation}
\mean{u} = \mean{u^{(0)}} + \mean{u^{(1)}}\epsilon + O(\epsilon^2).
\end{equation}
In other words, the reconstructed solution $\mean{u^{(0)}} + \mean{u^{(1)}}\epsilon$ is an approximation of order~$O(\epsilon^2)$ of the averaged exact solution $\mean{u}$.
\end{theo}

\begin{coro}[motion by slices]\label{cor-motion-by-slices}
Under the same assumptions than the above~\cref{th-asymptotic-diffusion}, let us assume in addition that~{$\flux = \nu\nabla\unknown + \boldsymbol m$}, where~{$\nu > 0$} and~{$\boldsymbol m$} is a vector field that does not depend on~$\nabla\unknown$. If~$\boldsymbol m$ expands as~{$\boldsymbol m = \boldsymbol m^{(0)} + \boldsymbol m^{(1)}\epsilon + \boldsymbol O(\epsilon^2)$}, {$q^{(-1)} = 0$} and~{$\nu = O(1)$}, then
\begin{equation*}
\pdd{\unknown^{(0)}}{x_3^{(1)}} = 0. \label{eq-asymptotic-dt-block-0}
\end{equation*}
\end{coro}

\begin{proof}
By construction, \mbox{$\flux[]_3^{(-1)} = \nu\textpdd{\unknown^{(0)}}{x_3^{(1)}}$}. The result follows then from~\eqref{eq-asymptotic-vertical-stress--1}.
\end{proof}

\begin{rema}[no-slip condition]\label{rem-lubrication}
More generally, if~$\boldsymbol m$ expands as~$\boldsymbol m = \boldsymbol m^{(-1)}/\epsilon + \boldsymbol m^{(0)} + \boldsymbol m^{(1)}\epsilon + \boldsymbol m^{(2)}\epsilon^2 + \boldsymbol O(\epsilon^3)$ and~\mbox{$\nu = O(\epsilon^n)$}, with~\mbox{$n \in \{0, \ldots, 3\}$}, then
\begin{equation}
\pdd{\unknown^{(0)}}{x_3^{(1)}} = \inv{\nu}\paren{q^{(n-1)} - m_3^{(n-1)}},\; \text{that is } \mean{\unknown^{(0)}} = \unknown^{(0)}(x_3 = 0) + \dfrac{\height^{(1)}}{2\nu}q^{(n-1)} - \bigmean{\int_0^{x_3} \inv{\nu}m_3^{(n-1)} \dx_3}.
\end{equation}
Instead of the Neumann boundary condition~\eqref{eq-asymptotic-diffusion-bc}, we could consider a Dirichlet condition on~$u$. In that case, we do no more have information on the normal flux~\mbox{$\flux \cdot \boldsymbol n$} on the substrate. This could be a difficulty since this term is involved in the right-hand-side of~\eqref{eq-asymptotic-diffusion-eps}. Note that a popular case of such Dirichlet condition is the no-slip condition on the substrate, which the starting point of lubrication shallow models. In that case, there is no more motion by slices, but this asymptotic expansion still allows to obtain an explicit expression of the zeroth order depth average of~$\unknown$. Moreover, from~\eqref{eq-asymptotic-vertical-flux} and under assumption~$\mathcal{H}^{(k)}$ (see~\cref{lem-asymptotic-diffusion-substrate-flux}), it is possible to extend the result to any $k \in \{-1,0,1\}$ as
\begin{equation}
\mean{\unknown^{(k+1)}} = \unknown^{(k+1)}(x_3 = 0) + \dfrac{\height^{(1)}}{2\nu}q^{(n+k)} - \bigmean{\int_0^{x_3} \inv{\nu}m_3^{(n+k)} \dx_3}, \label{eq-asymptotic-dt-block-k}
\end{equation}
See e.g.~\cite[\S3.1]{BalCra-1999} or~\cite{BerSarSmu-2014} for practical examples of lubrication problems.
\end{rema}

\begin{coro}[asymptotic analysis of the pressure]\label{cor-asymptotic-diffusion-pressure}
For any~\mbox{$i\in\{1,2,3\}$}, assume the~$i$-th component of a vector field~$\ubold$ is solution of a problem of the form~\eqref{eq-asymptotic-diffusion} with vector flux~$(\boldsymbol\sigma_{ij})_{1 \leqslant j \leqslant 3}$ and operator~$\mathcal{A}_i$. Assume in addition that these three problems satisfy the expansion~\eqref{eq-asymptotic-diffusion-expansion}. Assume also that the tensor flux writes~\mbox{$\boldsymbol\sigma = \boldsymbol\tau - p\boldsymbol I$}, where~$\boldsymbol\tau$ is a tensor field and~$p$ is a scalar field, called the pressure, both of them being of order at most~\mbox{$O(1/\epsilon)$}. Then the depth average of the pressure~$p$ is given by
\begin{align*}
\mean{p^{(-1)}} &= \mean{\tau_{33}^{(-1)}} - q^{(-1)},
\\
\mean{p^{(k+1)}} &= \mean{\tau_{33}^{(k+1)}} - q^{(k+1)} - \bigmean{\int_{x_3}^{\height}\paren{\Div_s\boldsymbol \tau_{3s}^{(k)} - \mathcal{A}_3^{(k)}}\dx_3} - \boldsymbol\tau_{3s}^{(k)}(x_3=\height) \cdot \nabla_s\height^{(1)},
\end{align*}
for any $k \in \{-1,0,1\}$.
\end{coro}

\begin{proof}
The first result is immediate from~\eqref{eq-asymptotic-vertical-stress--1}. From~\eqref{eq-asymptotic-diffusion-identification}, we get
\begin{equation*}
\pdd{p^{(k+1)}}{x_3} = \pdd{\tau_{33}^{(k+1)}}{x_3} + \Div_s\boldsymbol \tau_{3s}^{(k)} - \mathcal{A}_3^{(k)},
\end{equation*}
and then, by integrating between~$x_3$ and~$\height$
\begin{equation*}
p^{(k+1)}(x_3=\height) - p^{(k+1)} = \tau_{33}^{(k+1)}(x_3=\height) - \tau_{33}^{(k+1)} + \int_{x_3}^{\height}\paren{\Div_s\boldsymbol \tau_{3s}^{(k)} - \mathcal{A}_3^{(k)}}\dx_3.
\end{equation*}
From the Neumann boundary condition~\eqref{eq-asymptotic-diffusion-identification-cl-freesurface-k}, we obtain the relation
\begin{equation*}
\tau_{33}^{(k+1)}(x_3=\height) - p^{(k+1)}(x_3=\height) = \boldsymbol\tau_{3s}^{(k)}(x_3=\height) \cdot \nabla_s\height^{(1)} + q^{(k+1)}.
\end{equation*}
A depth integration allows then to conclude.
\end{proof}

\begin{rema}[simplifications]
If~$\boldsymbol \tau_{3s}^{(k)}$ is independent upon~$x_3$, we have the much simpler result
\begin{equation*}
\mean{p^{(k+1)}} = \mean{\tau_{33}^{(k+1)}} - q^{(k+1)} - \Div_s\paren{\height^{(1)}\boldsymbol \tau_{3s}^{(k)}} + \half[\height^{(1)}]\Div_s\boldsymbol \tau_{3s}^{(k)} + \bigmean{\int_{x_3}^{\height}\mathcal{A}_3^{(k)}\dx_3}.
\end{equation*}
Similarly, if $\mathcal{A}_3^{(k)}$ is also independent of $x_3$, we have the even simpler result
\begin{equation*}
\mean{p^{(k+1)}} = \mean{\tau_{33}^{(k+1)}} - q^{(k+1)} - \Div_s\paren{\height^{(1)}\boldsymbol \tau_{3s}^{(k)}} + \half[\height^{(1)}]\Div_s\boldsymbol \tau_{3s}^{(k)} + \half[\height^{(1)}]\mathcal{A}_3^{(k)}.
\end{equation*}
\end{rema}

\subsection{Time derivatives}

The aim of this section is to derive the depth average of the time derivatives that we introduced in \cref{subsec-setting-equations}, namely the material derivative $\textmaterialdd{}$ and the vector and tensor Gordon--Schowalter derivatives $\textgsdd[a]{}$, and express them as functions of their equivalents in the plane $\Omega$. We will do so under the assumptions that lead to motion by slices (see~\cref{cor-motion-by-slices}):
\begin{equation}
\pdd{\velocity_s^{(0)}}{x_3^{(1)}} = \boldsymbol 0, \label{eq-asymptotic-dt-block}
\end{equation}
and the fact that the velocity field~$\velocity$ is divergence-free (incompressibility) and satisfies to the no-penetration constraint on the substrate. To exploit this, we will first need to perform the asymptotic analysis of the divergence-free condition.

\subsubsection{Divergence-free condition}\label{sec-dt-incompr}

Let us start by the asymptotic analysis of the divergence-free condition
\begin{equation}
\Div\,\velocity = 0, \label{eq-asymptotic-pb-mass}
\end{equation}
under the boundary conditions
\begin{subequations}\label{eq-asymptotic-pb-mass-cl}
\begin{alignat}{2}
\velocity \cdot \boldsymbol n = 0 \iff \velocity[]_3 &= 0 &\text{ on }& \substrate, \label{eq-asymptotic-pb-mass-cl-slip}
\\
\materialdd[s]{\height} &= \velocity[]_3 &\text{ on }& \freesurface(t), \label{eq-asymptotic-pb-mass-cl-freesurface}
\end{alignat}
\end{subequations}
where we noted $\materialdd[s]{} \coloneqq \pddt{} +
\advection[s]{\velocity_s}$ the material derivative advected by the velocity $\velocity_s$.
\begin{prop}[averaged vertical velocity]\label{prop-asymptotic-mass-v3}
Assume that $\velocity$ is divergence-free and that the property of motion by slices~\eqref{eq-asymptotic-dt-block} and the boundary conditions~\eqref{eq-asymptotic-pb-mass-cl} are satisfied. Then the vertical velocity $\velocity[]_3$ is given by:
\begin{equation}
\velocity[]_3 = - \Div_s\paren{\velocity_s^{(0)}}x_3 - \mean[x_3]{\Div_s\paren{\velocity_s^{(1)}}}\epsilon + O(\epsilon^3),\; \text{that is}\; \velocity[]_3^{(1)} = - \Div_s\paren{\velocity_s^{(0)}}x_3^{(1)}.
\end{equation}
\end{prop}
\begin{proof}
Note that $\Div\,\velocity = \Div_s\velocity_s + \textpdd{\velocity[]_3}{x_3} = 0$. Therefore, the result follows directly from an integration of the equation~\eqref{eq-asymptotic-pb-mass} using the hypothesis of motion by slices and the boundary condition~\eqref{eq-asymptotic-pb-mass-cl-slip}.
\end{proof}
\begin{coro}[averaged free surface evolution]\label{cor-asymptotic-mass-height}
Under the same conditions as in the above proposition, the height $\height$ is solution of the equation
\begin{equation}
\pddt{\height} + \Div_s\paren{\height\velocity_s^{(0)}} + \Div_s\paren{\height\mean{\velocity_s^{(1)}}}\epsilon = O(\epsilon^3). \label{eq-asymptotic-mass-H}
\end{equation}
\end{coro}
\begin{proof}
Immediate from~\cref{prop-asymptotic-mass-v3}, \cref{cor-asymptotic-leibniz-nabla} (see \cref{sec-appendix}) and boundary condition~\eqref{eq-asymptotic-pb-mass-cl-freesurface}.
\end{proof}

By identification, we immediately obtain the equation satisified by the height at order $O(\epsilon)$:
\begin{equation}
\pdd{\height^{(1)}}{t} + \Div_s\paren{\height^{(1)}\velocity_s^{(0)}} = 0. \label{eq-asymptotic-mass-H-1}
\end{equation}
Solving~\eqref{eq-asymptotic-mass-H-1} instead of~\eqref{eq-asymptotic-mass-H} leads to an error of order $O(\epsilon)$ since\footnote{As a reminder, the height was scaled by the characteristic length $L$, not by the characteristic height $H$, as mentioned in \cref{sec-framework-asymptotic}.} $\height/\epsilon = \height^{(1)} + O(\epsilon)$.


\subsubsection{Lagrangian derivative}

The depth average of the Lagrangian derivative is expressed as a function of the planar material derivative, advected by the velocity $\velocity_s^{(0)}$, noted
\[
\materialdd[s][(0)]{\varphi} \coloneqq \pddt{\varphi} +
\advection[s]{\velocity_s^{(0)}}\varphi,
\]
for any scalar field $\varphi$ defined in $\Omega\times \R_+$. It extends componentwise to any vector or tensor field. Thanks to the new mass conservation law~\eqref{eq-asymptotic-mass-H}, we deduce the
\begin{prop}[averaged Lagrangian derivative]\label{prop-asymptotic-time-derivative-material}
Assume that $\velocity$ is divergence-free and that the property of motion by slices~\eqref{eq-asymptotic-dt-block} and the boundary conditions~\eqref{eq-asymptotic-pb-mass-cl} are satisfied. Then, any scalar field $\unknown$ defined in $\body \times \R_+$ satisfies the relation
\begin{equation}
\bigmean{\materialdd{\unknown}} = \materialdd[s][(0)]{\mean{\unknown}} + O(\epsilon),
\end{equation}
Moreover, if we assume $\cov(\nabla_s\unknown, \velocity_s^{(1)}) = O(\epsilon)$, then we have a better approximation:
\begin{equation}
\bigmean{\materialdd{\unknown}} = \materialdd[s][(0)]{\mean{\unknown}} + \advection[s]{\mean{\velocity_s^{(1)}}}\mean{\unknown}\epsilon + O(\epsilon^2),
\end{equation}
\end{prop}

\begin{proof}
We show this property in the second case, the first result being independent of the assumption $\cov(\nabla_s\unknown, \velocity_s^{(1)}) = O(\epsilon)$. The calculation is performed in two steps, by noting that $\textmaterialdd{} = \textmaterialdd[s]{} + \velocity[]_3\textpdd{}{x_3} = \textmaterialdd[s]{} + \velocity[]_3^{(1)}\textpdd{}{x_3^{(1)}} + O(\epsilon)$ and by using the linearity of the operator~$\mean{\,\cdot\,}$.
\begin{enumerate}
\item Two direct applications of the Leibniz formula (\cref{cor-asymptotic-leibniz-L} in \cref{sec-appendix}) and one of the formula~\eqref{eq-asymptotic-mean-cov-konig-huygens} gives
\begin{multline*}
\height\bigmean{\materialdd[s]{\unknown}} = \materialdd[s][(0)]{}\paren{\height\mean{\unknown}} - \unknown(x_3 = \height)\materialdd[s][(0)]{\height} + \height\mean{\advection[s]{\velocity_s^{(1)}}\unknown}\epsilon + O(\epsilon^3) \\
= \height\paren{\materialdd[s][(0)]{\mean{\unknown}} + \advection[s]{\mean{\velocity_s^{(1)}}}\mean{\unknown}\epsilon} + \height\cov(\nabla_s\unknown, \velocity_s^{(1)})\epsilon \\
+ \paren{\mean{\unknown} - \unknown(x_3 = \height)}\paren{\materialdd[s][(0)]{\height} + \advection[s]{\mean{\velocity_s^{(1)}}}\height\epsilon} + O(\epsilon^3).
\end{multline*}
\item A direct application of the appendix formula~\eqref{eq-asymptotic-leibniz-kpp}, after expansion of the vertical velocity using the~\cref{prop-asymptotic-mass-v3} and use of the hypothesis of motion by slices, gives
\begin{align*}
\height\bigmean{\velocity[]_3\pdd{\unknown}{x_3}} &= - \paren{\height\Div_s\velocity_s^{(0)} + \height\mean{\Div_s\velocity_s^{(1)}}\epsilon}\bigmean{x_3\pdd{\unknown}{x_3}} + O(\epsilon^3) \\
&= \paren{\height\Div_s\velocity_s^{(0)} + \height\mean{\Div_s\velocity_s^{(1)}}\epsilon} \paren{\mean{\unknown} - \unknown(x_3 = \height)} + O(\epsilon^3).
\end{align*}
\end{enumerate}
The sum of the two above results and the application of equality~\eqref{eq-asymptotic-mass-H} allow to conclude.
\end{proof}

This result still holds true when replacing $\unknown$ with a vector or a tensor field, in particular with the velocity field $\velocity$.

\subsubsection{Gordon--Schowalter derivative of a vector}

First, the depth average of the Gordon--Schowalter derivative of a vector is expressed as a function of its plane equivalent, with respect to the zeroth-order planar velocity $\velocity_{s}^{(0)}$:
\begin{equation}
\gsdd[a][(0)]{\boldsymbol{w}_s} \coloneqq \materialdd[s][(0)]{\boldsymbol{w}_s} - W_s(\velocity_{s}^{(0)}) \cdot \boldsymbol{w}_s - aD_s(\velocity_{s}^{(0)}) \cdot \boldsymbol{w}_s,
\end{equation}
for any planar vector $\boldsymbol{w}_s$ defined in $\Omega\times\R_+$.

\begin{theo}[Gordon--Schowalter derivative of a vector]\label{thm-asymptotic-mean-gsdd-vect}
Let $\ubold$ be a vector field defined in $\body \times \R_+$ such that $\var(\ubold) = O(\epsilon)$ (closure assumption). Assume that $\velocity$ is divergence-free and that the property of motion by slices~\eqref{eq-asymptotic-dt-block} and the boundary conditions~\eqref{eq-asymptotic-pb-mass-cl} are satisfied. Then the depth-average of the Gordon--Schowalter derivative of $\ubold$ is given by
\begin{subequations}
\begin{multline}
\paren{\bigmean{\gsdd[a]{\ubold}}}_s = \gsdd[a][(0)]{\mean{\ubold_s}} - \half[1 + a]\bigmean{\pdd{\velocity_{s}^{(1)}}{x_3^{(1)}}u_3} - \left(\bigmean{\brackets*{W_s\paren{\velocity_{s}^{(1)}} + aD_s\paren{\velocity_{s}^{(1)}}}\cdot\ubold_s} \vphantom{\bigmean{\pdd{\velocity_{s}^{(2)}}{x_3^{(1)}}u_3}}\right.\\
\left.{}+ \half[1+a]\bigmean{\pdd{\velocity_{s}^{(2)}}{x_3^{(1)}}u_3} + \dfrac{(1 - a)\height^{(1)}}{4}\mean{u_3}\nabla_s\Div_s\velocity_s^{(0)}\right)\epsilon + \boldsymbol O(\epsilon^2),
\end{multline}
\begin{multline}
\paren{\bigmean{\gsdd[a]{\ubold}}}_3 = \materialdd[s][(0)]{\mean{u_3}} + \half[1-a]\bigmean{\pdd{\velocity_{s}^{(1)}}{x_3^{(1)}}\cdot \ubold_s} + a\Div_s(\velocity_s^{(0)})\mean{u_3} \\
+ \paren{\half[1-a]\bigmean{\pdd{\velocity_{s}^{(2)}}{x_3^{(1)}}\cdot \ubold_s} + \dfrac{(1 + a)\height^{(1)}}{4}\nabla_s\Div_s\velocity_s^{(0)} \cdot \mean{\ubold_s} + a\mean{\Div_s(\velocity_s^{(1)})u_3}}\epsilon + O(\epsilon^2).
\end{multline}
\end{subequations}
\end{theo}

\begin{proof}
Recall that
\begin{equation*}
\gsdd[a]{\ubold} \coloneqq \materialdd{\ubold} - W(\velocity) \cdot \ubold - aD(\velocity) \cdot \ubold.
\end{equation*}
As in the previous demonstration, we will treat each of the terms separately, and we use the previous result for the first term. The expansions~\eqref{eq-asymptotic-expansion}, \eqref{eq-asymptotic-expansion-coords} and~\eqref{eq-asymptotic-expansion-diff} applied to the velocity field, and the motion by slices assumption lead to
%\begin{subequations}
\begin{gather*}
2W(\velocity) =
\begin{pmatrix}
2W_s(\velocity_s^{(0)}) + 2W_s(\velocity_s^{(1)})\epsilon & \pdd{\velocity_s^{(1)}}{x_3^{(1)}} + \paren{\pdd{\velocity_s^{(2)}}{x_3^{(1)}} - \nabla_s\velocity[]_3^{(1)}}\epsilon \\
\paren{\nabla_s\velocity[]_3^{(1)} - \pdd{\velocity_s^{(2)}}{x_3^{(1)}}}\epsilon - \pdd{\velocity_s^{(1)}}{x_3^{(1)}} & 0
\end{pmatrix} + \boldsymbol O(\epsilon^2), \\
2D(\velocity) =
\begin{pmatrix}
2D_s(\velocity_s^{(0)}) + 2D_s(\velocity_s^{(1)})\epsilon &  \pdd{\velocity_s^{(1)}}{x_3^{(1)}} + \paren{\pdd{\velocity_s^{(2)}}{x_3^{(1)}} + \nabla_s\velocity[]_3^{(1)}}\epsilon \\
\paren{\nabla_s\velocity[]_3^{(1)} + \pdd{\velocity_s^{(2)}}{x_3^{(1)}}}\epsilon + \pdd{\velocity_s^{(1)}}{x_3^{(1)}} & 2\pdd{\velocity[]_3^{(1)}}{x_3^{(1)}} + 2\pdd{\velocity[]_3^{(2)}}{x_3^{(1)}}\epsilon
\end{pmatrix} + \boldsymbol O(\epsilon^2).
\end{gather*}
%\end{subequations}
Therefore, the second and third terms in the definition of the Gordon--Schowalter derivative are
%\begin{subequations}
\begin{gather*}
W(\velocity) \cdot \ubold =
\begin{pmatrix}
W_s(\velocity_s^{(0)}) \cdot \ubold_s + \half\pdd{\velocity_s^{(1)}}{x_3^{(1)}}u_3 + \paren{W_s(\velocity_s^{(1)}) \cdot \ubold_s + \pdd{\velocity_s^{(2)}}{x_3^{(1)}}\half[u_3] - \nabla_s\velocity[]_3^{(1)}\half[u_3]}\epsilon + \boldsymbol O(\epsilon^2) \\
- \half\pdd{\velocity_s^{(1)}}{x_3^{(1)}} \cdot \ubold_s + \half[\epsilon] \paren{\nabla_s\velocity[]_3^{(1)} - \pdd{\velocity_s^{(2)}}{x_3^{(1)}}} \cdot \ubold_s + O(\epsilon^2)
\end{pmatrix}, \\
D(\velocity) \cdot \ubold =
\begin{pmatrix}
D_s(\velocity_s^{(0)}) \cdot \ubold_s + \half\pdd{\velocity_s^{(1)}}{x_3^{(1)}}u_3 + \paren{D_s(\velocity_s^{(1)}) \cdot \ubold_s + \pdd{\velocity_s^{(2)}}{x_3^{(1)}}\half[u_3] + \nabla_s\velocity[]_3^{(1)}\half[u_3]}\epsilon + \boldsymbol O(\epsilon^2) \\
\half\pdd{\velocity_s^{(1)}}{x_3^{(1)}} \cdot \ubold_s + \pdd{\velocity[]_3^{(1)}}{x_3^{(1)}}u_3 + \paren{\half\brackets*{\nabla_s\velocity[]_3^{(1)} + \pdd{\velocity_s^{(2)}}{x_3^{(1)}}} \cdot \ubold_s + \pdd{\velocity[]_3^{(2)}}{x_3^{(1)}}u_3}\epsilon + O(\epsilon^2)
\end{pmatrix},
\end{gather*}
%\end{subequations}
where we used the~\cref{prop-asymptotic-mass-v3}. Depth averaging is easily done using the motion by slices assumption. However, it is necessary to treat the terms consisting in a product of two terms, a priori dependent on the coordinate $x_3$, carefully. We deal with the term $u_3\nabla_s\velocity[]_3^{(1)} = -x_3^{(1)}\nabla_s\Div_s\paren{\velocity_{s}^{(0)}}u_3$ for the example, whose average can be calculated through $\mean{x_3^{(1)}u_3}$. By successively applying the formulas~\eqref{eq-asymptotic-mean-cov-konig-huygens} and~\eqref{eq-asymptotic-mean-examples}, we get
\begin{equation*}
\mean{x_3^{(1)}u_3} = \half[\height^{(1)}]\mean{u_3} + \cov(x_3^{(1)}, u_3) \approx \half[\height^{(1)}]\mean{u_3}.
\end{equation*}
The closure assumption allows then to prove the approximation since, from the result~\eqref{eq-asymptotic-mean-cauchy-schwarz} and the formula~\eqref{eq-asymptotic-mean-examples},
\begin{equation*}
\abs*{\cov(x_3^{(1)}, u_3)}^2 \leqslant \var(x_3^{(1)})\var(u_3) = \dfrac{\height^{(1)^2}}{12}\var(u_3) = O(\epsilon).
\end{equation*}
The linearity of $\mean{\,\cdot\,}$ and the~\cref{prop-asymptotic-time-derivative-material} on the asymptotic analysis of the Lagrangian derivative finally lead to the conclusion.
\end{proof}

\begin{rema}[average of a product]
The expression of the depth-average of the Gordon--Schowalter derivative involves averages of products, which is not desirable for closing the reduced problem: we would prefer instead a product of averaged quantities. The simplification of those products is actually case-dependent: one at least of the two quantities involved by the product should be independent upon~$x_3$. For instance, either the first-order planar velocity~$\velocity_{s}^{(1)}$ or the vector field~$\boldsymbol{u}$ involved by the Gordon--Schowalter derivative should be independent upon~$x_3$.
\end{rema}

%\end{document}

\subsubsection{Gordon--Schowalter derivative of a tensor}

Then, the depth average of the Gordon--Schowalter derivative of a tensor is expressed as a function of its plane equivalent, with respect to the zeroth-order planar velocity $\velocity_{s}^{(0)}$:
\begin{equation}
\gsdd[a][(0)]{\boldsymbol\sigma_{ss}} \coloneqq \materialdd[s][(0)]{{\boldsymbol\sigma_{ss}}} - 2\sym\brackets*{W_s(\velocity_{s}^{(0)}) \cdot {\boldsymbol\sigma_{ss}} - aD_s(\velocity_{s}^{(0)}) \cdot {\boldsymbol\sigma_{ss}}},
\end{equation}
for any planar tensor $\boldsymbol\sigma_{ss}$ in $\Omega\times\R_+$.

\begin{theo}[Gordon--Schowalter derivative of a tensor]\label{thm-asymptotic-mean-gsdd-tens}
Let $\boldsymbol\tau$ be a symmetric tensor field defined in $\body \times \R_+$, whose zeroth order vertical variation with respect to $\epsilon$ is negligible, that is $\var(\boldsymbol\tau^{(0)}) = O(\epsilon)$ (closure assumption). Assume that $\velocity$ is divergence-free and that the property of motion by slices~\eqref{eq-asymptotic-dt-block} and the boundary conditions~\eqref{eq-asymptotic-pb-mass-cl} are satisfied. Then the depth-average of the Gordon--Schowalter derivative of $\boldsymbol\tau$ is given by
\begin{subequations}
\begin{multline}
\paren{\bigmean{\gsdd[a]{\boldsymbol\tau}}}_{ss} = \gsdd[a][(0)]{\mean{\boldsymbol\tau_{ss}}} - (1+a)\sym\paren{\bigmean{\pdd{\velocity_{s}^{(1)}}{x_3^{(1)}} \otimes \boldsymbol\tau_{3s}}} \\
- \brackets*{\sym\left(\bigmean{\braces*{W_s\paren{\velocity_s^{(1)}} + a D_s\paren{\velocity_s^{(1)}}} \cdot \tau_{ss}} + (1+a)\bigmean{\pdd{\velocity_{s}^{(2)}}{x_3^{(1)}} \otimes \boldsymbol\tau_{3s}}\right. \\
\left.\vphantom{\bigmean{\pdd{\velocity_{s}^{(2)}}{x_3^{(1)}} \otimes \boldsymbol\tau_{3s}}}+ \half[(1-a)\height^{(1)}]\nabla_s\Div_s\velocity_s^{(0)} \otimes \mean{\boldsymbol\tau_{3s}}\right)}\epsilon + \boldsymbol O(\epsilon^2),
\end{multline}
\begin{multline}
\paren{\bigmean{\gsdd[a]{\boldsymbol\tau}}}_{s 3} = \gsdd[a][(0)]{\mean{\boldsymbol\tau_{s3}}} + \half[1-a]\bigmean{\pdd{\velocity_s^{(1)}}{x_3^{(1)}} \cdot \boldsymbol\tau_{ss}} - \half[1+a]\bigmean{\pdd{\velocity_s^{(1)}}{x_3^{(1)}}\tau_{33}} + a\Div_s\paren{\velocity_s^{(0)}}\mean{\boldsymbol\tau_{s3}} \\
+ \braces*{-\bigmean{\brackets*{W_s\paren{\velocity_s^{(1)}} + aD_s\paren{\velocity_s^{(1)}}} \cdot \boldsymbol\tau_{s3}} + \dfrac{\height^{(1)}}{4}\brackets*{(1+a)\mean{\boldsymbol\tau_{ss}} - (1-a)\mean{\tau_{33}}I_{ss}} \cdot \nabla_s\Div_s\velocity_s^{(0)} \\
+ \half[1-a]\bigmean{\boldsymbol\tau_{ss} \cdot \pdd{\velocity_s^{(2)}}{x_3^{(1)}}} - \half[1+a]\bigmean{\pdd{\velocity_s^{(2)}}{x_3^{(1)}}\tau_{33}} + a \mean{\Div_s\paren{\velocity_s^{(1)}}\boldsymbol\tau_{3s}}}\epsilon + \boldsymbol O(\epsilon^2),
\end{multline}
\begin{multline}
\paren{\bigmean{\gsdd[a]{\boldsymbol\tau}}}_{33} = \materialdd[s][(0)]{\mean{\tau_{33}}} + (1-a)\bigmean{\pdd{\velocity_s^{(1)}}{x_3^{(1)}} \cdot \boldsymbol\tau_{s3}} + 2a\Div_s\paren{\velocity_s^{(0)}}\mean{\tau_{33}} \\
+ \brackets*{\bigmean{\paren{\half[(1+a)\height^{(1)}]\nabla_s\Div_s\velocity_s^{(0)} + (1-a)\pdd{\velocity_s^{(2)}}{x_3^{(1)}}} \cdot \boldsymbol\tau_{s3}} + 2a\mean{\Div_s\paren{\velocity_s^{(1)}}\tau_{33}}}\epsilon + O(\epsilon^2),
\end{multline}
\end{subequations}
\end{theo}

\begin{proof}
By noticing that we just need to develop the expression of the tensor product $(W(\velocity) + a D(\velocity)) \cdot \boldsymbol\tau$ and to use the decomposition
\begin{equation*}
2\sym(\boldsymbol \xi \cdot \boldsymbol \zeta) =
\begin{pmatrix}
2\sym(\boldsymbol\xi_{ss}\cdot\boldsymbol\zeta_{ss} + \boldsymbol\xi_{s3}\otimes\boldsymbol\zeta_{3s})\ &\ \boldsymbol\xi_{ss}\cdot\boldsymbol\zeta_{s3} + \boldsymbol\xi_{3s}\cdot\boldsymbol\zeta_{ss} + \boldsymbol\xi_{s3}\zeta_{33} + \xi_{33}\boldsymbol\zeta_{3s} \\
" & 2(\boldsymbol\xi_{3s}\cdot\boldsymbol\zeta_{s3} + \xi_{33}\zeta_{33})
\end{pmatrix},
\end{equation*}
we can return to the proof of the previous theorem and conclude.
\end{proof}

%\end{document}
\section{Examples}\label{sec-examples}

\subsection{From Navier--Stokes equations to viscous shallow water equations}

\subsubsection{Problem statement}

The constitutive equation of an incompressible Newtonian fluid write:
\begin{equation}
\stress = 2\eta D(\velocity) - p \boldsymbol{I}, \label{eq-examples-ns-stress}
\end{equation}
where~$\stress$ is the symmetric Cauchy stress tensor and~$p$ a contribution to the pressure. The constant material parameter is the viscosity~\mbox{$\eta > 0$}. These constitutive equations are coupled with standard conservation of mass~\eqref{eq-examples-ns-mass} and momentum~\eqref{eq-examples-ns-momentum} equations. Assuming that the density~$\rho$ is constant, i.e. the fluid is incompressible, the conservation equations write:
\begin{subequations}
\begin{empheq}{align}
\Div\,\velocity &= 0, \label{eq-examples-ns-mass}
\\
\rho\materialdd{\velocity} - \Div\,\stress &= \rho\boldsymbol g, \label{eq-examples-ns-momentum}
\end{empheq}
\end{subequations}
where~$\boldsymbol{g} = - g\boldsymbol{\hat{u}}$ is the constant gravity vector, with $\boldsymbol{\hat{u}}$ some given unit vector. The previous set of equations is closed by initial and boundary conditions. For the latter, we use, on the free surface, the conditions~\eqref{eq-setting-bc-h} with no surface tension ($\gamma = 0$), and, on the substrate, the conditions~\eqref{eq-setting-bc-0} with no friction ($k = 0$).

\subsubsection{Asymptotic analysis}

With $\Sigma \coloneqq \eta U / L$ as characteristic stress, the dimensionless problem writes:

\begin{enonce}{Problem}\label{pb-example-NS}
Find $\velocity$ defined in $\Omega$ such that
\begin{subequations}
\begin{empheq}{align}
\Reop \materialdd{\velocity} - \Div\stress &= -\dfrac{\Reop}{\Fr^2}\boldsymbol{\hat{u}}, \\
\Div\,\velocity &= 0,
\end{empheq}
\end{subequations}
with $\stress = 2D(\velocity) - p\boldsymbol I$, where\/ $\Reop \coloneqq \rho U L / \eta$ is the Reynolds number and $\Fr \coloneqq U / \sqrt{g L}$ is the Froude number.
\end{enonce}

The same notation is used for dimensioned and dimensionless field for convenience.

We would like to use the framework developed in \cref{sec-asymptotic-diffusion} to get a shallow formulation of the above Navier--Stokes problem. To do so, $\mathcal{A}(\unknown)$ is identified with $\Reop\textmaterialdd{\velocity[]_i} + (\Reop/\Fr^2)\hat{u}_i$ and $\flux$ with $\paren{\stress_{ij}}_{1 \leqslant j \leqslant 3}$, for $i$ successively taken in $\{1,2,3\}$. Then, the condition of motion by slices is satisfied, up to order 1 at least, for each component of the velocity, thanks to~\cref{cor-motion-by-slices}. As a consequence, the rate of deformation expands with respect to $\epsilon$ as
\begin{equation}
2D(\velocity) =
\begin{pmatrix}
2D_s(\velocity_s^{(0)}) & \boldsymbol 0 \\
\boldsymbol 0 & -2\Div_s\velocity_s^{(0)}
\end{pmatrix} +
\begin{pmatrix}
2D_s(\velocity_s^{(1)}) & \pdd{\velocity_s^{(2)}}{x_3^{(1)}} - \nabla_s\Div_s\velocity_s^{(0)}x_3^{(1)} \\
\pdd{\velocity_s^{(2)}}{x_3^{(1)}} - \nabla_s\Div_s\velocity_s^{{(0)}^\top}x_3^{(1)} & -2\Div_s\velocity_s^{(1)}
\end{pmatrix}\epsilon + \boldsymbol O(\epsilon^2), \label{eq-examples-Dv-expansion}
\end{equation}
where we used~\cref{prop-asymptotic-mass-v3}.

Following the literature, we make the fundamental assumption that $\Reop / \Fr^2 = O(1/\epsilon)$. Therefore, due to the form of the deformation rate tensor~\eqref{eq-examples-Dv-expansion}, applying the~\cref{cor-asymptotic-diffusion-pressure} leads to
\begin{equation}
p^{(-1)} = 0,\; \mean{p^{(0)}} = -2\Div_s(\velocity_s^{(0)}) + \half[\hat{u}_3]\dfrac{\epsilon \Reop}{\Fr^2}\height^{(1)}.
\end{equation}
We can now apply~\cref{th-asymptotic-diffusion},
\cref{prop-asymptotic-time-derivative-material} and
\cref{cor-asymptotic-mass-height} to get a shallow formulation of the above Navier--Stokes problem.

\begin{enonce}{Problem}\label{pb-example-NS-asymptotic}
Find $\velocity_s^{(0)}$ and $\height^{(1)}$ defined in $\Omega$ such that
\begin{subequations}\label{eq-example-NS-asymptotic}
\begin{empheq}{align}
&\pddt\height^{(1)} + \Div_s\paren{\height^{(1)}\velocity^{(0)}} = 0, \\
&\Reop\height^{(1)}\materialdd[s][(0)]{\velocity_s^{(0)}} - \Div_s\paren{\height^{(1)}\mean{\stress_{ss}^{(0)}}} = -\dfrac{\epsilon \Reop}{\Fr^2}\height^{(1)}\boldsymbol{\hat{u}}_s, \\
&\mean{\stress_{ss}^{(0)}} = 2D_s(\velocity_s^{(0)}) + 2\Div_s(\velocity_s^{(0)})\boldsymbol I_{ss} - \half[\hat{u}_3]\dfrac{\epsilon \Reop}{\Fr^2}\height^{(1)}\boldsymbol I_{ss}.
\end{empheq}
\end{subequations}
\end{enonce}

This reduced model was first obtained in 2006 by Marche~\cite{Mar-2006-sv} after substantial formal calculations. In comparison, our approach requires less than one page.
Solving \cref{pb-example-NS-asymptotic} instead of \cref{pb-example-NS} leads to an error of order $O(\epsilon)$ on velocity and height (cf.~\cref{th-asymptotic-diffusion} or end of \cref{sec-dt-incompr}).

\subsection{Maxwell viscoelastic model}

We use the constitutive equations of a Maxwell viscoelastic fluid~\eqref{eq-maxwell} along with initial and boundary conditions. For the latter, we use, on the free surface, the conditions~\eqref{eq-setting-bc-h} with no surface tension ($\gamma = 0$), and, on the substrate, the conditions~\eqref{eq-setting-bc-0} with no friction ($k = 0$).

With $\Sigma \coloneqq \eta U / L$ again as characteristic stress, the dimensionless problem writes:

\begin{enonce}{Problem}\label{pb-example-Maxwell}
Find $\velocity$ and $\boldsymbol{\tau}$ defined in $\Omega$ such that
\begin{subequations}
\begin{empheq}{align}
\Reop \materialdd{\velocity} - \Div\stress &= -\dfrac{\Reop}{\Fr^2}\boldsymbol{\hat{u}}, \\
\Div\,\velocity &= 0, \\
\We\gsdd[1]{\boldsymbol{\tau}} + \boldsymbol{\tau} &= 2D(\velocity), \label{eq-maxwell-evolution-adim}
\end{empheq}
\end{subequations}
with $\stress = \boldsymbol{\tau} - p\boldsymbol I$, where $\We \coloneqq \lambda V / L$ is the Weissenberg number.
\end{enonce}

The same notation is used for dimensioned and dimensionless field for convenience.

\subsubsection{Sufficient conditions for the motion by slices}

With the help of our formalism, we can derive a new sufficient condition for the motion by slices in the context of the Maxwell model:
\begin{coro}[motion by slices for the Maxwell model]\label{cor-asymptotic-gsdd-motion-slices}
Let $\boldsymbol\tau$ be the elastic stress defined in the Maxwell model~\eqref{eq-maxwell}, and assume that it expands as
\mbox{$\boldsymbol\tau = \boldsymbol\tau^{(0)} +
\boldsymbol\tau^{(1)}\epsilon + \boldsymbol O(\epsilon^2)$} and that the Weissenberg number $\We$ is of order $O(1)$, then
\begin{equation*}
\pdd{\ubold_s^{(0)}}{x_3^{(1)}} =
\pdd{\ubold_s^{(1)}}{x_3^{(1)}} = \boldsymbol 0.
\end{equation*}
\end{coro}
\begin{proof}
Let us write the components $s3$ and $33$ of the equation~\eqref{eq-maxwell-evolution-adim}:
\begin{subequations}
\begin{empheq}{align}
\We\paren{\materialdd{\boldsymbol\tau_{s3}} - \tau_{33}\pdd{\velocity_s}{x_3} - \paren{\nabla_s\velocity_s + \pdd{\velocity[]_3}{x_3}\boldsymbol I} \cdot \boldsymbol\tau_{s3} - \boldsymbol\tau_{ss} \cdot \nabla_s\velocity[]_3} + \boldsymbol\tau_{s3} &= \pdd{\velocity_s}{x_3} + \nabla_s\velocity[]_3, \label{eq-asymptotic-gsdd-motion-slices-s3}
\\
\We\paren{\materialdd{\tau_{33}} - 2\nabla_s\velocity[]_3 \cdot \boldsymbol\tau_{s3} - 2\pdd{\velocity[]_3}{x_3}\tau_{33}} + \tau_{33} &= 2\pdd{\velocity[]_3}{x_3}. \label{eq-asymptotic-gsdd-motion-slices-33}
\end{empheq} By assumption, we know the only term of order $O(1/\epsilon)$ in~\eqref{eq-asymptotic-gsdd-motion-slices-s3} satisfies
\begin{equation}
-\We\tau_{33}^{(0)}\pdd{\velocity_s^{(0)}}{x_3^{(1)}} = \pdd{\velocity_s^{(0)}}{x_3^{(1)}}. \label{eq-asymptotic-gsdd-motion-slices-s3--1}
\end{equation}
\end{subequations}
Now, $\We = 0$ leads to the result. Assume the opposite and that $\tau_{33}^{(0)} = -1/\We$. Then, from equation~\eqref{eq-asymptotic-gsdd-motion-slices-33}, we have
\begin{equation}
\tau_{33}^{(0)} = -\dfrac{1}{\We} = 0
\end{equation}
since $- 2\We\tau_{33}^{(0)} = 2$ and $2\nabla_s\velocity[]_3 \cdot \boldsymbol\tau_{s3}$ is at most of order $O(\epsilon)$. This leads to a contradiction: $\tau_{33}^{(0)}$ cannot be uniformly constant equal to $-1 / \We$. We finally conclude from~\eqref{eq-asymptotic-gsdd-motion-slices-s3--1} that
\mbox{$\textpdd{\velocity_s^{(0)}}{x_3^{(1)}} = \boldsymbol 0$}. Note that if $\tau_{33}^{(0)} = 0$ is assumed, then the proof is still valid.

Hypothesis $\mathcal{H}^{(-1)}$ is satisfied then, from~\cref{lem-asymptotic-diffusion-substrate-flux}, there holds $\boldsymbol\sigma_{s3}^{(0)} = \boldsymbol 0$, which leads to $\boldsymbol\tau_{s3}^{(0)} = \boldsymbol 0$, and then, from~\eqref{eq-asymptotic-gsdd-motion-slices-s3}
\begin{equation}
\paren{\We\tau_{33}^{(0)} + 1}\pdd{\velocity_s^{(1)}}{x_3^{(1)}} = \boldsymbol 0,
\end{equation}
which leads to the second result, since $\tau_{33}^{(0)}$ cannot be uniformly constant equal to $-1 / \We$.
\end{proof}

\subsubsection{Asymptotic analysis}

In what follows, we assume that $\boldsymbol{\tau}$ asymptotic expansion is of the form $\boldsymbol{\tau} = \boldsymbol\tau^{(0)} + \boldsymbol\tau^{(1)}\epsilon + \boldsymbol\tau^{(2)}\epsilon^2 + \boldsymbol O(\epsilon^3)$, whereas $\Reop/\Fr^2$ is of order $O(1/\epsilon)$, as in the previous example. \cref{cor-asymptotic-gsdd-motion-slices} ensures that the property of motion by slices is satisfied.


\paragraph{Momentum equation -- Asymptotic.}\ 
The asymptotic analysis of the momentum equation has already been done in the previous example. Similar to what was done to deduce the pressure in the previous example, the zeroth order of the averaged pressure is given by
\mbox{$\mean{p^{(0)}} = \mean{\tau_{{33}}^{(0)}} + \hat{u}_3\epsilon \Reop/(2\Fr^2)h^{(1)}$}, whence
\begin{equation}
\mean{\stress^{(0)}_{ss}} =
\mean{\boldsymbol{\tau}_{{ss}}^{(0)}} - \mean{\tau_{{33}}^{(0)}}\boldsymbol I_{ss} - \half[\hat{u}_3]\dfrac{\epsilon \Reop}{\Fr^2}\height^{(1)}\boldsymbol I_{ss}.
\end{equation}
In addition, assumption~$\mathcal{H}^{(-1)}$ (see~\cref{lem-asymptotic-diffusion-substrate-flux}) is satisfied for both planar velocity components, then
\begin{equation}
\stress^{(0)}_{s3} = \boldsymbol{\tau}_{{s3}}^{(0)} = \boldsymbol 0. \label{eq-example-maxwell-stress-s3-0}
\end{equation}

\paragraph{Maxwell equation -- Asymptotic.} Taking into account the very last result~\eqref{eq-example-maxwell-stress-s3-0} and applying~\cref{thm-asymptotic-mean-gsdd-tens} and~\cref{prop-asymptotic-mass-v3} lead to the shallow formulation of the Maxwell equation~\eqref{eq-maxwell}:
\begin{subequations}
\begin{align}
\We\gsdd[1][(0)]{\mean{\boldsymbol{\tau}_{{ss}}^{(0)}}} + \mean{\boldsymbol{\tau}_{{ss}}^{(0)}} &= 2 D_s(\velocity_s^{(0)}), \\
\We\paren{\materialdd[s][(0)]{\mean{\tau_{{33}}^{(0)}}} + 2\Div_s(\velocity_s^{(0)})\mean{\tau_{{33}}^{(0)}}} + \mean{\tau_{{33}}^{(0)}} &= -2 \Div_s(\velocity_s^{(0)}). \label{eq-examples-oldroyd-asymptotic-maxwell-33}
\end{align}
\end{subequations}
\begin{rema}
From a microscopic perspective, it is known~\cite[p.~90]{Bird-1987} that the elastic stress can be expressed in terms of the conformation tensor $\boldsymbol c \propto \braket{\boldsymbol\ell \otimes \boldsymbol\ell}$ defined as an average over all possible directions $\boldsymbol\ell$ of molecules constituting the considered fluid:
\begin{equation}
\boldsymbol{\tau} = \dfrac{\eta}{\lambda}\paren{\boldsymbol c - \boldsymbol I}.
\end{equation}
On the reasonable assumption that $\boldsymbol\ell$ expands with respect to $\epsilon$ as $\boldsymbol\ell =
\begin{smallpmatrix} \boldsymbol O(1) \\ O(\epsilon)
\end{smallpmatrix}$, we deduce a possible expansion for the conformation tensor as $\boldsymbol c =
\begin{smallpmatrix} \boldsymbol O(1) & \boldsymbol O(\epsilon) \\ \boldsymbol O(\epsilon) & O(\epsilon^2)
\end{smallpmatrix}$.
Thus, at zeroth order, the component $33$ of the elastic stress reads $\tau_{{33}}^{(0)} = -1/\We$, but it would lead to the contradiction $-1/\We = 0$, according to equation~\eqref{eq-examples-oldroyd-asymptotic-maxwell-33}. If this expansion does not work for $\boldsymbol c$, does it at least work for $\boldsymbol{\tau}$? Assume $\tau_{{33}}^{(0)} = 0$, then $\Div_s(\velocity_s^{(0)}) = 0$ from equation~\eqref{eq-examples-oldroyd-asymptotic-maxwell-33}. If this does not lead to a contradiction, it implies a strong physical constraint, which does not hold in most situations. We therefore understand that the order with respect to $\epsilon$ of a given quantity (here $\boldsymbol{\tau}$) must be componentwise at least of the order of the second member (here $D(\velocity)$).
\end{rema}

\paragraph{Summary.} The shallow formulation of the Maxwell problem writes

\begin{enonce}{Problem}\label{pb-example-Maxwell-asymptotic}
Find $\velocity_s^{(0)}$, $\height^{(1)}$, $\mean{\boldsymbol{\tau}_{{ss}}^{(0)}}$ and $\mean{\tau_{{33}}^{(0)}}$ defined in $\Omega$ such that
\begin{subequations}
\begin{empheq}{align}
&\pddt\height^{(1)} + \Div_s\paren{\height^{(1)}\velocity^{(0)}} = 0, \\
&\Reop\height^{(1)}\materialdd[s][(0)]{\velocity_s^{(0)}} - \Div_s\paren{\height^{(1)}\mean{\stress_{ss}^{(0)}}} = -\dfrac{\epsilon \Reop}{\Fr^2}\height^{(1)}\boldsymbol{\hat{u}}_s, \\
&\mean{\stress^{(0)}_{ss}} =
\mean{\boldsymbol{\tau}_{{ss}}^{(0)}} - \mean{\tau_{{33}}^{(0)}}\boldsymbol I_{ss} - \half[\hat{u}_3]\dfrac{\epsilon \Reop}{\Fr^2}\height^{(1)}\boldsymbol I_{ss}, \\
&\We\gsdd[1][(0)]{\mean{\boldsymbol{\tau}_{{ss}}^{(0)}}} + \mean{\boldsymbol{\tau}_{{ss}}^{(0)}} = 2 D_s(\velocity_s^{(0)}), \\
&\We\paren{\materialdd[s][(0)]{\mean{\tau_{{33}}^{(0)}}} + 2\Div_s(\velocity_s^{(0)})\mean{\tau_{{33}}^{(0)}}} + \mean{\tau_{{33}}^{(0)}} = -2 \Div_s(\velocity_s^{(0)}).
\end{empheq}
\end{subequations}
\end{enonce}

It is consistent with what has been derived in~\cite{BouBoy-2013}. When $\We = 0$, it reduces to the same asymptotic~\eqref{eq-example-NS-asymptotic} as the Navier--Stokes equation.
Solving \cref{pb-example-Maxwell-asymptotic} instead of \cref{pb-example-Maxwell} leads to an error of order $O(\epsilon)$ on velocity, height and all other averaged quantities (cf.~\cref{th-asymptotic-diffusion} or end of \cref{sec-dt-incompr}).

\section{Conclusion}

A new framework for asymptotic analysis has been presented in this paper. We have proposed a set of tools for performing thin-layer approximations of free-surface fluid models with flat topography and no-penetration boundary condition. The motion by slices is no more an arbitrary assumption, and minimal assumptions for obtaining it have been pointed out. Moreover, for first time, to the best of our knowledge, the first-order expansion of the Gordon--Schowalter derivatives has been developed. The proposed framework makes it possible to clarify the necessary assumptions and to simplify the asymptotic analysis of a wide variety of partial differential equations from continuum physics, whether they are scalar, vector or tensor, with objective derivative and/or diffusive, reactive or advective terms. This framework has been illustrated by two applications: the Navier--Stokes and the Maxwell models. In perspective, a generalisation of this work would be to consider here an arbitrary topography, following the work made in~\cite{BouWes-2004} for the Saint-Venant system.

\appendix
\section{Leibniz formulas}\label{sec-appendix}

\begin{subequations}
\begin{prop}[Leibniz's Integral Rule]
Let $\mathcal{X}$ and $I$ be two intervals of $\R$ and \mbox{$f \colon (x, t) \in \mathcal{X} \times I \mapsto f(x, t) \in \R^\spdim$}, $\spdim \in \mathbb{N}\backslash\{0\}$, a sufficiently smooth function. Let also $a$ and $b$ two differentiable functions defined from $\mathcal{X}$ to $I$. Then the parameterized integral function $F$, defined in $\mathcal{X}$ by
\begin{equation*}
F(x) = \int_{a(x)}^{b(x)}f(x, t)\dt,
\end{equation*}
is differentiable and
\begin{equation*}
F'(x) = \int_{a(x)}^{b(x)}\pdd{f}{x}\dt + f(x, b(x))b'(x) - f(x, a(x))a'(x).
\end{equation*}
\end{prop}

\begin{coro}[Leibniz's Integral Rule for depth average]\label{cor-asymptotic-leibniz-L}
Let $\mathcal{L}$ be a first order differential linear operator, depending only on partial derivatives $\partial_t$, $\partial_{x_1}$ and $\partial_{x_2}$. For any sufficiently smooth scalar field $\varphi$, the following formula is satisfied:
\begin{equation}
\mathcal{L}(\height\mean{\varphi}) = \height\mean{\mathcal{L}\varphi} + \varphi(x_3 = \height)\mathcal{L}\height.
\end{equation}
\end{coro}

\begin{coro}\label{cor-asymptotic-leibniz-nabla}
For any sufficiently smooth vector field $\ubold$, the following formulas are satisfied:
\begin{align}
\Div_s(\height\mean{\ubold}) &= \height\mean{\Div_s\ubold} + \ubold_s(x_3 = \height)\cdot\nabla_s\height, \\
\Delta_s(\height\mean{\ubold}) &= \height\mean{\Delta_s\ubold} + 2\nabla_s\ubold_s(x_3 = \height)\cdot\nabla_s\height + \ubold_s(x_3 = \height)\Delta_s\height, \\
\height\mean{\Delta_s\ubold} &= \Div_s(\height\mean{\nabla_s\ubold}) - \nabla_s\ubold_s(x_3 = \height)\cdot\nabla_s\height.
\end{align}
\end{coro}

We close this section by giving a last relation that will be useful:
\begin{equation}\label{eq-asymptotic-leibniz-kpp}
\mean{x_3\partial_{x_3}\varphi} + \mean{\varphi} = \varphi(x_3 = \height).
\end{equation}
\end{subequations}


\subsection*{Acknowledgements} 

This work has taken place in the context of an interdisciplinary collaboration with two physicists, Hélène Delanoë-Ayari (Institut Lumière Matière, Univ. Lyon 1) and François Graner (Laboratoire Matière et Systèmes Complexes, Univ. Paris 7), on the modelling of biological tissues. We warmly thank them for fruitful discussions. This work is supported by the CNRS interdisciplinary mission (MITI).

\bibliographystyle{crplain}
\bibliography{CRMATH_Shourick_20230365}
\end{document}