User:Egm6322.s09.lapetina/project

 See my comments below. Egm6322.s09 14:29, 1 May 2009 (UTC)

=Coastal Hydrodynamic Modeling in Two and Three Dimensions=

=Introduction= My research focuses on two and three dimensional modeling in coastal areas. This is a cutting edge issue in fluid mechanics. In this article we will examine only a few of the many ways partial differential equations are critical to the development of a robust ocean model. Here, particular emphasis will be given to the basic equations, as well as grid generation, an under-appreciated topic.

The report will begin by looking at the continuity and momentum equations, and what simplifications are made for coastal models. Boundary and initial conditions will be discussed, as well as means of solving these equations numerically. A few schemes for solving these will be discussed.

Next, the report will examine the issue of curvilinear grid generation. Structured grid models require the use of orthogonal or curvilinear grids. In the coastal zone, curvilinear grids are required to offer adaptability to surface features. However, their generation is very non-trivial. Formation requires the use of coordinate transformations and highly non-linear partial differential equations.

The last section will discuss my particular areas of interest in the development of a coastal model, touching on the partial differential equations aspect there.

=Modeling With PDEs=

The Typical Equations of Fluid Motion in Oceans
To lead the discussion, the two dimensional equations of motion will be discussed.

Virtually all problems of coastal modeling begin with the continuity equation and momentum equation:

$$\frac{\partial \rho}{\partial t}+\nabla \cdot(\rho \vec v)=0$$

$$\frac{\partial}{\partial t}(\rho \vec v)+ \nabla \cdot (\rho \vec v \vec v)= - \nabla P + \nabla \cdot \mathbf{\tau} + \rho \mathbf{F}$$

where $$\nabla$$ is the differential operator, $$\rho$$ is the density of the fluid, seawater, $$\vec v$$ is the velocity vector, $$P$$ is pressure, $$\tau$$ is the viscous stress tensor, and $$\mathbf{F}$$ are the body accelerations, such as gravity and Coriolis effects.

Proof of Dimensional Homogeneity

In the classroom it was stressed that all good scientists and engineers check the dimensional homogeneity of their equations. Here would be an appropriate place to do so. Using $$L$$ for length, $$T$$ for time, and $$M$$ for mass, the continuity equation is dimensionally expressed as:

$$\frac{\frac{M}{L^3}}{T}=\frac{1}{L} \frac{M}{L^3} \frac{L}{T} $$

which obviously reduces to:

$$\frac{M}{T L^3}$$

Likewise, the dimensions of the momentum equation can be expressed as:

$$\frac {1}{T}\frac{M}{L^3} \frac{L}{T} + \frac{1}{L} \frac{M}{L^3} \frac{L^2}{T^2}=\frac{1}{L} \frac{M \frac{L}{T^2}}{L^2} +\frac{1}{L} \frac{M \frac{L}{T^2}}{L^2}+ \frac{M}{L^3} \frac{L}{T^2}$$

$$\frac{M}{L^2 T^2}+ \frac{M}{L^2 T^2} = \frac{M}{L^2 T^2} +\frac{M}{L^2 T^2}$$.

Clearly both of these equations work out appropriately.

Standard Simplifications
Typically, simplifications and assumptions are made to this tensor equation based upon the physical system of the coastal ocean. Constant density is assumed using the Boussinesq Approximation, removing its time varying term. The Boussinesq Approximation is a bit more complex than simply doing this, and worth a look in other fluid systems like clouds. Due to scaling, we will assume that vertical currents are much smaller than horizontal, eliminating vertical advection. Additionally, the only body forces considered in the horizontal planes are Coriolis effects. Eddy viscosity closure is used to simplify the turbulence terms. With these simplifications, the momentum equations reduce to:

$$\frac {\partial u}{\partial t} + u \frac{\partial u}{\partial x} +u \frac{\partial u}{\partial y} + u \frac {\partial u}{\partial z}+ fu = -\frac{1}{\rho} \frac{\partial P}{\partial x}+ \frac {\partial}{\partial x} \left [ A_x\frac{\partial u}{\partial x} \right ] +\frac {\partial}{\partial y} \left [ A_y \frac{\partial u}{\partial y} \right ]+\frac {\partial}{\partial z} \left [ A_z \frac{\partial u}{\partial z}\right ]$$

$$\frac {\partial v}{\partial t} + v \frac{\partial v}{\partial x} +v \frac{\partial v}{\partial y} + v \frac{\partial v}{\partial z}+ fv = -\frac{1}{\rho} \frac{\partial P}{\partial y}+ \frac {\partial}{\partial x} \left [ A_x\frac{\partial v}{\partial x}\right ] +\frac {\partial}{\partial y} \left [ A_y \frac{\partial v}{\partial y}\right ] +\frac {\partial}{\partial z} \left [ A_z \frac{\partial v}{\partial z}\right ]$$

and

$$\frac {1}{\rho}\frac{\partial P}{\partial z}+g=0$$



The hydrostatic approximation is used in the open ocean. It is rooted in the fact that for the ocean, the depth of the ocean $$H$$ is much less than $$L$$, the length scales of coastal and ocean models. This further reduces the equation by removing horizontal viscous eddy coefficients. Thereby the typical two-dimensional ocean model (including advection) operates with the following equations of motion in the $$x$$ and $$y$$ directions:

$$\frac{Du}{Dt}-fv=-\frac{1}{\rho} \frac{\partial P}{\partial x} + \frac{\partial}{\partial z} \left [ A_v \frac{\partial u}{\partial z} \right ]$$

and

$$\frac{Dv}{Dt}-fu=-\frac{1}{\rho}\frac{\partial P}{\partial y} + \frac{\partial}{\partial z} \left [ A_v \frac{\partial v}{\partial z}\right ]$$

The terms represent inertial terms including advection, Coriolis forces, pressure gradients, and vertical eddy viscosity, respectively.

These equation is a second order, non-linear partial differential equations. The order is two because the last friction terms are second derivatives, and the equation is non-linear because the advective terms are functions of the momentum itself. Assuming a constant vertical eddy viscosity, the $$u$$ equation can be re-expressed as:

$$\frac {\partial u}{\partial t}+u \frac{\partial u}{\partial x} + u \frac{\partial u}{\partial y} -fv =-\frac{1}{\rho} \frac{\partial P}{\partial x} + A_v \frac{\partial^2 u}{\partial z^2}$$

$$A_v \frac{\partial^2 u}{\partial z^2}=-u \frac{\partial u}{\partial x} - u \frac{\partial u}{\partial y} +fv - \frac{1}{\rho} \frac{\partial P}{\partial x} -\frac{\partial u}{\partial t}$$

These equations are of more than 2 independent variables, so classification does not fall into the hyperbolic, elliptic, or parabolic. They are non-linear.

 When you solve these nonlinear PDEs by the Newton-Raphson method, you would need to linearize them. Once linearized, these PDEs look more like a parabolic equation. Egm6322.s09 16:27, 1 May 2009 (UTC)

Gaspard Gustave de Coriolis 1792-1843

Meteorologists and oceanographers are forever indebted to Coriolis, a French engineer and mathematician. Educated at the École Polytechnique, he was a contemporary of Navier and Cauchy. In addition to first describing his namesake effect, he did substantial work on the definition of "work". He died at the age of 51 in Paris.

 Could link to this Wikipedia article on Coriolis. Egm6322.s09 16:27, 1 May 2009 (UTC)

Boundary Conditions
Boundary conditions applied in ocean modeling vary depending upon the model purpose, scale, and goal. A model which includes flow over vertical walls will have different boundary conditions than one measuring saltwater intrusion. Here we will discuss some of the boundary conditions applied to two and three dimensional ocean models.

In describing water wave equations, Dean and Dalrymple provide the three boundary conditions for the three-dimensional momentum equation on the open ocean.

Assuming a distance $$H$$ between mean water level and the bottom of the ocean, at $$z=-H$$, $$\frac{\partial w}{\partial z}=0$$ for all $$x,y$$. Physically, this means that water cannot flow through the bottom of the ocean. It is the bottom boundary condition

The pressure term above $$\frac{1}{\rho} \frac{\partial P}{\partial x} $$ contains several components, including ambient pressure, baroclinic pressure gradients due to salinity gradients, and barotropic pressure gradients due to differential elevations. This latter term is of primary concern here, and introduces boundary conditions. In the equations above, the barotropic pressure term can be expressed as:

$$P_{barotropic_x}= \frac {1}{\rho} \frac{\partial \zeta }{\partial x}$$

$$P_{barotropic_y}= \frac {1}{\rho} \frac{\partial \zeta }{\partial y}$$

where $$\zeta$$ is the free surface of water, and has an average value of zero.

Like in the bottom boundary condition, the vertical velocity, $$w$$ cannot break the free surface. Therefore, at the elevation $$z=\zeta$$, for all $$x$$ and $$y$$,

$$w= \frac{\partial \zeta}{\partial t}+ u \frac{\partial \zeta}{\partial x} +v \frac{\partial \zeta}{ \partial y}$$.

This is known as the kinematic free-surface boundary condition.

A third boundary condition, of special importance in coastal models used for hurricanes, is the dynamic free surface boundary condition. It states that on free surfaces such as the ocean, Bernoulli's Equation applies, and the surface must respond to changes in pressure by physically moving.

In our models, we apply other boundary conditions in the $$x$$ and $$y$$ directions. Tidal boundary conditions are applied using the principle of superposition. Each phase of the tide is added linearly and as time progresses, the amplitude of water at the open boundary changes in a modified sinusoid. This, along with other inputs, such as wind, rivers, and salinity, causes momentum within the coastal model.

Methods to Solve
Typically the momentum and continuity equations are grouped to solve for wave propagation within the coastal area. Each cell has a value of $$\zeta$$ and as many currents as there are dimensions associated with it. Solutions for the entire set of differential equations are of two types, iterative and exact solutions. One of my favorite methods of solving the many differential equations is the Conjugate Gradient Method. It is an exact solution found by an iteration, but it provides very accurate answers after only a few iterations. It utilizes both coordinate transformations and conics. Other iterative methods include the Successive Over Relaxation Method, Jacobi Method, and Gauss Seidel.

Conversion from Cartesian to Curvilinear Coordinates
Accurate modeling of coastal areas requires the use of either curvilinear or unstructured grids. Here we will discuss curvilinear grids. use of them requires their generation, which is not a simple concept. As depicted on the right, boundary points ($$A$$ through $$F$$) on the physical plane are chosen. These are typically done at obvious points within the physical domain. Additionally, evenly spaced $$(\xi,\eta)$$ points are chosen on the computational plane. The number of points on the physical plane between points $$A$$ through $$F$$ is also chosen by the user. $$(\xi, \eta)$$ intersection points are integer values, making them easy to determine. However, converting those $$(\xi, \eta)$$ points to $$x,y$$ points on the physical plane is much more complicated. But, once those values are known, important computational values of $$\Delta x$$ and $$\Delta y$$ values can be determined for modeling.



To use curvilinear coordinates, it is critical in grid generation to minimize the skewness of all cells in the transformed grid. We want to minimize the quantity of the second derivative of the curves in the physical plane.

The corner points $$A$$ through $$F$$ have set position both planes, and the number of points between them on both planes is known. The distance between points on the computational plane is in units of 1. We want to find $$(x,y)$$ points corresponding to integer $$(\xi,\eta)$$ values. To do this, we want to express the equations:

$$\eta_{xx}+\eta_{yy}=0$$ and

$$\xi_{xx}+\xi_{yy}=0$$ in terms of:

$$x_{\eta \eta}, x_{\xi\xi}, x_{\xi\eta}, y_{\eta \eta}, y_{\xi\xi}, y_{\xi\eta}, x_\xi, y_\xi, x_\eta$$ and $$y_\eta$$.

This way, if there are no inflection points between the associated $$\eta$$ and $$\xi$$ as functions of $$x$$ and $$y$$, the physical grid will have minimal curvature.

We know the Jacobian to change from $$(x,y)$$ to $$\xi,\eta$$:

$$\mathbf{J} = \left [ \frac {\partial {x_i}}{\partial {\bar x_j}} \right ]=

\begin{bmatrix} x_\xi   &  \  x_\eta   \\ y_\xi    & \ y_\eta \end{bmatrix}

$$

 If $$\displaystyle (x,y) \equiv (x_1, x_2)$$ represent the old coordinates, and $$\displaystyle (\xi, \eta) \equiv (\xi_1, \xi_2)$$, the new coordinates, then the Jacobian matrix for the coordinate transformation $$\displaystyle (x_1, x_2) \mapsto (\xi_1 , \xi_2)$$ is usually defined as $$\displaystyle \mathbf J = \left[ \partial \xi_i / \partial x_j \right]$$, as done in class. So what you wrote would be $$\displaystyle \mathbf J^{-1} = \left[ \partial x_i / \partial \xi_j \right]$$, the inverse of the Jacobian matrix, i.e., the Jacobian matrix of the inverse map $$\displaystyle (\xi_1, \xi_2) \mapsto (x_1 , x_2)$$. Egm6322.s09 16:27, 1 May 2009 (UTC)

and the determinant of this matrix, expressed as $$J$$ (non-boldface),

$$J=x_\xi y_\eta -x_\eta y_\xi$$.

Likewise, to transform in the opposite direction, $$\mathbf{J_{\alpha}} = \left [ \frac {\partial {x_i}}{\partial {\bar x_j}} \right ]=

\begin{bmatrix} \xi_x   &  \  \xi_y   \\ \eta_x    & \ \eta_y \end{bmatrix}

$$.

To develop equations for the equivalence of

$$\eta_{xx}+\eta_{yy}=0$$ and

$$\xi_{xx}+\xi_{yy}=0$$ in terms of:

$$x, y$$, and their derivatives with respect to $$\xi$$ and $$\eta$$, we can use the relationship that:

$$

\begin{bmatrix} \xi_x   &  \  \xi_y   \\ \eta_x    & \ \eta_y \end{bmatrix}

\begin{bmatrix} x_\xi   &  \  x_\eta   \\ y_\xi    & \ y_\eta \end{bmatrix}

=

\mathbf{I_2}

$$ This leaves us with four equations:

$$x_\xi \xi_x + y_\xi \xi_y =1$$

$$x_\eta \xi_x +y_\eta \xi_y =0$$

$$x_\xi \eta_x +y_\xi \eta_y=0$$

and

$$x_\eta \eta_x +y_\eta \eta_y=1$$.

We want to put the derivatives of $$\xi$$ and $$\eta$$ with respect to $$x$$ and $$y$$ in terms of derivatives of $$x$$ and $$y$$ with respect to $$\xi$$ and $$\eta$$, and the determinant of $$\mathbf{J}$$.

Therefore:

$$x_\xi=J \eta_y \; \; y_\xi=-J \eta_x$$

$$x_\eta=-J \xi_y \; \; y_\eta=J \xi_x$$.

We can re-express our original equations as:

$$(\eta_x)_x +(\eta_y)_y=0$$ and

$$(\xi_x)_x +(\xi_y)_y=0$$.

Using substitution:

$$(\frac{-y_\xi}{J})_x +(\frac{x_\xi}{J})_y=0$$ and

$$(\frac{y_\eta}{J})_x +(\frac{-x_\eta}{J})_y=0$$.

Here will will split the two equations, and solve the top one first.

Using the chain rule, we find:

$$\xi_x (\frac{-y_\xi}{J})_\xi + \eta_x (\frac{-y_\xi}{J})_\eta + \xi_y (\frac{x_\xi}{J})_\xi + \eta_y (\frac{x_\xi}{J})_\eta=0$$

Applying the quotient rule and substitution, we arrive at:

$$\left [ \frac{y_\eta}{J} \left ( \frac{J(-y_{\xi\xi})+ y_\xi J_\xi}{J^2} \right ) + \frac{-y_\xi}{J} \left ( \frac{J(-y_{\xi\eta})+ y_\xi J_\eta}{J^2} \right ) + \frac{-x_\eta}{J} \left ( \frac{J(x_{\xi\xi})- x_\xi J_\xi}{J^2} \right ) + \frac{x_\xi}{J} \left ( \frac{J(x_{\xi\eta})- x_\xi J_\eta}{J^2} \right ) \right ]=0$$

If there is a coordinate transformation, $$J$$ is non-zero, so it can be removed from the calculation.

Derivatives of the determinant of the Jacobian are:

$$J_{\xi}= x_{\xi\xi} y_\eta+ y_{\eta \xi}x_\xi - x_{\eta \xi} y_\xi - y_{\xi \xi} x_\eta$$ and

$$J_{\eta}= x_{\xi\eta} y_\eta+ y_{\eta \eta}x_\xi - x_{\eta \eta} y_\xi - y_{\xi \eta} x_\eta$$

Substituting this leaves us with:

$$ y_\eta ( -y_{\xi \xi} x_\xi y_\eta +y_{\xi \xi} x_\eta y_\xi +x_{\xi\xi} y_\eta y_\xi +y_{\eta \xi}x_\xi y_\xi- x_{\eta \xi} y^2_{\xi} -y_{\xi\xi} x_\eta y_\xi )   -   $$

$$ y_\xi ( -y_{\xi \eta} x_\xi y_\eta +y_{\xi \eta} x_\eta y_\xi +x_{\xi \eta} y_\eta y_\xi +y_{\eta \eta}x_\xi y_\xi- x_{\eta \eta} y^2_{\xi} - y_{\xi\eta} x_\eta y_\xi )  - $$ $$

x_\eta ( x_{\xi \xi} x_\xi y_\eta -x_{\xi \xi} x_\eta y_\xi - x_{\xi \xi} y_\eta x_\xi - y_{\eta \xi} x^2_\xi + x_{\eta \xi} y_{\xi} x_\xi + y_{\xi\xi} x_\eta y_\xi )  +   $$

$$ x_\xi ( x_{\xi \eta} x_\xi y_\eta - x_{\xi \eta} x_\eta y_\xi - x_{\xi \eta} x_\xi y_\eta - y_{\eta \eta}x^2_\xi + x_{\eta \eta} y_{\xi} x_\xi + y_{\xi\eta} x_\eta y_\xi ) =0$$.

It is important to note here that all the terms are composed of derivatives of $$x$$ and $$y$$ with respect to $$\xi$$ and $$\eta$$. This is the form we want for our final expression.

To simplify this expression, we group the terms by the second derivatives:

$$y_{\xi\xi} \left [ -y^2_\eta x_\xi - x_\eta y_\eta y_\xi + x_\eta y_\eta y_\xi - x^2_\eta x_\xi \right ] + $$

$$x_{\xi \xi} \left [ y^2_\eta y_\xi - x_\eta y_\eta x_\xi + x_\eta y_\eta x_\xi - x^2_\eta y_\xi \right ] +$$

$$y_{\xi \eta} \left [ x_\xi y_\xi y_\eta +x_\xi y_\xi y_\eta - x_\eta y^2_\xi + x_\eta y^2_\xi + x^2_\xi x_\eta + x^2_\xi x_\eta \right ] +$$

$$x_{\xi \eta} \left [ -y_\xi x_\xi x_\eta - y_\xi x_\xi x_\eta - y_\eta y^2_\xi - y_\eta y^2_\xi + x^2_\xi y_\eta - x^2_\xi y_\eta \right ] $$ $$ y_{\eta \eta} \left [ -x^3_\xi -x_\xi y^2_\xi \right ] + x_{\eta \eta} \left [ y^3_\xi +y_\xi x^2_\xi \right ]=0$$.

We can factor out a $$-x_\xi$$ from three of the terms, and a $$y_\xi$$ from the other three. After canceling out terms, this leaves us with:

$$-x_\xi \left [ y_{\xi\xi} \left [ y^2_\eta + x^2_\eta \right ] + y_{\xi \eta} \left [ -2 y_\xi y_\eta - 2x_\xi x_\eta \right ]  +  y_{\eta \eta} \left [ x^2_\xi +y^2_\xi \right ]  \right ]  $$

$$ + y_\xi \left [ x_{\xi \xi} \left [ y^2_\eta + x^2_\eta  \right ] + x_{\xi \eta} \left [ -2 x_\xi x_\eta - 2 y_\eta y_\xi \right ]  + x_{\eta \eta} \left [ y^2_\xi +x^2_\xi \right ] \right ]=0 $$

Solving the other equation:

$$\xi_{xx}+\xi_{yy}=0$$

using the same substitutions and groupings yields a similar equation:

$$x_\eta \left [ y_{\xi\xi} \left [ y^2_\eta + x^2_\eta \right ] + y_{\xi \eta} \left [ -2 y_\xi y_\eta  - 2x_\xi x_\eta \right ]  +  y_{\eta \eta} \left [ x^2_\xi +y^2_\xi \right ]  \right ]  $$

$$ -y_\eta \left [ x_{\xi \xi} \left [ y^2_\eta + x^2_\eta  \right ] + x_{\xi \eta} \left [ -2 x_\xi x_\eta - 2 y_\eta y_\xi \right ]  + x_{\eta \eta} \left [ y^2_\xi +x^2_\xi \right ] \right ]=0 $$

From the definition of the non-zero Jacobian,

$$J=x_\xi y_\eta -x_\eta y_\xi$$

neither $$x_\xi $$ and $$y_\xi$$ nor $$y_\eta$$ and $$ x_\eta $$ are zero if a transformation occurs.

Thus, if we declare:

$$\bar a=y^2_\eta + x^2_\eta $$

$$\bar b= y_\xi y_\eta + x_\xi x_\eta $$

and

$$\bar c=y^2_\xi +x^2_\xi $$

we are left with:

$$\bar a x_{\xi \xi} -2 \bar b x_{\xi \eta} +\bar c x_{\eta \eta} =0$$

and

$$\bar a y_{\xi \xi} -2 \bar b y_{\xi \eta} +\bar c y_{\eta \eta} =0$$

as the equations for coordinate transformation.

Contravariants, Covariants, and Fluxes
Once grids are generated, the task of expressing the continuity and momentum equations in curvilinear coordinates is far from complete. To have conservative, perpendicular fluxes between grid cells, the use of covariant or contravariant vectors is required. Good documentation on the derivation of these is available here.

=Conclusions=

Turbulence in Modeling
One of my particular areas of interest is wetland-hurricane interaction, where not all these assumptions are valid. Friction losses must be estimated from the three-dimensional Reynolds stresses. These losses take forms dependent upon the conditions of the flow (initial conditions and boundary conditions). For flow through the lower half of a vegetative canopy, varying drags can be applied to the Reynolds Averaged Navier Stokes Equations. However, in flows near the peak of vegetation or over inundated vegetation, momentum exchange becomes much more complex, and three-dimensional considerations are required. The expression of the friction terms for wetland regions is one of my main areas of research interest.

Regardless of how turbulence is treated, it generally appears in the tensor term: $$ \nabla \cdot \mathbf{\tau}$$

from the Navier-Stokes Equations. It is topical to point out that this expression in tensor form is indicative of turbulence, it retains its properties through coordinate transformations, and acts in all dimensions.

Discussion of Elegance and Meaning of Curvilinear Grid Generation Problem
These two equations can be expressed as:

$$

\left \lfloor \partial_\xi \ \partial_\eta \right \rfloor

\begin{bmatrix} \bar a \ -\bar b \\ -\bar b \ \bar c \end{bmatrix}

\begin{Bmatrix}

\partial_\xi x \\ \partial_\eta x

\end{Bmatrix}

=

\left \lfloor \partial_\xi \ \partial_\eta \right \rfloor

\mathbf{\bar A}

\begin{Bmatrix}

\partial_\xi x \\ \partial_\eta x

\end{Bmatrix}

=0$$

and

$$

\left \lfloor \partial_\xi \ \partial_\eta \right \rfloor

\begin{bmatrix} \bar a \ -\bar b \\ -\bar b \ \bar c \end{bmatrix}

\begin{Bmatrix}

\partial_\xi y \\ \partial_\eta y

\end{Bmatrix}

= \left \lfloor \partial_\xi \ \partial_\eta \right \rfloor

\mathbf{\bar A}

\begin{Bmatrix}

\partial_\xi y \\ \partial_\eta y

\end{Bmatrix}

=0$$.

The difficulty of solving this problem rapidly is the fact that $$\bar a$$, $$\bar b$$, and $$\bar c$$, components of the matrix $$\mathbf{\bar A}$$ are functions of $$x$$ and $$y$$. Expressing them as terms of the Jacobian:

$$\mathbf{J} = \left [ \frac {\partial {x_i}}{\partial {\bar x_j}} \right ]=

\begin{bmatrix} x_\xi   &  \  x_\eta   \\ y_\xi    & \ y_\eta \end{bmatrix}=

\begin{bmatrix} J_{11}   &  \  J_{12}   \\ J_{21}    & \ J_{22} \end{bmatrix}

$$

shows $$\mathbf{\bar A}$$ to contain:

$$ \begin{bmatrix} J_{22}^2 +J_{21}^2   &  \ -(J_{21}J_{22}+J_{11}J_{12})  \\ -(J_{21}J_{22}+J_{11}J_{12})  & \ J_{21}^2 +J_{11}^2 \end{bmatrix} $$.

The complexity of solving the problem of minimizing the curvatures while generating curvilinear grids is demonstrated in this matrix. It is not readily decomposed, reduced, or inverted. Innovation in this field would greatly enhance grid generation speeds and potentially accuracies.

Example of grids generated using this transformation are here and here.

=Concluding Remarks=

This paper has made evident the importance of partial differential equations in all aspects of coastal modeling, not merely in solving equations. Partial differential equations are needed even to build grids on which to use a model, and advances are needed to solve these equations faster. Expression of existing problems in forms such as tensors and matrices may provide opportunities for insights not achieved by looking at the polynomial expressions of problems.

=Signature= --Egm6322.s09.lapetina 04:06, 30 April 2009 (UTC)

=References=