Numerical grids used in a coastal ocean model with breaking wave effects

In coastal ocean modeling, traditional single-block rectangular (Cartesian) grids have been most commonly used for their simplicity. In many cases, these grids may be not well suited (even at very high resolutions) for regions with complicated physical fields, open boundaries, coastlines, and bottom bathymetry. The numerical curvilinear nearly orthogonal/orthogonal, single/multi-block coastline-following grids for the Mediterranean Sea, Monterey Bay and the South China Sea (SCS) are presented. These grids can be used in coastal ocean modeling to enhance model numerical solutions and save computer resources by giving better treatment of regions with high gradients such as areas of complicated coastlines and steep slopes of shelf breaks, complicated bottom topography, open boundaries, and multi-scale physical phenomenon. Grid generation techniques are used to designed these grids. This kind of grids can also easily increase horizontal resolutions in the subregion of the model domain, without increasing the computational expense, with a higher resolution over the entire domain.A three dimensional coastal ocean model with breaking wave effects is also presented and applied. The ocean system is a primitive equation modeling system with grid generation routines and a turbulent closure which is capable of taking surface breaking wave effects into account. The system also includes a grid package which allows model numerical grids to be coupled with the ocean model. The model code is written for multi-block grids, but a single-block grid is used for the South China Sea (SCS). The model with breaking wave effects and a grid of 121 × 121 grid points are used to simulate the winter circulation of the SCS as an example. The model output of the 60-day run shows the observed upwelling locations in the sea surface salinity field.     

Domain imbedding methods for the Stokes equations

Summary We study direct and iterative domain imbedding methods for the Stokes equations on certain non-rectangular domains in two space dimensions. We analyze a continuous analog of numerical domain imbedding for bounded, smooth domains, and give an example of a simple numerical algorithm suggested by the continuous analysis. This algorithms is applicable for simply connected domains which can be covered by rectangular grids, with uniformly spaced grid lines in at least one coordinate direction. We also discuss a related FFT-based fast solver for Stokes problems with physical boundary conditions on rectangles, and present some numerical results.     

The effects of interplay between the rotation and shoaling for a solitary wave on variable topography

This paper presents specific features of solitary wave dynamics within the framework of the Ostrovsky equation with variable coefficients in relation to surface and internal waves in a rotating ocean with a variable bottom topography. For solitary waves moving toward the beach, the terminal decay caused by the rotation effect can be suppressed by the shoaling effect. Two basic examples of a bottom profile are analyzed in detail and supported by direct numerical modeling. One of them is a constant‐slope bottom and the other is a specific bottom profile providing a constant amplitude solitary wave. Estimates with real oceanic parameters show that the predicted effects of stable soliton dynamics in a coastal zone can occur, in particular, for internal waves.     

Stiffness characteristic and the applications of the method of lines on nonuniform grids

The method of lines is investigated for the numerical solution of the stream-function-and-vorticity form of the Navier-Stokes equations on nonuniform grids. Stiffness characteristics of a linear one-dimensional model equation are examined to establish the feasibility of applying the method to the vorticity equation in two dimensions. The governing equations are transformed from the physical domain with a highly variable grid to a computational domain with a uniform grid. The method of lines is used to solve only the vorticity equation, and the successive-over relaxation technique is used to solve the stream-function equation. It is observed that the transformed governing equations become stiffer with increased concentration of grid points and also as the number of grid points increases. It is also observed that the differencing technique affects the stiffness characteristics. The use of forward differencing is not feasible, and backward differencing is preferable to central differencing for high Reynolds numbers. The results of specific applications for the solution of flow in curved-wall diffusers and a driven cavity demonstrate that the method of lines under certain circumstances is feasible for the numerical solution of physical problems on domains covered with variable grids.     

矩形到任意多边形区域的Schwarz-Christoffel变换数值解法

运用Schwarz-Christoffel变换方法，建立多边形区域到带状区域共形映射数学模型.对于模型中的约束条件和奇异积分问题，根据Riemann(黎曼)原理，建立复参数与实参数互逆变换，消除非线性系统的约束条件；经过合理积分路径的确定，模型中的奇异积分转化为Gauss-Jacobi(高斯 雅可比)型积分；采用Levenberg-Marquardt算法对非线性系统模型进行求解.根据第一类椭圆函数性质，建立了矩形区域到带状区域共形映射数学模型，通过复参数椭圆函数的计算，得到矩形边界与带状区域边界的关系.最后，对8点对称多边形区域与27点不规则条带状区域计算，将不规则封闭区域边界映射到矩形区域边界，矩形区域内的正交网格，通过变换之后在多边形区域内依然满足正交性，为研究不规则区域到规则区域映射的数值计算奠定基础.     

Supersonic flow past a flat lattice of cylindrical rods

Two-dimensional supersonic laminar ideal gas flows past a regular flat lattice of identical circular cylinders lying in a plane perpendicular to the free-stream velocity are numerically simulated. The flows are computed by applying a multiblock numerical technique with local boundary-fitted curvilinear grids that have finite regions overlapping the global rectangular grid covering the entire computational domain. Viscous boundary layers are resolved on the local grids by applying the Navier–Stokes equations, while the aerodynamic interference of shock wave structures occurring between the lattice elements is described by the Euler equations. In the overlapping grid regions, the functions are interpolated to the grid interfaces. The regimes of supersonic lattice flow are classified. The parameter ranges in which the steady flow around the lattice is not unique are detected, and the mechanisms of hysteresis phenomena are examined.     

Nonlinear Wave Equations for Oceanic Internal Solitary Waves

In the coastal ocean, the interaction of barotropic tidal currents with topographic features such as the continental shelf, sills in narrow straits, and bottom ridges are often observed to generate large amplitude, horizontally propagating internal solitary waves. These are long nonlinear waves and hence can be modeled by equations of the Korteweg–de Vries type. Typically they occur in regions of variable bottom topography, with the consequence that the appropriate nonlinear evolution equation has variable coefficients. Further, as these waves can be long‐lived it is necessary to take account of the effects of the Earth's background rotation. We review this family of model evolution equations and some of their pertinent solutions, obtained both asymptotically and numerically.     

The influence of shelf zone topography and coastline geometry on coastal trapped waves

The results of numerical experiments with a model of coastal trapped waves are presented to identify two important features for regional modeling of the interaction of a shelf zone with open ocean. First, a wave train of this type can be formed by wind action at a considerable distance from the place of impact. The waves propagate along a coastline without significant loss of energy, provided that the coastline and shelf zone topography have no features comparable to the Rossby radius. However, the waves lose energy while passing over capes and submarine canyons and when the shelf width decreases. For regional modeling, remote generation of waves must be thoroughly investigated and taken into account. The other feature is that the propagating waves can use part of energy to form density anomalies on the shelf by raising intermediate waters from the adjacent offshore areas of the open ocean. Thus, coastal trapped waves can carry wind energy from wind action areas to other coastal areas to form density anomalies and other types of motion.     

Study of storm surge due to Typhoon Linda (1997) in the Gulf of Thailand using a three dimensional ocean model

Numerical integrations using the three dimensional ocean model based on the princeton ocean model (POM) were applied for the study of both sea level elevation and ocean circulation patterns forced by the wind fields during typhoons that moved over the Gulf of Thailand (GoT). The simulation concerned a case of Typhoon Linda which occurred during November 1-4, 1997. Typhoon Linda was one of the worst storms that passed the Gulf of Thailand and hit the southern coastal provinces of Thailand on November 3, 1997. It caused flooding and a strong wind covering large areas of agriculture and fisheries, which destroyed households, utilities and even human lives. The model is the time-dependent, primitive equation, Cartesian coordinates in a horizontal and sigma coordinate in the vertical. The model grid has 37 × 97 orthogonal curvilinear grid points in the horizontal, with variable spacing from 2 km near the head of the GoT to 55 km at the eastern boundary, with 10 sigma levels in the vertical conforming to a realistic bottom topography. Open boundary conditions are determined by using radiation conditions, and the sea surface elevation is prescribed from the archiving, validation and interpretation of satellite oceanographic data (AVISO). The initial condition is determined from the spin up phase of the first model run, which was executed by using wind stress calculated from climatological monthly mean wind, restoring-type surface heat and salt and climatological monthly mean freshwater flux. The model was run in spin up phase until an ocean model reached an equilibrium state under the applied force. A spatially variable wind field taken from the European Centre for Medium-Range Weather Forecasts (ECMWF) is used to compute the wind stress directly from the velocity fluctuations. Comparison of tendency between the sea surface elevations from model and the observed significant wave heights of moored buoys in the Gulf of Thailand under Seawatch project is investigated. The model predicts the sea level elevation up to 68.5 cm at the Cha-Am area located in the north of where the typhoon strands to the shore. Results of sea level elevation show that there is an area of peak set-up in the upper gulf, particularly in the western coast, and the effects of the storm surge are small at the lower gulf. During the entire period of this study, the surge in the gulf was induced by the northeasterly wind blowing over it.     

Variational method for adaptive mesh generation

A variational method is suggested for generating adaptive grids composed of hexahedral cells. The method is based on the minimization of a functional written on a manifold in a space whose variables are usual spatial coordinates in a physical domain and the components of a monitor vector function. A grid is constructed in the manifold, and its projection onto the physical domain yields an adaptive grid. Examples of adaptive grid generation are given.     

Numerical modelling of free-surface shallow flows over irregular topography with complex geometry

A well-balanced Godunov-type finite volume algorithm is developed for modelling free-surface shallow flows over irregular topography with complex geometry. The algorithm is based on a new formulation of the classical shallow water equations in hyperbolic conservation form. Unstructured triangular grids are used to achieve the adaptability of the grid to the geometry of the problem and to facilitate localised refinement. The numerical fluxes are calculated using HLLC approximate Riemann solver, and the MUSCL-Hancock predictor–corrector scheme is adopted to achieve the second-order accuracy both in space and in time where the solutions are continuous, and to achieve high-resolution results where the solutions are discontinuous. The novelties of the algorithm include preserving well-balanced property without any additional correction terms and the wet/dry front treatments. The good performance of the algorithm is demonstrated by comparing numerical and theoretical results of several benchmark problems, including the preservation of still water over a two-dimensional hump, the idealised dam-break flow over a frictionless flat rectangular channel, the circular dam-break, and the shock wave from oblique wall. Besides, two laboratory dam-break cases are used for model validation. Furthermore, a practical application related to dam-break flood wave propagation over highly irregular topography with complex geometry is presented. The results show that the algorithm can correctly account for free-surface shallow flows with respect to its effectiveness and robustness thus has bright application prospects.     

A technology for 3D grid generation

New results concerning the technology for adaptive 3D grid generation are obtained. This technology is based on the numerical solution of the inverted one-, two-, and three-dimensional Beltrami and diffusion equations with respect to the monitor metric. The one- and two-dimensional equations are used to generate grids on the edges and the faces of the domain, respectively. The three-dimensional equations are used to generate a grid in the interior of the domain. Examples of adaptive 3D hexahedral and prismatic grids are demonstrated.     

Modeling fully developed laminar flow in a helical duct with rectangular cross section and finite pitch

The mass and momentum transport equations are written in an orthogonal coordinate system using Germano’s transformation to model a laminar flow in a helical duct with a rectangular cross section and finite pitch. The system of governing equations are discretized and solved by the finite-volume numerical method. The three dimensional domain is reduced to a two dimensional slab of cells, orthogonal to the main flow direction, enforcing the fully developed state for 2π/(τ · d_h) >> 1 where τ and d_h representing the duct’s centerline torsion and its hydraulic diameter. This approximation and the use of an orthogonal grid allow a great simplification on the numerical procedure. Comparisons of the numerical solution against experimental data are drawn to assess the accuracy of the approximation.     

A survey of dynamically-adaptive grids in the numerical solution of partial differential equations

The construction of dynamically-adaptive curvilinear coordinate systems based on numerical grid generation and the use thereof in the numerical solution of partial differential equations is surveyed, and correlations are made among the various approaches. These adaptive grids are coupled with the physical solution being done on the grid so that the grid points continually move in the course of the solution in order to resolve developing gradients, or higher variations, in the solution. Particular attention is given to systems using elliptic grid generation based on variational principles. It is noted that dynamic grid adaption can remove the oscillations common when strong gradients occur on fixed grids, and that it appears that when the grid adapts to the solution most numerical solution algorithms work well. Particular applications in computational fluid dynamics and heat transfer are noted.     

Generation of structured difference grids in two-dimensional nonconvex domains using mappings

Generation of structured difference grids in two-dimensional nonconvex domains is considered using a mapping of a parametric domain with a given nondegenerate grid onto a physical domain. For that purpose, a harmonic mapping is first used, which is a diffeomorphism under certain conditions due to Rado’s theorem. Although the harmonic mapping is a diffeomorphism, its discrete implementation can produce degenerate grids in nonconvex domains with highly curved boundaries. It is shown that the degeneration occurs due to approximation errors. To control the coordinate lines of the grid, an additional mapping is used and universal elliptic differential equations are solved. This makes it possible to generate a nondegenerate grid with cells of a prescribed shape.     

Well-posedness analysis and numerical implementation of a linearized two-dimensional bottom sediment transport problem

A two-dimensional linearized model of coastal sediment transport due to the action of waves is studied. Up till now, one-dimensional sediment transport models have been used. The model under study makes allowance for complicated bottom relief, the porosity of the bottom sediment, the size and density of sediment particles, gravity, wave-generated shear stress, and other factors. For the corresponding initial–boundary value problem the uniqueness of a solution is proved, and an a priori estimate for the solution norm is obtained depending on integral estimates of the right-hand side, boundary conditions, and the norm of the initial condition. A conservative difference scheme with weights is constructed that approximates the continuous initial–boundary value problem. Sufficient conditions for the stability of the scheme, which impose constraints on its time step, are given. Numerical experiments for test problems of bottom sediment transport and bottom relief transformation are performed. The numerical results agree with actual physical experiments.     

Saturated-unsaturated groundwater modeling using 3D Richards equation with a coordinate transform of nonorthogonal grids

Saturated-unsaturated flow under a complex terrain is usually solved using the Richards equation. Finite difference or finite volume methods are commonly employed for discretization of Richards equation because of simplicity of coding. Complex subsurface boundary geometries lead to nonorthogonal grids in curvilinear grid systems, which leads to difficulty in discretization and mesh generation. This paper develops a vertical coordinate transform, enabling a computational domain regular in the vertical direction. As a result, the grid of curvilinear surfaces can be successfully transformed to a computational grid that allows solution of the Richards equation with efficient computation and simpler coding. The anisotropic Richards equation in the Cartesian coordinate system is transformed to the equation in the arbitrary coordinate system and further expressed as a form appropriate to the orthogonal coordinate system. The generalized third boundary condition is transformed to a form suited to the orthogonal coordinate system. The finite volume method is used to solve the Richards equation in the orthogonal coordinate system. Four cases are used to validate the present orthogonal coordinate system. The computational results from the orthogonal coordinate system are in good agreement with the results from Ansys Fluent solved in a Cartesian coordinate system for the subsurface flow case. For the coupled case of hill slopes, a good agreement between the computational results and the experimental data is obtained. The present results for V-titled catchment and slab case accord well with the results obtained from HydroGeoSphere and PAWS. The present algorithm can improve grid generation for solution of Richards equation in a hydrological model for a complex domain.     

On grid generation for numerical models of geophysical fluid dynamics

A simple geometric condition that defines the class of classical (stereographic, conic and cylindrical) conformal mappings from a sphere onto a plane is derived. The problem of optimization of computational grid for spherical domains is solved in an entire class of conformal mappings on spherical (geodesic) disk. The characteristics of computational grids of classical mappings are compared for different spherical radii of geodesic disk. For a rectangular computational domain, the optimization problem is solved in the class of classical mappings and respective area of the spherical domain is evaluated.     

Finite difference simulation of elastic wave propagation through 3D heterogeneous multiscale media based on locally refined grids

To simulate the interaction of seismic waves with microheterogeneities (like cavernous/fractured reservoirs), a finite difference technique based on grids locally refined in time and space is used. These grids are used because the scales of heterogeneities in the reference medium and in the reservoir are different. Parallel computations based on domain decomposition of the target area into elementary subdomains in both the reference medium (a coarse grid) and the reservoir (a fine grid) are performed. Each subdomain is assigned to a specific processor unit, which forms two groups: one for the reference medium, and the other for the reservoir. The data exchange between the groups within a processor unit is performed by non-blocking iSend/iReceive MPI commands. The data exchange between the two groups is performed simultaneously with coupling the coarse and a fine grids, and is controlled by a specially chosen processor unit. The results of a numerical simulation for a realistic model of fracture corridors are presented and discussed.     

Repeated Richardson extrapolation applied to the two-dimensional Laplace equation using triangular and square grids

The focus of this work is to verify the efficiency of the Repeated Richardson Extrapolation (RRE) to reduce the discretization error in a triangular grid and to compare the result to the one obtained for a square grid for the two-dimensional Laplace equation. Two different geometries were employed: the first one, a unitary square domain, was discretized into a square or triangular grid; and the second, a half square triangle, was discretized into a triangular grid. The methodology employed used the following conditions: the finite volume method, uniform grids, second-order accurate approximations, several variables of interest, Dirichlet boundary conditions, grids with up to 16,777,216 nodes for the square domain and up to 2097,152 nodes for the half square triangle domain, multigrid method, double precision, up to eleven Richardson extrapolations for the first domain and up to ten Richardson extrapolations for the second domain. It was verified that (1) RRE is efficient in reducing the discretization error in a triangular grid, achieving an effective order of approximately 11 for all the variables of interest for the first geometry; (2) for the same number of nodes and with or without RRE, the discretization error is smaller in a square grid than in a triangular grid; and (3) the magnitude of the numerical error reduction depends on, among other factors, the variable of interest and the domain geometry.