\documentclass[twoside]{article}
\usepackage{amssymb} % used for R in Real numbers
\usepackage{epsfig}
\pagestyle{myheadings}
\setcounter{page}{149}
\markboth{\hfil  Governing Equations of Fluid Mechanics \hfil}%
{\hfil Swungho Lee \&  Bharat K. Soni \hfil}
\begin{document}
\title{\vspace{-1in}\parbox{\linewidth}{\footnotesize\noindent
{\sc  Differential Equations and Computational Simulations III}\newline
J. Graef, R. Shivaji, B. Soni, \& J. Zhu (Editors)\newline
Electronic Journal of Differential Equations, Conference~01, 1997, pp. 149--157. \newline
ISSN: 1072-6691. URL: http://ejde.math.swt.edu or http://ejde.math.unt.edu
\newline ftp  ejde.math.swt.edu  (login: ftp)}
 \vspace{\bigskipamount} \\
  Governing Equations of Fluid Mechanics  in
  Physical Curvilinear Coordinate System
\thanks{ {\em 1991 Mathematics Subject Classifications:} 65N30, 35Q30. 
 \hfil\break\indent
{\em Key words and phrases:} Physical curvilinear coordinate systems,
\hfil\break\indent  incompressible flows, transformations.
\hfil\break\indent
\copyright 1998 Southwest Texas State University  and University of
North Texas. \hfil\break\indent
Published November 12, 1998.} }

\date{}
\author{Swungho Lee \&  Bharat K. Soni}
\maketitle

\begin{abstract}
This paper presents the development of unsteady three-dimensional
incompressible Navier-Stokes and Reynolds-averaged Navier-Stokes equations
in an unsteady physical curvilinear coordinate system.  It is demonstrated
that the numerical simulations based on governing equations in  a physical
curvilinear coordinate system are less mesh sensitive than other schemes.
\end{abstract}

\section*{Introduction}

The ultimate goal of the current research is to develop a numerical flow
simulation system applicable to unsteady three dimensional incompressible
 Navier-Stokes equations that is accurate and efficient in view of CPU time
taken for convergence. Patel~et~al.~\cite{pat} note that,
\begin{quote}If the geometry is highly curved and skewness of angles between
the velocity components and the faces of the computational cells are large,
an approach that transforms only the independent coordinate variables in
the equations representing the transport of mass and momentum may lead
to increased numerical diffusion.\end{quote}
This fact provides the motivation for the development of governing equations
in physical curvilinear component form.  In this coordinate system,
 components of the velocity have the same direction as that of the coordinate
lines and have physical values.  The physical curvilinear-components form
of the velocity was first introduced by Truesdell~\cite{tru}.
Demirdzic~et~al.~\cite{dem} derived the physical curvilinear-components
form in nonorthogonal coordinates for Reynolds-averaged Navier-Stokes
equations and the equations of the $k - \epsilon$ turbulence model.
In their derivation, the equations of the Cartesian tensor forms were
transformed directly into physical curvilinear-component forms by a two step
procedure. Takizawa~et~al.~\cite{tak} used this form to simulate a
two-dimensional free-surface problem using a different concept for the
connection coefficients. However, in this approach, practical applications
have been limited to two-dimensional problems because of the large storage
requirements of geometric tensors and connection coefficients, as well as
the numerical error associated with the evaluation of these geometric
tensors and connection coefficients.

This work explores an approach which is different from
Demirdzic~et~al~\cite{dem} in that the partial differential equations with
the coordinate-free vector form are transformed into the physical curvilinear
coordinate system using general transformation laws.  The resulting unsteady
three-dimensional incompressible viscous equations based on the  unsteady
physical curvilinear coordinate system derived in this research are validated
by numerically simulating  the  laminar three-dimensional lid-driven cavity
and the free surface turbulent flows.It is demonstrated that these numerical
simulations are less mesh sensitive.

\section*{Basic transformations}

In the following analysis, $x_i$ are Cartesian coordinates, $\xi^i$ are
curvilinear coordinates and $\xi^{(i)}$ are physical curvilinear coordinates.
 First, consider the general transformation law under the two changes of the
coordinates system $x_i$ to $\xi^i$ and $\xi^i$ to $\xi^{(i)}$.  The
relationship between the coordinates $x_i$ and the coordinates $\xi^i$ can be
expressed as follows:
\begin{equation}
\xi^i=\xi^i(x_j, t)
\end{equation}
The relationship between the coordinates $\xi^i$ and the coordinates
$\xi^{(i)}$ can be expressed as:
\begin{equation}
\xi^{(i)}
= \xi^{(i)}(\xi^i)\,, \textrm{ where }\Delta\xi^{(i)}
= \sqrt{\mathrm{g}_{ii}}\Delta\xi^i
\end{equation}
$\sqrt{\mathrm{g}_{ii}}$ are evaluated at $\xi^k$ = constant and $k \neq i$.
 $\xi^{(i)}$ resemble the coordinate stretching in each direction of $\xi^i$.
 In view of transforming the coordinates from $x_i$ to $\xi^i$ and $\xi^i$ to
$\xi^{(i)}$, the vector $d\vec r$ can be written as:
\begin{eqnarray}
d\vec r
& = & {\frac{\partial\vec r}{\partial\xi^i}}d\xi^i
= \vec a_id\xi^i\qquad\textrm{(sum on i)}  \\
& = & {\frac{\partial\vec r}{\partial\xi^{(i)}}}d\xi^{(i)}
= \vec a_{(i)}d\xi^{(i)}\,, \textrm{ where } \vec a_{(i)}
= {\frac{1}{\sqrt{\mathrm{g}_{ii}}}}\vec a_i
\end{eqnarray}
$\vec a_i$ are covariant base vectors in the curvilinear coordinate system
and $\vec a_{(j)}$ are covariant base vectors in the physical curvilinear
 coordinate system. The relationships between each coordinate system for
each covariant and contravariant metric tensors are written as:
\begin{eqnarray}
\vec a_{(i)} \cdot \vec a_{(j)}
= \mathrm{g}_{(ij)}
& = & {\frac{1}{{\sqrt{\mathrm{g}_{ii}}}{\sqrt{\mathrm{g}_{jj}}}}}
\mathrm{g}_{ij} \,,
\textrm{ where }
\vec a_{i}\cdot \vec a_{j}
= \mathrm{g}_{ij} \\
\vec a^{(i)}\cdot \vec a^{(j)}
= \mathrm{g}^{(ij)}
& = & \sqrt{\mathrm{g}_{ii}}\sqrt{\mathrm{g}_{jj}}\mathrm{g}^{ij} \,,
\textrm{ where } \vec a^{i}\cdot \vec a^{j}
= \mathrm{g}^{ij}\,, \quad \vec a^{(i)}
= \sqrt{\mathrm{g}_{ii}}\vec a^{i}
\end{eqnarray}
$\mathrm{g}_{(ij)}$ and $\mathrm{g}^{(ij)}$ are the physical covariant and
the physical contravariant metric tensors, respectively.  The physical
covariant metric tensor $\mathrm{g}_{(ij)}$ are equal to one if the
subscripts $i$ and $j$ are the same.  To obtain the divergence, gradient
and Laplacian operators of a vector in the physical curvilinear coordinate
system, one starts from the covariant and the contravariant derivatives of
the base vectors.  Before obtaining the divergence, however, one needs to
define the physical curvilinear components of the velocity vector.
These components can be defined as the magnitude of the $i^\mathrm{th}$
component projected onto the $i^\mathrm{th}$ physical curvilinear coordinate
direction,
\begin{equation}
u^{(i)} = \vec u \cdot \vec a_{(i)}\,.
\end{equation}
Here $u^{(i)}$ represents the physical curvilinear component of the velocity
vector and has physical values.  In the Cartesian coordinate system,
$u^{(i)}$ are identical to the physical components of the velocity $u{(i)}$.

The derivatives of the covariant base vector in the physical curvilinear
coordinate system are obtained by taking the derivative of the covariant
base vector in the curvilinear coordinate system.
\begin{eqnarray} \label{six}
{\frac{\partial\vec a_{(i)}}{\partial\xi^{(j)}}}
& = & [{\frac{\sqrt{\mathrm{g}_{kk}}}{{\sqrt{\mathrm{g}_{ii}}}
{\sqrt{\mathrm{g}_{jj}}}}}\Gamma_{ij}^{k}-\delta_{i}^{k}
{\frac{\mathrm{g}_{km}}{{\mathrm{g}_{ii}}{\sqrt{\mathrm{g}_{jj}}}}}
\Gamma_{ij}^{m}\,]\vec a_{(k)} \quad
\textrm{(sum on k and m)} \nonumber \\
& = & \Gamma_{(ij)}^{(k)}\vec a_{(k)} \,\,
\textrm{(sum on k)} \,,
\end{eqnarray}
where
$$ \Gamma_{(ij)}^{(k)}
= {\frac{\sqrt{\mathrm{g}_{kk}}}{{\sqrt{\mathrm{g}_{ii\!}}}
{\sqrt{\mathrm{g}_{jj\!}}}}}\Gamma_{ij}^{k}-\delta_{i}^{k}
{\frac{\mathrm{g}_{km}}{{\mathrm{g}_{ii\!}}{\sqrt{\mathrm{g}_{jj\!}}}}}
\Gamma_{ij}^{m}\,.
$$
Similarly,  the derivatives of the contravariant base vector in the
physical curvilinear coordinate system are obtained by taking the
derivative of the contravariant base vector.
\begin{equation}
{\frac{\partial \vec a^{(i)}}{\partial \xi^{(j)}}}
= -\Gamma_{(kj)}^{(i)}\vec a^{(k)} \quad
\textrm{(sum on k)}
\end{equation}
Here the Christoffel and physical counterparts of the Christoffel symbols
have the following properties:
\begin{equation}
\Gamma_{jk}^{i}
= \Gamma_{kj}^{i} \,\,
\textrm{and} \,\,
\Gamma_{(jk)}^{(i)}
\neq\Gamma_{(kj)}^{(i)}
\end{equation}

Using the general transformation laws for a scalar $\phi$, the gradient can
be written as follows:
\begin{equation}
\nabla\phi
= {\frac{\partial\phi}{\partial\xi^{(i)}}}\vec a^{(i)}
= \mathrm{g}^{(ik)} {\frac{\partial\phi}{\partial\xi^{(k)}}}\vec a_{(i)}
\quad
\textrm{(sum on i and k)}\,,
\end{equation}
where $\vec a^{(i)} = \mathrm{g}^{(ik)}\vec a_{(k)}$.
Also, the gradient of a vector $\vec u$ can be expressed using equation  as:
\begin{eqnarray*}
\nabla\vec u
& = & {\frac{\partial\vec u}{\partial\xi^{(i)}}}\vec a^{(i)}
= {\frac{\partial(u^{(k)}\vec a_{(k)})}{\partial\xi^{(i)}}}\vec a^{(i)}
\quad \textrm{(sum on i and k)}  \\
& = & u_{(,i)}^{(k)}\vec a_{(k)}\vec a^{(i)}
= \mathrm{g}^{(ij)}u_{(,i)}^{(k)}\vec a_{(k)}\vec a_{(j)} \,,
\end{eqnarray*}
where
$$ u_{(,i)}^{(k)} =
{\frac{\partial u^{(k)}}{\partial\xi^{(i)}}} + u^{(m)}\Gamma_{(im)}^{(k)}\,.
$$
The quantity $u_{(,j)}^{(k)}$ is called the covariant derivative of the
physical curvilinear components of a vector $\vec u$. One can easily evaluate
the divergence of a vector $\vec u$ using equation (\ref{six})  as:
\begin{eqnarray}
\nabla\cdot\vec u
& = & {\frac{\partial\vec u}{\partial\xi^{(i)}}}\cdot\vec a^{(i)}
= {\frac{\partial u^{(k)}\vec a_{(k)}}{\partial\xi^{(i)}}}\cdot\vec a^{(i)}
\quad \textrm{(sum on i and k)} \nonumber \\
& = & u_{(,i)}^{(i)}
= {\frac{\sqrt{\mathrm{g}_{ii}}}{{J}}} {\frac{\partial}{\partial\xi^{(i)}}}
({\frac{{J}}{\sqrt{\mathrm{g}_{ii}}}}u^{(i)}) \quad \textrm{(sum on i)}
\end{eqnarray}
The Laplacian can be evaluated by the divergence of the gradient of a vector
$\vec u$,
\begin{eqnarray*}
\nabla^{2} \vec u
&\! = &\! \nabla \cdot \hat{T}
= {\frac{\partial \hat{T}}{\partial\xi^{(j)}}} \cdot \vec{a}^{(j)}
\,, \textrm{ where } \quad
\hat{T} = \nabla_{\vec u}
= u_{(i,)}^{(k)} \vec a_{(k)} \vec a^{(i)}  \\
&\! = &\! \textrm{g}^{(jk)}{} [u_{(,j)}^{(m)} \Gamma_{(m\!k)}^{(i)}
-{} u_{(,m)}^{(i)} \Gamma_{(j\!k)}^{(m)}
+ {\frac{\partial u_{(,j)}^{(i)}}{\partial \xi^{(k)}}}] \vec a_{(i)}
\, \textrm{(sum on i, j, k and m)}
\end{eqnarray*}

The time derivative of a scalar or a vector F is given as
\mbox{Warsi}~\cite{war}:
\begin{eqnarray} \label{13a}
{\frac{\partial F}{\partial t}} |_{x_i}
& = & [{\frac{\partial F}{\partial \tau}} + {\frac{\partial F}{\partial
\xi^{(i)}}}
{\frac{\partial \xi^{(i)}}{\partial t}}]|_{\xi^(i)}  \\ \label{13b}
{\frac{\partial F}{\partial \tau}} |_{\xi^(i)}
& = & [{\frac{\partial F}{\partial t}} + {\frac{\partial F}{\partial x_{i}}}
{\frac{\partial x_{i}}{\partial \tau}}] |_{x_i}
\end{eqnarray}
Here $\tau$ represents the time in an unsteady physical curvilinear coordinate
system. The grid speed can be easily evaluated by replacing $F$ with
$\xi^{(i)}$ in equation (\ref{13b}).
\begin{equation}
w^{(i)}  =  {\frac{\partial \xi^{(i)}}{\partial t}}
= - {\frac{\partial \xi^{(i)}}{\partial x_{l}}} {\frac{\partial x_{l}}{\partial
\tau}}
= - {\frac{{\sqrt{\mathrm {g}}_{ii}}}{J}} b_{l}^{i} {\frac{\partial
x_{l}}{\partial \tau}} \quad
\textrm{(sum on l)}
\end{equation}
where
$$ b_{l}^{i}
  =  {\frac{\partial x_{m}}{\partial \xi^{j}}} {\frac{\partial x_{n}}
{\partial \xi^{k}}}- {\frac{\partial x_{n}}{\partial \xi_{j}}}
{\frac{\partial x_{m}}{\partial \xi^{k}}}
$$
with $i, j, k$, and $l, m, n$ in each cyclic order.

The divergence of the grid speed vector $\vec w$ is written as:
$$
\nabla \cdot \vec w
 =  {\frac{{\sqrt{\mathrm {g}}_{ii}}}{J}} {\frac{\partial({\frac{J}
{\sqrt{\mathrm {g}}_{ii}}} w^{(i)})}{\partial \xi^{(i)}}}
 =  - {\frac{l}{J}} {\frac{\partial J}{\partial \tau}}
$$

\section*{Governing equations in the unsteady \\ physical curvilinear
coordinate system}
By replacing $F$ with $\vec u$ in equation (\ref{13a}), one can get the
equations into an unsteady coordinate system.  The vector form of
incompressible Reynolds-averaged Navier-Stokes equations, with the body
force in unsteady coordinates system, is given by equation (\ref{16}).
\begin{equation} \label{16}
{\frac{\partial \vec u}{\partial \tau}} + (\nabla \vec u) \cdot \vec v
= - \nabla P + \nu_{E} \nabla^{2} \vec u + [\nabla \vec u +
(\nabla \vec u)^{T}] \cdot \nabla \nu_{E}\,,
\end{equation}
where
$$\vec v = \vec u + \vec w,{} P
= \rho + {\frac{z}{Fn^{2}}} + {\frac{2}{3}}k \, \textrm{and }\, \nu_{E}
= {\frac{l}{R_{\mathrm{eff}}}}
= {\frac{l}{R_{\mathrm{e}}}} + \nu_{t}\,.
$$
Here $\vec w$ is a grid speed vector, $Fn$ is a Froude number, $P$ is total
pressure, and
$\rho$ is a static pressure. $\nu_{t}$ and $R_{\mathrm{eff}}$ represent the
eddy viscosity and the effective Reynolds number, respectively.

The procedure for the transformation of the incompressible Reynolds-averaged
Navier-Stokes equations, based on an unsteady physical curvilinear
coordinate
system, is now introduced using the derivations for the gradient, divergence
operator, Laplacian, and time derivative.  An unsteady physical curvilinear
component form of the continuity and the Reynolds-averaged Navier-Stokes
equations can be presented as:
\begin{eqnarray}
&&{\frac{\partial}{\partial \xi^{i}}} ({\frac{J}{\sqrt\mathrm{g}_{ii}}}
u^{(i)}) \ = \ 0 \nonumber \\
\lefteqn{ R_{\mathrm{eff}} {\frac{\partial u^{(i)}}{\partial\tau}} +
R_{\mathrm{eff}} v^{(j)} ({\frac{\partial u^{{(i)}}}{\partial \xi^{(j)}}} +
u^{(k)} \Gamma_{(kj)}^{(i)} ) }  \label{17b} \\
&=&
- R_{\mathrm{eff}}\, [\mathrm{g}^{(ij)} {\frac{\partial P}{\partial \xi^{(j)}}}
- {\frac{\partial \nu_{E}}{\partial \xi^{(j)}}}
[\mathrm{g}^{(jk)} {\frac{\partial u^{(i)}}{\partial\xi^{k}}}
+ \mathrm{g}^{(jk)} \Gamma_{(lk)}^{(i)} u^{(l)}
+ \mathrm{g}^{(ik)} {\frac{\partial u^{(i)}}{\partial\xi^{k}}} \nonumber \\
&&+ \mathrm{g}^{(ik)}  \Gamma_{(lk)}^{(j)} u^{(l)}]\,]
+ \mathrm{g}^{(jk)} [{\frac{\partial^{2} u^{(i)}}{\partial \xi^{(k)} \partial
\xi^{(j)}}}
+ {\frac{\partial (u^{(l)} \Gamma_{(lj)}^{(i)})}{\partial \xi^{(k)}}}
+ \Gamma_{(lk)}^{(i)} u_{,(j)}^{(l)}
- \Gamma_{(jk)}^{(l)} u_{,(l)}^{(i)}] \nonumber
\end{eqnarray}

Equation (\ref{17b}), which has a nonconservative form, is rearranged into
the standard form for the use of the finite analytic method \cite{swun}.
Using the stretched coordinates $\xi^{i}\ast$, the 12-point finite analytic
discretization scheme based on the local nonuniform grid spacing \cite{swun}
is employed.   The stretched coordinates $\xi^{i}\ast$ are defined as:
\begin{displaymath}
\xi^{i} = {\sqrt{\mathrm{g}^{ii}}} \xi^{i}\ast  \textrm{ or } \xi^{i}\ast =
{\frac{l}{\sqrt{\mathrm{g}^{ii}}}} \xi^{i}
\end{displaymath}
The first equation shows the relation between the curvilinear coordinates
and the stretched coordinates. The $\sqrt{\mathrm{g}^{ii}}$ are evaluated
at $\xi^{k}$ = constant and $k\neq i$.

\section*{Results and discussions}

An unsteady three-dimensional incompressible flow solver based on the physical
curvilinear coordinate system has been developed \cite{swun}. The 12-point
finite analytic scheme with enhanced kinematic boundary condition and
numerical approach was utilized in this development. The detailed
discussions on the numerical scheme can be found in \cite{swun}.
The results of the  following two test cases to validate the pertinent
numerical simulations are presented here.

\subsection*{Lid-driven three-dimensional cavity flow}

The velocity components on the wall are zero, except on the moving wall with
a velocity of  1.   The computations are performed on a grid consisting of
\mbox{$16\!\times\!16\!\times\!16$} grid points.  The dimensionless time
step is taken from 0.01 to 2. The matrix that consists of the coefficients
resulting from the finite analytic method was solved using the GMRES
(Generalized Minimal RESidual) method \cite{bal}.  To
obtain the solution for the steady state, only one iteration per time step
is used.  All computations are performed using a relaxation factor of 1.
Figure 1 shows that the rate of the convergence depends on the size of the
time step in the range from 0.01 to 2.  The computations lead to a fully
converged solution within fewer than 200 iterations.
\begin{figure}%[hbtp]
\epsfig{file=fig1.eps, width=4.5in}
\caption{Convergence Histories for Different Time Steps}
\end{figure}

\begin{figure}%[hbtp]
\epsfig{file=fig2.eps, width=4.5in}
\caption{Comparison of the Centerline u-Velocity Profile}
\end{figure}

Figure 2 shows the comparison of the center line u-velocity profiles at the
different Reynolds numbers.  The computations were performed on both uniform
and nonuniform grids. The time increment is set to $\Delta t = 1$.

\mbox{Peric~\cite{per}} has reported that if the angle between the two
coordinate lines is greater than $135^\circ$ or less than $45^\circ$, then
the pressure correction equation does not converge at all, or the convergence
rate is too slow.  \mbox{Cho and Chung~\cite{cho}} used a new treatment method
for nonorthogonal terms in the pressure-correction equation in order to enlarge
the ranges for convergence and found that the smaller the angle, the narrower
the region of relaxation factor.  In the present research, the computations are
performed at several inclined angles \\  ($90^\circ, 60^\circ, 45^\circ,
30^\circ, 15^\circ,10^\circ, \textrm{and} 5^\circ$) to check the influence on
the  rate of the convergence due to the grid skewness.
As mentioned, all computations of the  three-dimensional cavity flow are
performed
using a relaxation factor of 1.  The solution always converged, even for very
small
inclined angles,  but more iterations were required for convergence for small
inclined angles.  Figure 3 shows the convergence histories in various inclined
angles from 90 to 5 degrees.

\begin{figure}%[hbt]
\epsfig{file=fig3.eps,width=4.5in}
\caption{Convergence Histories in Various Inclined Angles (Degree) }
\end{figure}

\subsection*{Ship flow with the free surface}

It has been shown that the present code is less mesh sensitive and converges
well
even at the large grid skewness for the three-dimensional cavity flow.  For
the next case, three-dimensional moving free surface turbulent flow was
simulated.
 The upper boundaries move arbitrarily with the flow, and the grid in the
computational
domain is generated at every time step until the solution of the steady state
is obtained.
\begin{table} [!hbp] \caption{Grid Dependence Tests for the Wigley Hull}
\label{table 1}
\begin{center}
\begin{tabular}{|lccc|}  \hline
& I & II & III \\ \hline
Grid Points & $125\times35\times34$ & $125\times40\times40$ &
$125\times50\times48$ \\
Total Nodes & 148,750 & 200,000 & 300,000 \\
Time Increment & 0.005 & 0.005 & 0.005 \\
Total time steps & 400 & 400 & 400 \\
Total CPU (hours) & 34.91 & 46.56 & 60.68 \\
Reynolds Number & $1.0\times10^6$ & $1.0\times10^6$ & $1.0\times10^6$ \\
Froude Number & 0.289 & 0.289 & 0.289 \\ \hline
\end{tabular}
\end{center}
\end{table}

Table~\ref{table 1} shows the information for these computations.
The computations were performed with three different grids under the same
conditions. The Baldwin-Lomax
turbulence model \cite{bal} was used to calculate the eddy viscosity in the
turbulent
flow.   The Froude number and the Reynolds number used in the experiment
\cite{coop}
are 0.289 and $3.3 \times 10^6$, respectively.  The wave elevations on the hull
surface and a perspective view of the wave elevation is shown in Figure 4.  A
deviation
of wave profiles is observed in the bow region, while better agreement is seen
toward the stern. The bow peak is not captured properly due to the large
spacing of ($\Delta \textrm{x}$) in a region of relatively high gradient.   The
residue remains around $10^{-4}$after $t = 0.5$.

\begin{figure}%[hbtp]
\epsfig{file=fig4.eps,width=4.5in}
\caption{The Wave Elevations on the Hull Surface for different grid sizes and a
Perspective view of the wave elevation}
\end{figure}


The comparison of these numerical simulations with the results reported in
open literature have shown very good agreement. It is demonstrated, especially
in the cavity flow simulation, that the numerical simulations involving a
physical
curvilinear coordinate system are less mesh sensitive.

\begin{thebibliography}{99}

\bibitem{pat} Patel, V.~C., Chen, H.~C. and Ju,~S.,
``Solutions of the Fully-Elliptic Reynolds-Averaged Navier-Stokes Equations
and Comparisons with Experiments,'' IIHR Report Number~323, The University
of Iowa, 1988.

\bibitem{tru} Truesdell, C., ``The Physical components of vectors and
Tensors,''Journal of  Applied Math. Mech., Vol~33, pp~345--356, 1953.

\bibitem{dem} Demirdzic, I., Gosman, A.~D., Issa, R.~I., and Peric, M.,
``A Calculation Procedure for Turbulent Flow in Complex Geometries,''
Computer \& Fluids Vol.~15, No.~3, pp~251--273, 1987.

\bibitem{tak} Takizawa, A., Koshizuka, S. and Kondo, S.,
``Generalization of Physical Component
Boundary Fitted Co-ordinate (PCBFC) Method for the Analysis of Free-Surface
Flow,'' International Journal for Numerical Methods in Fluids, Vol.~15,
pp~1213-1237, 1992.

\bibitem{war} Warsi, Z.~U.~A., ``Conservation Form of the Navier-Stokes
Equations in General Nonsteady Coordinates,'' AIAA Journal,
Vol.~19, No.~2, 1981.

\bibitem{swun} Swungho Lee, ``Three-Dimensional Incompressible Viscous
Solutions Based on the Physical Curvilinear Coordinate System,''
Dissertation, Department of Computational
Engineering, Mississippi State University, December 1997.

\bibitem{per} Peric, M., ``Analysis of Pressure Velocity Coupling on
Non-Staggered Grids,'' Numerical Heat Transfer, Part~B, Vol.~17, pp~63--82,
1990.

\bibitem{cho} Cho, M.~J. and Chung M.~K., ``New Treatment of Nonorthogonal
Terms in the Pressure-Correction Equation,'' Numerical Heat Transfer,
Part~B, Vol.~26, pp~133--145, 1994.

\bibitem{bal} Baldwin, B.~S. and Lomax, H., ``Thin Layer Approximation and
Algebraic Model for Separated Turbulent Flows,'' AIAA 78--257, 1978.

\bibitem{coop} ``Cooperative Experiments on Wigley Parabolic Models in
Japan,'' 17th ITTC Resistance Committee Report, 1983.

\bibitem{chen} Chen, C.~J., Bravo, R.~H., Chen, H.~C. and Xu, Z.,
``Accurate Discretization of Incompressible Three-Dimensional
Navier-Stokes Equations,'' Numerical Heat Transfer, Part~B, Vol.~27,
pp~371--392, 1995.

\bibitem{ros} Rosenfeld, M. and Kwak, D., ``A Fractional Step Solution
Method for the Unsteady Incompressible Navier-Stokes Equations in
Generalized Coordinate Systems,'' Journal of Computational Physics,
Vol.~94, pp~102--137, 1991.

\bibitem{kar} Karki, K.~C., Sathyamurthy, P.~S. and Patankar, S.~V.,
``Performance of a Multigrid Method with an Improved Discretization Scheme
 for Three-Dimensional Fluid Flow Calculations,''
Numerical Hear Transfer, Part~B, Vol.~29, pp~275--288, 1996.

\end{thebibliography}
\bigskip

\noindent{\sc Swungho Lee} \\
Hydrodynamics Research Dept. \\
Hyundai Maritime Research Institute \\
1, Cheonha-Dong, Dong-Ku, Ulsan, Korea \\
E-mail address: slee@hhi.co.kr
  \medskip

\noindent{\sc Bharat K. Soni} \\
Department of Aerospace Engineering \\
Mississippi State University \\
Mississippi State, MS  39210   USA \\
E-mail address:  bsoni@erc.msstate.edu


\end{document}