\documentstyle[11pt]{article}
\textwidth=156mm
\textheight=240mm
\topmargin=-1.7cm
\oddsidemargin=0.0mm
\evensidemargin=0.0mm
\pagestyle{empty}

\begin{document}

\noindent
\begin{center}
{\large\bf NUMERICAL SIMULATIONS OF SELF-ORGANIZATION IN
2-D CIRCULAR SHEAR FLOWS\\}
\end{center}

\vspace{5mm}
\noindent
K. Bergeron$^1$, E.A. Coutsias$^1$,
J.P. Lynov, \underline{A.H. Nielsen.}\\

\vspace{1ex}
\noindent
Association EURATOM - Ris\o$\;$National Laboratory,
Optics and Fluid Dynamics Department, \\
DK-4000 Roskilde, Denmark 

\vspace{1ex}
\noindent
$^1$Department of Mathematics and Statistics,
University of New Mexico, Albuquerque,\\ NM 87131, USA

\vspace{5mm}
\noindent
A fluid submitted to a velocity shear can be unstable, and may form a line
of vortices which is in itself unstable and evolves through vortex-pairing
processes. Such a shear flow is studied numerically in an annulus geometry.
We have performed direct numerical simulations of the experiment of 
Rabaud-Couder$^1$ (and more resently Chomaz et al.$^2$). Here a flow is
enclosed in a vertical cylinder of radius $R_{2}$ and very small height e. 
A disk of radius $R_{1}$ at the centre of each end, rotates at angular velocity 
$\Omega_{1}$, while the rest of the cell rotates at $\Omega_{2}$.
Viscous friction on the containing walls forces the fluid into solid rotations
at $\Omega_{1}$ in the central part and $\Omega_{2}$ in the external annulus
respectively, seperated by a shear zone. Chomaz et al. have demonstrated that 
for low differential rotation (i.e. low $\Delta \Omega = \mid \Omega_{2} -
\Omega_{1} \mid$) the flow will be stable while for increasing $\Delta \Omega$ 
a Pitchfork bifurcation point is reached where the flow becomes unstable and a necklace
of vortices is observed. By increasing $\Delta \Omega$, the size of the vortices will 
grow and a series of transitions towards patterns with even smaller numbers of 
vortices will take place. Arguing that the Ekman layers at the top and bottom bounding 
walls dominate the flow, they proposed a model in which the flow at the mid-plane is 
essentially two dimensional, with the motion of the vessel manifesting itself as 
a body force.

We have modelled the experiments by the Navier-Stokes equation:
\begin{equation}
\frac{\partial \omega}{\partial t} + \vec{u} \cdot \nabla \omega  = 
\nu \nabla^2 \omega + \nabla \times \vec{F} \; ,\label{eq:eq1}
\end{equation}
combined with the Poisson equation
\begin{equation}
\nabla^{2} \psi  =  - \omega  \; , \label{eq:eq2}
\end{equation}
where $\vec{u}(x,y) = \nabla \psi \times \vec{z}$ is the velocity, 
$\omega (x,y)$ the scalar vorticity, $\psi (x,y)$ the stream function, $\nu$ the viscosity, and $\vec{F}$
represent the body force due to the friction on horizontal walls. 

The Navier-Stokes equation, given by eqn.~\ref{eq:eq1} - ~\ref{eq:eq2}, are solved numerically in an annulus geometry employing a fully dealiased, spectral scheme based on Fourier-Chebyshev expansions. This implies very high resolution of boundary layers. The flow field is subject to no-slip conditions, so that it matches the wall velocity at the boundaries. Precomputed solvability constrains are enforced on each Fourier mode of $\omega$. A more detailed description of how these solvability constrains are implemented in the numerical scheme is given by Coutsias and Lynov$^3$. 


The forcing term has to be modelled in order to obtain a closed set of equations. Assuming a parabolic vertical dependence characteristic of a Poiseuille flow between two planes we obtain:
\begin{equation}
\vec{F} \approx \frac{8 \nu}{e^2} ( \vec{u}^* - \vec{u} ) , \label{eq:eq3}
\end{equation}
where $\vec{u}^*$ is the boundary velocity on top or bottom walls. In the 
experiments this velocity is a discontinuous function of the radial distance to rotation axis. Such a discontinuous function is, however, not possible to use in numerical simulations.  Rabaud-Couder$^2$ has shown that the typical length of the shear is of the order of the cell height e. We have thus choosen a 
$\omega^* \vec{z} = \nabla \times \vec{u}^*$ which is constant in two radial regions with a sharp transition region of the scale length e. In order to obtain a steady state solution correcsponding to this forcing we insert the prescribed $\omega^*$ in equation ~\ref{eq:eq1} which transform it to the Helmholtzs equation:
\begin{equation}
(1 - \frac{e^2}{8} \nabla^{2} )\bar{\omega} = \omega^* \; .
\end{equation}
This $\bar{\omega}$ is used as an initial condition for the code.



\vspace{6mm}
\noindent
1. M. Rabaud and Y. Couder, J. Fluid Mech., {\bf 136} (1983) 291-319.

\noindent
2. J.M. Chomaz, R. Rabaud, C. Basdevant and Y. Couder, J. Fluid Mech.,
{\bf 187} (1988) 115-140.

\noindent
3. E.A. Coutsias and J.P. Lynov, Physica D., {\bf 51} (1991) 482-497.


\end{document}

