\input{aefm-macro} 2-8-selec-iso.tex\\ \\ %%%%%%%%%%%%%%%%%%%%%%%%%% \setcounter{chapter}{2} %\renewcommand{\baselinestretch}{2} \def\f{\frac} \def\erf{\rm erf} \def\p{\partial} \def\be{\begin{equation}} \def\ee{\end{equation}} \def\a{\alpha} \def\n{\nabla} \def\z{\zeta} \def\la{\langle} \def\ra{\rangle} \def\lb{\left[} \def\rb{\right]} \def\lcb{\left\{} \def\rcb{\right\}} \def\lp{\left(} \def\rp{\right)} \def\rhu{\rightharpoonup} \def\o{\over} \def\ep{\epsilon} \def\gl{\stackrel{>}{\scriptstyle <}} \def\lg{\stackrel{<}{\scriptstyle >}} \def\lap{\nabla^2} \def\no{\noindent} \def\h{\hat{}} \def\hs{\hspace{0.3 cm}} \def\definition{\stackrel{\rm def}{=}} \def\erf{{\rm erf}} \def\erfc{{\rm erfc}} \def\pvi{{\int\mskip -33 mu - \quad}} \def\Real{\mbox{Re}\,} \def\Imag{\mbox{Im}\,} \def\pvint{{\int\!\!\!\!\!\!-}} \def\2int{\int\!\!\!\int} \def\3int{\int\!\!\!\int\!\!\!\int} \def\sgn{\mbox{sgn}} \def\bfq{{\bf q}} \def\bfn{{\bf n}} \def\bfx{{\bf x}} \def\theequation{\thesection.\arabic{equation}} \pagestyle{myheadings} %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \setcounter{section}{7} \section{Selective withdrawal in an isothermal stratifed fluid} \setcounter{equation}{0} [References]: \\ Gelhar \& Mascolo, 1966\\ R.C. Y. Koh, 1966 {\it J. Fluid Mechanics}, {\bf 24}, pp. 555-575.\\ Brooks, N. H., \& Koh, R. C. Y., Selective withdrawal from density stratified reservoirs. {\it J Hydraulics ASCE}, HY4, July 1969. 1369-1400.\\ Ivey, G. N. \\ Monosmith, et. al. \\ Nuclear power plants are often located near a large lake so that cold water can be withdrawn to cool the engines. If the lake is thermally stratified so that there is apprciable temperature gradient vertically. Cold water can be withdrawn from lake bottom, and warm water from the power plant can be returned to the top of lake. Due to the fact that stratification supresses vertical motion, water motion hence thermal mixing should be limited to a vertically thin layer not far away from the intake. In this section we treat the slow and steady flow of an isothermal but stratified fluid, due perhaps to salinity variation, into a two- dimensional line sink. The molecular diffusivity of salt is ignored. These simplifications permits a analytical solution which is easy to examine the physical implicaions. For a slow flow of a viscous fluid flowing in a vertically stratified fluid. we expect, in light of Yih's theorem, that the motion of the fluid should be confined within a thin layer. \subsection{Estimation of scales} Governing equation for continuity : \be u_x+w_z=0 \ee Ignoring mass diffusion, the incompressiblity condition reads \be u \rho_x + w \rho_z = 0 \ee \be \rho \lp uu_x + wu+z \rp = -p_x + \mu \lp u_{xx} + u_{zz} \rp \ee \be \rho \lp uw_x + ww_z \rp = - p_z - \rho g + \mu \lp w_{xx} + w_{zz} \rp \ee Anticipating a thin boundary layer due to vertical suppression, we introduce the following normalization. \be u= Uu', ~w = \f{U\delta}{L} w' , ~ \rho = \rho_0 \rho' , ~x = Lx' , ~ z = \delta z' \ee then \be u'_x + w'_z= 0 \ee \be u'\rho'_x + w'\rho'_z = 0 \ee \be \f{\rho_0U^2}{L} \lb \rho' \lp u'u'_x + w'u'_z \rp \rb = - \f{P}{L} \f{\p p'}{\p x'} + \f {\mu U}{\delta^2} \lp u'_{xx} \f{\delta^2}{L^2} + u'_{zz} \rp \ee \be \f{\rho U^2}{L} \f{\delta}{L} \lb \rho' \lp u'w'_x + w'w'_z \rp \rb = - \f{P}{\delta} \f{\p p'}{\p z'} - \rho_0 g\rho' + \mu \f{U}{\delta^2} \f{\delta}{L} \lb w'_{xx} \f{\delta^2}{L^2} + w'_{zz} \rb \ee Here $\delta $ is the boundary layer thickness, $ L\sim x $ is the horizontal length scale, and $U $ is the horizontal velocity scale which must be of the order $\sim Q/\delta$. After dividing the $x,y$ momentum equations by $\mu U\delta^2$, and noting \be \f{\rho_0 U^2}{L} \,\, \f{\delta^2}{\mu U} = \lp \f{\rho_0 U \delta}{\mu} \rp \f{\delta}{L} \equiv R \f{\delta}{L} \ee we get \be R \f{\delta}{L} \lb \rho' \lp u'u'_x + w'w'_z \rp \rb = - \lp \f{P\delta}{\mu U} \f{\delta}{L} \rp p'_x + \lb u'_{xx} \lp \f{\delta}{L} \rp^2 + u'_{zz} \rb \ee Assume \be R= \f{\rho_0 U \delta}{\mu}=O(1), ~\mbox{and}~ \f{\delta}{L} \ll 1. \ee We must have \be \f{P\delta}{\mu U} \f{\delta}{L} = O(1)\ee so that the viscous stress is forced by the pressure gradient. Thus we can take the pressure scale to be \be P = \mu \f{UL}{\delta^2} \label{press-scale}\ee To the leading order the $x$ momentum equation is simply \be 0 = -p_x +\mu u_{zz} \ee Now the $z$ momentum equation, \be \f{\rho_0 U^2}{L} \f{\delta}{L} \lb \rho' \lp u'w'_x + w'w'_z \rp \rb = - \f{P}{L} \f{1}{\delta /L} \f{\p p'}{\p z'} - \rho_0 g \rho' + \f{\mu U}{\delta^2} \f{\delta}{L} \lb w'_{xx} \f{\delta^2}{L^2} + w'_{zz} \rb . \ee Dividing by $P/L$ and using (\ref{press-scale}), \be R \lp \f{\delta}{L} \rp \lb \rho' \lp w'w'_x + w'w'_z \rp \rb = - \f{1}{\lp \f{\delta}{L} \rp} \f{\p p'}{\p z} - \f{\rho_0 gL}{P} \rho' + \lp \f{\delta}{L} \rp \lb w'_{xx} \lp \f{\delta}{L} \rp^2 + w'_{zz} \rb \ee i.e., \be R \lp \f{\delta}{L} \rp^2 \lb \rho' \lp w'w'_x + w'w'_z \rp \rb = - \f{\p p'}{\p z'} - \f{\rho_0 g \delta}{P} \rho' + \lp \f{\delta}{L} \rp^2 \lb w'_{xx} \lp \f{\delta}{L} \rp^2 + w'_{zz} \rb . \ee Since gravity must be important, we must have \be \f{\rho_0 g \delta}{P} = O(1),\label{buoyancy}\ee The $z$ momentum equation reduces to \be 0= - \f{\p p}{\p z} - \rho g \ee meaning that pressure is hydrostatic. Now (\ref{buoyancy}) implies \be \f{\rho_0 g \delta^3}{\mu U L} =O(1) \ee since $Q = O(U \delta ) $. It follows that \be \f{\rho_0 g}{\mu Q} \f{\delta^4}{L} = O(1) \ee i.e., \be \delta^4 \sim \f{\mu Q L}{\rho_0 g} \ee Since $L \sim x$, we have \be \delta \sim\lp \f{\mu Qx}{\rho_0 g}\rp ^{1/4}\ee Thus by a mere scale estimate, we not only achieved (confirmed) a simplification, but also conclude that the boundary layer thickness increases with $x^{1/4}$. \subsection{Approximate equations} We summarize the the approximate equations stating continuity, \be u_x + w_z= 0 \label{Eq:42.15} \ee incompressibility, \be u\rho_x + w \rho_z=0 \label{Eq:42.16} \ee horizontal momentum balance, \be 0 = -p_x + \mu \, u_{zz} \label{Eq:42.17} \ee and vertical momentum balance \be 0 = -p_z - \rho g . \label{Eq:42.18} \ee Let the stream function $\psi$ be defined such that \be u = \psi_z, \qquad w = -\psi_x \label{streamfunction}\ee Equation (\ref{Eq:42.16}) means that in the direction of the local velocity, the density is constant. Thus the the density can only be a function of the stream function, i.e., \be \rho = \rho (\psi ) \label{Eq:42.density}\ee Let us check that (\ref{Eq:42.density}) implies (\ref{Eq:42.16}). Clearly \[ \rho_x = \f{d\rho}{d\psi} \psi_x, \quad \rho_z = \f{d\rho}{d\psi} \psi_z,\] so that \[ \f{d\rho}{d \psi} = \f{ \rho_x}{ \psi_x}=\f{\rho_z}{\psi_z} \] Eq. (\ref{Eq:42.16}) follows by using (\ref{streamfunction}). Eliminating $p$ between Eqn. (\ref{Eq:42.17}) and Eqn. (\ref{Eq:42.18}) we get \[ \mu \, u_{zzz} + g\rho_x = 0, \] which can be rewritten in terms of the $\psi$, \be \nu \, \psi_{zzzz} + \f{g}{\rho} \, \f{d\rho}{d \psi} \, \psi_x = 0 . \label{Eq:42.19'} \ee Clearly this is a boundary layer equation. For a weakly stratifed fluid, we approximate $\rho$ by a constant in the boundary layer. The boundary conditions are: \be u, \quad z \to \pm \infty. \label{Eq:42.20a} \ee which implies that \[ u_z, \, u_{zz}, u_{zzz}, ... \to 0, \quad z \to \pm \infty \] It follows from Eqn. (\ref{Eq:42.19'}) that \be \psi_x = w \to 0 \qquad z \rightarrow \pm \infty \label{Eq:42.20b} \ee if $d\rho/d\psi \neq 0$ as $z \to\pm \infty$. Integrating Eqn. (\ref{Eq:42.15}) from $z = -\infty$ to $z = \infty$, and taking into account Eqn. (\ref{Eq:42.20b}), we get \[ \f{\p}{\p x} \, \int^\infty_{-\infty} \, u \, dz = 0 . \] Therefore, \be \int^\infty_{-\infty} \, u \, dz = - Q \; \mbox{(const.)} \label{Eq:42.21} \ee where $Q$ denotes the steady discharge. In terms of $\psi$, the boundary conditions Eqns. (\ref{Eq:42.20a}), (\ref{Eq:42.20b}) and the integral constraint Eqn. (\ref{Eq:42.21}) can be expressed as \be \psi_x, \, \psi_z, \, \psi_{zz} \rightarrow 0 \qquad (z \to \pm \infty) \label{Eq:42.22} \ee and \be \psi(x,\infty) - \psi(x,-\infty) = -Q . \label{Eq:42.23} \ee {\bf Remark}: Unlike the 2-D laminar jet problem discussed before, the momentum flux rate is relatively unimportant in the present case of slow flow. \subsection{Similarity solution for linear stratification} In the present case $u<0$ hence $\psi$ increases as $z$ decreases. Therefore $ d \rho/d\psi$ and $d \rho/dz$ have opposite signs. For a stably stratified fluid, $\rho $ increases as $z$ decreases. Hence $d \rho/d \psi >0$. In the special case where \[ \f{g}{\rho}\f{d\rho}{d \psi} = \mbox{const.} \, \equiv \nu \, c > 0 \] (\ref{Eq:42.19'}) becomes linear \be \psi_{zzzz} + c \psi_x = 0. \label{Eq:42.24} \ee Let us look for a one-parameter transformation \[x =\lambda^\alpha x^* ,z =\lambda^\beta z^* , \psi=\lambda^\gamma\psi^* \] such that Eqns. (\ref{Eq:42.22}) - (\ref{Eq:42.24}) are invariant. Since Eqn. (\ref{Eq:42.22}) is homogeneous, it gives no information about the choice of $\alpha, \beta, \gamma$. From Eqns. (\ref{Eq:42.23}) and (\ref{Eq:42.24}), invariance of the boundary value problem requires that \[ \alpha = 4 \beta, \qquad \gamma = 0 . \] Therefore, a similarity solution can be found in the form \be \f{\psi}{Q} = f(\eta), \qquad \eta = \lp \f{c}{x} \rp^{1/4} z \label{Eq:42.25} \ee From Eqn. (\ref{Eq:42.25}), it follows that \be w = -\psi_x = \f{Q}{4} \, \f{\eta}{x} \, f'(\eta) \label{vertical}\ee \be u = \psi_z = \lp \f{c}{x} \rp^{1/4}Q \, f'(\eta), \label{horizontal}\ee \be \psi_{zz} = \lp \f{c}{x} \rp^{1/2}Q \, f''(\eta) \ee \be \psi_{zzzz} = \f{c}{x} Q \, f''''(\eta) \ee Eq(\ref{Eq:42.24}) then becomes \be 4 \, f'''' - \eta \, f' = 0 \label{Eq:42.26} \ee and Eqns. (\ref{Eq:42.22}) and (\ref{Eq:42.23}) become \be \eta \, f', \, f', \, f'' \rightarrow 0 \qquad (\eta \rightarrow \pm \infty) \label{Eq:42.27} \ee \be f(\infty) - f(-\infty) = -1 . \label{Eq:42.28} \ee {\bf Remarks}: (i) Since $u(x, 0)$ must be finite, $f'(0)$ is finite. It follows that \[ w(x,0) = \left. -\f{Q}{4} \, \f{\eta}{x} \, f'(\eta) \right|_{\eta = 0} = 0 . \] suggesting that $w$ is odd in $\eta$ and antisymmetric about $z = 0$, in the $x,z$ plane. Hence $f(\eta)$ is also anti-symmetric in $z$ about $z = 0$, i.e., $f(-\eta) = -f(\eta)$. Eqn. (\ref{Eq:42.28}) may be repleced by \be f(\infty) = -1/2 \label{Eq:42.28'} \ee It is always possible to take $f(0) = 0$, amounting to designating the axis as the streamline $\psi=0$. The anti-symmetry of $w$ also implies symmetry of $u$ in $z$ about $ z= 0$. These properties are the consequence of the assumption that $d\rho/d \psi $ is constant and cannot be expected in general because of gravity. (ii) The similarity transformation can also be derived heuristically by examining Eqns. (\ref{Eq:42.22}) - (\ref{Eq:42.24}). Specifically, Eqn. (\ref{Eq:42.24}) implies that \[ \f{1}{\delta^4} \sim \f{1}{x}, \quad \mbox{thus}\quad \delta \sim x^{1/4} \] where $\delta$ is the boundary layer thickness. Moreover, Eqn. (\ref{Eq:42.21}) gives $u \cdot \delta \sim O(1)$. Therefore, \[ u \sim \f{1}{\delta} \sim x^{1/4}, \quad \mbox{and} \quad {\psi \sim u \delta \sim x^0} . \] These results agree with (\ref{horizontal}). \subsection{Analytical Solution} Let \be g(\eta) = f'(\eta) . \label{Eq:42.29} \ee It follows from Eqn. (\ref{Eq:42.26}) that \be 4 \, g''' - \eta \, g = 0 . \label{Eq:42.30} \ee Furthermore, in view of Eqn. (\ref{Eq:42.27}), $g$ and its derivatives vanish at $\eta = \pm \infty$. Let us apply Fourier transform defined by \[ \hat{g}(k) = \f{1}{2\pi} \int^\infty_{-\infty} \, e^{-ik\eta} \, g(\eta) d\eta \, ; \] Eqn. (\ref{Eq:42.30}) is then converted to \[ -4 \, ik^3\hat{g} - i \, \f{d\hat{g}}{dk} = 0 \] The solution is \[ \hat{g} (k) = A \, e^{-k^4} , \] where $A$ is a coefficient. Taking the inverse Fourier transform of $\hat{g}$, we get \be g(\eta) = A \int^\infty_{-\infty} e^{-k^4 + ik\eta} \, dk = 2A \int^\infty_0 \, \cos \, k\eta \, e^{-k^4} \, dk. \label{Eq:42.31} \ee It is readily seen from Eqn. (\ref{Eq:42.31}) that $g(-\eta) = g(\eta)$, confirming the earlier arguments on the symmetry of $u$ and $f'(\eta)$. In view of Eqns. (\ref{Eq:42.28'}) and (\ref{Eq:42.29}), the constant $A$ must be chosen such that \[ 2A \int^\infty_0 \, d\eta \, \int^\infty_0 \, dk \, \cos \, k\eta \, e^{-k^4} =-\f{1}{2} \] It is shown in the Appendix that the double integral is $-\pi/2$ so that \[ A = -\f{1}{2\pi}. \] In particular, \[ g(0) = 2A \int^\infty_0 \, e^{-k^4} \, dk = \f{-\sqrt{2}}{4\Gamma (3/4)} = -0.2885 . \] The integral in Eqn. (\ref{Eq:42.31}) may be evaluated numerically. Alternatively, Eqn. (\ref{Eq:42.30}) is solved numerically by Chebyshev-polynomial method for unbounded domains (as, for example, J.P. Boyd (1992), {\bf J. Comp. Phys.}, \underline{45}, 43-79). Homework by Prof C. S. Wu, 1.63 class of 1995, now at U of Wisconsin, Madison) From Figure \ref{fig:selwith-yang} ( Prof. T.S. Yang, National Cheng Kung University, Taiwan (1.63 Class of 1995) it is seen that, away from the center line of the sink, there are regions in which $u > 0$. Under the assumption that $c > 0$, we have, $\p \rho/\p y > 0$ so that locally the fluid may be statically unstable. \setcounter{figure}{0} \begin{figure}[h] \begin{center} \includegraphics[scale=0.50]{dsf1} \end{center} {Velocity profile of a 2D slow flow into a line sink in a density-stratified fluid. Calculations by G. D. Lee, 2002} \label{fig:selwith-yang} \end{figure} \subsection*{Appendix: A double integral (Yang, 1995, Homework)} Consider the double integral \be I(n)=\int^\infty_0 \, dk \,e^{-k^n}\int^\infty_0 \, d\eta \, \cos \, k\eta \, \label{eq:doubleinetgral}\ee Note that \[ \int_{-\infty}^\infty \delta(k) e^{ik\eta} \, dk = 1\] By inverse transform we get the Fourier integral representation of the delta function: \[ \delta(k) = \f{1}{2\pi}\int_{-\infty}^\infty dk \,e^{-ik\eta}=\f{1}{\pi}\int_{0}^\infty d\eta \,\cos k\eta \] Now we use this result in (\ref{eq:doubleinetgral}) to get \be I(n) =\pi \int_0^\infty \delta(k) e^{-k^n} dk = \f{\pi}{2}\label{Yang-int}\ee It is interesting that the result is independent of $n$. Let us verify (\ref{Yang-int} ) by independent calculations for $n=1,2$: \[ I(1) =\int^\infty_0 \, dk \,e^{-k}\int^\infty_0 \, d\eta \, \cos \, k\eta \, \] Since \[ \int^\infty_0 \, dk \,e^{-k} \cos \, k\eta \, = \f{1}{1+\eta^2}\] and \[ \int_0^\infty d\eta \f{1}{1+\eta^2}= \left.\tan^{-1}\eta \right|_0^\infty = \f{\pi}{2}\] (\ref{Yang-int}) is proven. \[ I(2)=\int^\infty_0 \, dk \,e^{-k^2}\int^\infty_0 \, d\eta \, \cos \, k\eta \, \] Since \[ \int^\infty_0 \, dk \,e^{-k^2} \cos \, k\eta \, = \f{\sqrt{\pi}}{2}\, e^{-\eta^2/4}\] and \[\f{ \sqrt{\pi}}{2}\int_0^\infty d\eta \,e^{-\eta^2/4}= \f{ \sqrt{\pi}}{\sqrt{2}}\int_0^\infty d\eta \,e^{-\eta^2}= \f{\pi}{2}\] (\ref{Yang-int}) is proven again. \end{document}