\input{aefm-macro}
2-4spreadmud.tex 
\setcounter{chapter}{2}
\setcounter{section}{3}


\section{Spreading of a shallow mass on an incline}
\setcounter{equation}{0}

{\bf References:} \\
C. C. Mei, (1966), Nonlinear gravity waves in a thin sheet of viscous fluid, {\it J. Math. \& Phys.} 45, 482-496. \\
K. F. Liu \& C. C. Mei, 1989, Slow spreading of a sheet of  a Bingham fluid on an inclined plane, {\it J. Fluid Mech. } 207, 505-529.\\
X. Huang \& M. H. Garcia,1997, Asymnptotic solution for Bingham debris flows, {\it Debris -flow Hazards and Mitigation,} ASCE Proc.  pp.561-575.\\
 
The spreading of a finite mass of paint, paper pulp, mud or lava  on an inclined plane is of interest to a variety of industrial and  geological problems. For mud modeled as a Bingham-plastic non-Newtonian fluid, Liu \& Mei (1989) solved the equation similar to (2.3.22)  numerically. An analytical solution for Herschel-Bulkley fluid was given later by Huang \& Garcia (1997). We modify  their theories  for the simpler case of Newtonion fluid here. 


 The free surface of a thin layer is expected to flatten in time 
over most of the profile,  but it should steepen near the downstream front where the spatial derivative is much more
important   than elsewhere. Let us study the   problem by dividing the total fluid extent into two:  the far field not too
close to the steep 
front  and the near field around the front.  
\subsection{Far field away from the front}
Let the  total length of the thin layer be $L$ and the maximum depth be $H$, with $H/L\ll 1$. We  introduce the normalized variables suitable for the far field. 
\be h=Hh', \,\, x=Lx',\,\, t=Tt' \ee
where $T$ will be chosen to simplify the appearance of the final equation. Thus,
\[\f{H}{T} \f{\p h'}{\p t'} + \f{\rho g\sin \theta }{3 \mu}\f{H^3}{L} \f{\p}{\p x'} \lb \lp 1-\f{H}{L \tan\theta }\f{\p h'}{\p x'} \rp h'^3\rb= 0 \] 
Let us choose the time time scale to be
\be T=\f{L}{ {\rho g H^2\sin \theta}/{3\mu}}\label{timescale}\ee
where the denominator is the average fluid speed across the depth. The dimensionless equation for the farfield   reads,
  \be  \f{\p h'}{\p t'} +  \f{\p}{\p x'} \lb \lp 1-\f{H}{L \tan\theta }\f{\p h'}{\p x'} \rp h'^3\rb= 0 \label{norm-basic}\ee

Assume  \be \f{H}{L \tan\theta }\equiv \ep \ll 1\label{assumption}\ee
(\ref{norm-basic}) can then  be approximated  by the   hyperbolic equation,
\be \f{\p h'}{\p t'} +  \f{\p  h'^3}{\p x'}
=\f{\p h'}{\p
t'} +  3 h'^2   \f{\p  h' }{\p x'}  =0
\label{far}\ee
 It is of the class  called kinematic wave equation in
flood  hydrology.   Solution can be obtained by the theory of characteristics. 

In particular 
(\ref{far}) can be rewritten in the form
\be  \f{\p h}{\p
t}dt +    \f{\p  h }{\p x}dx  =dh=0
\label{char-form}\ee
if 
\be \f{dx'}{dt'} =  3h'^2 \label{char-curve}\ee 
 The  last two equations imply that $h'$ remains constant for all $t'$ along the characteristic curve $x'(t')$ defined by the differential equation (\ref{char-curve}).  
Moreover, if the initial profile is presecribed,
\be  h'(x', 0) = h'_o(x'),   \ee
then all characteristics are straight but of different slopes. 
The characteristic  originated from the initial point $x'=\xi'$ at $t'=0$ is the straight line
 \be  x(t,\xi) =  3h_o^2(\xi) t + \xi, \label{char}\ee
along which 
\be h'(x',t') = h'_o(\xi') \label{formalsol}\ee
 
In principle we can solve for $\xi' $ from (\ref{char}) in term of $x',t'$ and substitute the result into (\ref{formalsol}) to get $h'(x',t')$.  For any initial  hump, the characteristics at the front must intersect one another, impling multivaluedness of solution, which is phycially unacceptable. The solution can still be proceeded if a discontinuity, {\it shock},  is allowed at $x'=x'_s(t')$, as long as mass is conserved. 

As an example, consider a triangular  initial profile:
\be h'_o (x') = \lcb \begin{array} {cc}   x' & 0<x'<1\\
0, & x'<0, ~~ \& ~~ x'>1. \end{array} \right. \ee
The profile  has a shock front to begin with, which is  a mathematical idealization of course. 
From (\ref{char}) we get
\be  x' =  \lcb\begin{array}{cc} 3\xi'^2 t' +
\xi', & 0<\xi' <1,\\
\xi' , & \xi'<0, \xi' >1\end{array} \right.\label{char-ex} \ee
Within $0<\xi'<1$,   $\xi'$ can be solved in terms of $x',t'$, 
\be \xi' = \f{1}{6t'} \lp -1+\sqrt{1+12 x't'} \rp\ee
For any $t'>0$, 
\be h' = 0, \quad \mbox{ for}\quad   x' <0, ~~\& ~~x'>x'_s(t'), \ee
and \be h'=\xi'=\f{1}{6t'} \lp -1+\sqrt{1+12 x't'}\rp, \quad    \quad 0<x'<x'_s(t).\label{farfieldsolution}\ee
The shock front $x'_s(t')$   is unknown. 

 To locate the
shock front
$x'_s(t')$, we invoke  mass conservation,
\be  \mbox{initial volume}= \f{1}{2}=  \int_0^{x'_s(t)} h'\, dx' = \int_0^{h'_s}
h'\left.\f{dx'}{dh'}\right|_{t'}dh' \ee
From (\ref{char-ex}), 
\be x' = 3 h'^2 t' + h'\label{char-h}\ee
 so that \[ \left.\f{dx'}{dh'}\right|_{t'} = 6h't' +1\]
Thus
 \[ \f{1}{2} = \int_0^{h'_s} h'(6h't' +1) dh' = \f{h'^2_s}{2}(4h'_s t +1) \]
or \be t' = \f{1-h'^2_s}{4h'^3_s}\label{timeshock}\ee
This gives the shock height $h'_s$ implicitly as a function of $t'$.  To get the shock postion $x'_s(t')$ we put this result into  (\ref{char-h})
 \be x'_s= \f{3+h'^2_s}{4h'_s}\label{shockfront}\ee
Together (\ref{timeshock}) and (\ref{shockfront}) give  $x'_s(t')$, as plotted  in figure (\ref{fig:shockfront}).
\begin{figure}
\label{fig:shockfront}\begin{center}
\includegraphics[scale=0.5]{xs-hs-t.eps}\end{center}
\caption{Time variation of the location $x'_s(t') $ and depth $h'_s(t')$ of shock front.}
\end{figure}
To be needed later we note that the shock speed is given by 
\be 
\f{d x'_s}{dt'}=\f{d x'_s}{dh'}/\f{dt'}{dh'}=\f{1-3{h'_s}^{-2}}{4}\f{4{h'_s}^2}{1-3{h'_s}^{-2}}={h'_s}^2. \label{shockspeed}
\ee
 
%%%%%%%%%%%%%%%%

\subsection{Near field of downstream front} 
Clearly near the shock the local free surface slope cannot satisfy (\ref{assumption}). We define the longitudinal length scale for the near field  to be $L_s$ so that \be  \f{H}{L_s\tan\theta} = O(1), \quad \mbox{i.e.,} \quad\f{L_s}{L} =O(\ep)\ee
and renormalize the near field,
\be h=HH', ~~ t=Tt',~~\f{x-x_s(t)}{L_s}= X', ~~\mbox{or}\quad x'=x'_s+\ep X'\ee
  Thus the near field solution $H'$ is a function of  $X'(x',t') $ and $t'$. The derivatives are now related by the chain rule,
\be \f{\p H'(x',t')}{\p t'}  =    \f{\p H'}{\p t'} +\f{\p H'}{\p X'} \f{\p X'}{\p t'}    =  \f{\p H'}{\p t'}  -\f{1}{\ep}\f{d x'_s}{d t'} \f{\p H'}{\p X'}  \ee
\be \f{\p H'(x',t')}{\p x'} =  \f{\p H'}{\p X'}\f{\p X'}{\p x'} = \f{1}{\ep} \f{\p H'}{\p X'} \ee
 Eq. (\ref{norm-basic}) becomes
\be \f{\p H'}{\p t'} -\f{1}{\ep}\f{dx'_s}{dt'} \f{\p H'}{\p X'} + \f{1}{\ep}\f{\p }{\p X} \lb H'^3\lp 1-\f{\p H'}{\p X'} \rp \rb = 0 \ee
To the leading order  the near field is approximately governed by \be -\f{dx'_s}{dt'} \f{\p H'}{\p X'} + \f{\p }{\p X'} \lb H'^3\lp 1-\f{\p H'}{\p X'} \rp \rb = 0 \label{statwave}\ee
 which is the ordinary differential equation for a stationary wave, moving at the known shock speed ${d x'_s}/{dt'}={h'_s}^2$. Integrating with respect to $X'$ once, we get
\[ -{h'_s}^2H'+H'^3\lp 1-\f{d H'}{d X'}\rp = 0. \]
The constant of integration is taken to be zero so that the surface touches the dry bed where $H'=0$.  We may rewrite (\ref{statwave}) as 
\[ dX' = -\f{H'^2 dH'}{{h'_s}^2-H'^2} = dH' \lb 1-\f{h'_s}{2}\lp\f{1}{h'_s-H'} +\f{1}{h'_s+H'}\rp \rb \]
It follows upon   integration  that 
\be H'+\f{h'_s}{2} \log \lp\f{h'_s-H'}{h'_s+H'}\rp = X' -X'_* \label{nearfieldsolution}\ee
This is a smooth surface decreasing monotonically from $H' = h'_s$ at $X'=-\infty$   to $H'=0$ (dry bed) at  $X'=X'_*$. To determine the value of $X_*$, 
we   require that mass under the smooth profile be the same as that under the profile with a shock, as sketched in Figure \ref{fig:massconserv}. 
\begin{figure}
\label{fig:massconserv}\begin{center}
\includegraphics[scale=0.75]{matching.eps}\end{center}
\caption{Matching the near-field and far-field at the front by preserving mass.}
\end{figure}


\be \int_{-\infty}^0\lp h'_s - H' \rp dX' = \int_0^{X'_*} H'(X')  dX'  \label{area}\ee

  
Referring to Figure \ref{fig-massconserv}, the integrals can be equivalently carried out as 
  \be -\int_{H'(0)}^{h'_s}X'(H') dH' = \int_0^{H'(0)}X'(H')  dH' \label{matchingcondition} \ee
where $H'(0)$ is the height at $X'=0$, given by
\be H'(0)+\f{h'_s}{2} \log \lp\f{h'_s-H'(0)}{h'_s+H'(0)}\rp = -X'_* \label{profileat0}\ee
Substituting (\ref{nearfieldsolution}) into (\ref{matchingcondition}), we get, for the left-hand-side, 
\[  -\int_{H'(0)}^{h'_s} X'(H')dH' = X'_*\lp H'(0) -h'_s\rp +\f{1}{2}\lp {H'(0)}^2-{h'_s}^2 \rp 
-\f{h'_s}{2}\int_{H'(0)}^{h'_s}\log \lp \f{h'_s-H'}{h'_s+H'}\rp dH'\]
and for the right-hand side,
\[ 
 -\int_0^{H'(0)}X'(H')dH' = X'_*\lp H'(0) -0 \rp +\f{1}{2}\lp {H'(0)}^2-0 \rp 
-\f{h'_s}{2}\int^{H'(0)}_0\log \lp \f{h'_s-H'}{h'_s+H'}\rp dH'\]
Therefore,
\[ X'_s=-\f{1}{2}\lb \int_0^{h'_s}\log  \lp \f{h'_s-H'}{h'_s+H'}\rp dH' + h'_s \rb = -\f{h'_s}{2}\lb \int_0^1\log \lp \f{1-z}{1+z}\rp dz +1\rb\]
Using the facts
\[ \int_0^1\log (1-z)dz = -1, ~~~\int_0^1\log (1+z)dz = 2\log 2 -1\]
we get 
\be X'_*= \lp \f{1}{2}-\log 2\rp h'_s=0.19897 h'_s \label{mudfront}\ee
The near field solution is now complete.

Recall that the near and far fields coordinates  are related by
\be  x'=x'_s+\ep X'\label{nearandfar}\ee
in general and 
\be x'_*(t') = x'_s(t') + \ep \lp \f{1}{2}-\log 2\rp h'_s(t')\label{mudfrontfar}\ee in particular. 

We can get the uniformly valid  solution ($h_{U.V.}$) by adding the near and far field solutions and subtracting the common part. Thus, in terms of the far-field variables,
\be  h'_{U.V.} =h'+H'-h'_s\label{uniformsolution}\ee
where $h'$ is given by (\ref{farfieldsolution}) and $H'$  is given by (\ref{nearfieldsolution}), while the coordinates are related by (\ref{nearandfar}). 


Sample profiles of the free surface are shown for $ \ep = 0.1$  in Figure \ref{fig:UVsolution} for various instants.
\begin{figure}\begin{center}
\includegraphics[scale=0.75]{proevo.eps}\end{center}
\caption{Slow evolution of the free surface of a vicous fluid down a  dry incline. The  initial profile is triangular   with $\ep = 0.1$.  }
\label{fig:UVsolution}\end{figure}





 \end{document}