Nonhydrostatic Ocean Model Research Articles

A Numerical Simulation of Internal Wave Propagation on a Continental Slope and Its Influence on Sediment Transport

Significant current velocity near the sea bottom can be induced by internal waves, even for water a few hundred meters in depth. In this study, a nonhydrostatic ocean model was applied to simulate the generation and propagation of internal waves on the continental slope of the northern SCS. Based on the analyses of the vertical profiles of the currents, the propagation of internal waves along the continental slope can be categorized into six modes. The bed shear stress and the bedload transport were calculated to analyze the general characteristics of sediment transport along the continental slope of the northern SCS. Generally, there was no sediment transport on the sea bottom induced by the internal waves when the water depth was deeper than 650 m or shallower than 80 m. The downslope sediment transport dominated the slope at a water depth range of 200~650 m, while the upslope sediment transport dominated the slope at a water depth range of 80~200 m. The predicted directions of the bedload transport are coincident with the field observations of sand wave migration on the continental slope, which further confirms that the main cause of the generation and formation of sand waves on the continental slope of the northern SCS is the strong bottom current induced by the shoaling process of internal waves.

Nonhydrostatic simulation of internal lee wave-induced turbulence mixing: A comparative study of second-moment closure mixing schemes for an idealized marginal sea sill topography

To further understand the processes of wave-turbulence interactions under nonhydrostatic conditions and to gain insight into the selection of second-moment closure models in subsequent nonhydrostatic ocean models, a newly developed nonhydrostatic ocean model is coupled with the General Ocean Turbulence Model (GOTM) to form a robust version of the Marine Environment Research and Forecasting Model (MERF). On this basis, internal lee wave simulation with various turbulent mixing schemes is investigated for an idealized two-dimensional supercritical sea slope topography with high-resolution (100 m). First, the LES-type simulation based on a 3D Smagorinsky diffusivity scheme can reproduce well the characteristics of lee waves, serving as a benchmark for comparison. Second, the Mellor and Yamada (MY) scheme and the quasi-equilibrium scheme are able to effectively simulate the effects of internal waves on large-scale structures, while small-scale dynamic processes are partially suppressed by turbulent mixing. However, the K-ε model and the one-equation model further suppress the small-scale dynamic processes and weaken the propagation of internal lee waves due to overly strong turbulence mixing. Finally, the generation of a lee wave changes the instantaneous temperature field, and the use of the second-moment closure models converts more energy into irreversible mixing. Enhancing the model resolution can slightly improve the accuracy of the simulation in exchange for paying considerable computation.

A high-order spectral method for effective simulation of surface waves interacting with an internal wave of large amplitude

In this study, we propose a new method for simulating complex surface waves interacting with a large amplitude internal solitary wave. Our model is based on a high-order spectral method for surface waves with the bottom boundary conditions computed from an internal wave solver. The convergence of the model is first tested using an internal wave case without surface waves. We then perform simulations of a canonical nonlinear wave interaction case involving two surface wave components and a prescribed internal wave component. The initial wave energy growth rate of the resonant wave component is found to agree with the analytical value predicted by the perturbation theory. We also use the model to simulate the interaction between surface waves and a weakly nonlinear internal wave, and the result is nearly identical to that calculated using a two-layer model. Finally, we show the application of our model in a multiscale nested modeling framework, where the internal wave parameters are extracted from the mesoscale simulation using a nonhydrostatic ocean model. The evolution of the spatial variations of the surface roughness is captured and the surface wave orbital velocity also changes in space due to the surface motions induced by the internal wave. Our phase-resolved model provides a computationally efficient tool for simulating complex surface wave fields in the background of internal waves of large amplitudes.

MERF v3.0, a highly computationally efficient non-hydrostatic ocean model with implicit parallelism: Algorithms and validation experiments

Non-hydrostatic ocean models are widely required in coastal and estuarine studies but are limited in their application due to their very high computational expense. To alleviate this issue, this study developed a highly computationally efficient non-hydrostatic ocean model named MERF v3.0 (the third version of the Marine Environment Research and Forecasting model) in which a semi-implicit and variable layers (SIVL) scheme is proposed for dynamic pressure calculation. The SIVL scheme provides freedom to choose the ratio between the number of dynamic pressure layers and the number of velocity layers as required. Moreover, the model adopted numerical calculation framework OpenArray and PETSc (Portable, Extensible Toolkit for Scientific Computation) to achieve implicit parallelism and a concise code structure. Five classical numerical experiments were performed to validate the simulation capacity of this non-hydrostatic ocean model, including the standing wave in a closed basin, surface solitary wave, lock-exchange problem, periodic wave over a bar and tidally induced internal lee wave. The numerical simulation results are consistent with the analytical solutions and experimental data published in the literature. Benefitting from the advantages of OpenArray and PETSc, MERF v3.0 achieves an 87.8% parallel efficiency in the strong scaling test and an 84.1% parallel efficiency in the weak scaling test on 3584 processes, whereas the OpenArray limits the parallel efficiency of MERF v3.0 on more processes, which is beyond the scope of this article and will be further examined in our next study.

A simulation study on effects of suspended sediment through high riverine discharge on surface river plume and vertical water exchange

Rivers transport freshwater and suspended sediment matter (SSM) from land to coastal seas. In coastal seas termed as regions of freshwater influence (ROFIs), SSM is not only passively transported but also changes the density of ambient water and influences the physical characteristics especially in flood events, when a lot of SSM is supplied to the sea. Although the influence of SSM on the physical field in ROFIs would be significant, interactive physical processes, such as dynamics of river plumes and estuarine circulations, have hardly been investigated for hypopycnal plumes (i.e., the riverine sediment-freshwater is not denser than the seawater). In order to quantitatively estimate the interactive effects of SSM, we employ a non-hydrostatic ocean model with Lagrangian particles, which represents SSM and affects the density and buoyancy of ambient water. We use two experimental settings: (1) realistic simulations of the Tango Bay, Japan under the flooding of September 2013 and (2) idealized simulations for an open-bay ROFI. The former is conducted to assess whether the simulations could reproduce an actual event to some extent. The realistic simulations demonstrate that the choice of parameters such as SSM-particle size and composition is important for coastal simulations of flood events. The latter is conducted to understand the basic physics and to study the quantitative sensitivities of the physical processes to the riverine flux, composition, and particle size of SSM. The idealized simulations demonstrate that a large amount of riverine SSM affects the physical field in ROFIs through the following process: 1) horizontal density differences between nearshore and offshore waters are reduced as apparent density is increased close to the river mouth by riverine SSM, 2) the strength of vertical circulation is weakened by the reduced horizontal density difference, and 3) vertical water exchange between the surface and the subsurface layers decreases. The process in the control case of this study increases the relative amounts of surface freshwater in the river plume by 0.8% of the total riverine freshwater input. Sensitivity experiments with changing the parameters of SSM flux into the river, SSM composition, and SSM-particle diameter show that the percentage can be raised up to 2%. Meanwhile, the above-mentioned processes do not apply to extreme cases of small particle size and enormous SSM input wherein homopycnal (i.e., a similar density between the riverine sediment-freshwater and the ambient seawater) and hyperpycnal (i.e., the density of the riverine sediment-freshwater mixture exceeds the ambient seawater) plumes take place.

Monitoring and Modelling Bio-Physical Parameters for Hypoxia Hazard in a Coastal Sand Pit

Management of coastal areas requires monitoring and modeling of the anthropogenic drivers and the bio-physical processes affecting water quality. To assess the range of hydrographic conditions controlling oxygen distribution in the bottom layers of sand pits, a multi-year oceanographic survey has been conducted in a coastal area with several extraction pits. Hydrographic data including profiles of temperature, salinity and oxygen were collected and related to local wind conditions and circulation. Moreover, 1D and 3D high-resolution non-hydrostatic ocean models were used to describe turbulent mixing regimes and to obtain the range of wind speeds for which the critical anoxic conditions may occur. It is shown that wind speed appears to control the dynamics of oxygen concentrations, with oxygen depleted zones developing in a short time in low wind speed conditions. Moreover, the depth and the shape of the extraction pit contribute to decrease the mixing of the bottom layers and increase the water retention in the hole increasing the output and the persistence of oxygen depleted zones in the excavated area. The results of the numerical simulations show that the risk of hypoxia at the bottom of the sand pits is associated with higher temperatures and wind speed lower than 5 m/s, which is not infrequent during the summer season. However, the number of consecutive days of oxygen depletion can be considered lower than the danger threshold level assumed in the literature.

Partially implicit finite difference scheme for calculating dynamic pressure in a terrain-following coordinate non-hydrostatic ocean model

To improve the efficiency of the terrain-following σ-coordinate non-hydrostatic ocean model, a partially implicit finite difference (PIFD) scheme is proposed. By using explicit terms instead of implicit terms to discretize the parts of the vertical dynamic pressure gradient derived from the σ-coordinate transformation, the coefficient matrix of the discrete Poisson equation that the dynamic pressure satisfies can be simplified from 15 diagonals to 7 diagonals. The PIFD scheme is shown to run stably when it is applied to simulate five benchmark cases, namely, a standing wave in a basin, a surface solitary wave, a lock-exchange problem, a periodic wave over a bar and a tidally induced internal wave. Compared with the conventional fully implicit finite difference (FIFD) scheme, the PIFD scheme produces simulation results of equivalent accuracy at only 40–60% of the computational cost. The PIFD scheme demonstrates strong applicability and can be easily implemented in σ-coordinate ocean models.

A 3D unstructured non-hydrostatic ocean model for internal waves

A 3D non-hydrostatic model is developed to compute internal waves. A novel grid arrangement is incorporated in the model. This not only ensures the homogenous Dirichlet boundary condition for the non-hydrostatic pressure can be precisely and easily imposed but also renders the model relatively simple in its discretized form. The Perot scheme is employed to discretize horizontal advection terms in the horizontal momentum equations, which is based on staggered grids and has the conservative property. Based on previous water wave models, the main works of the present paper are to (1) utilize a semi-implicit, fractional step algorithm to solve the Navier-Stokes equations (NSE); (2) develop a second-order flux-limiter method satisfying the max–min property; (3) incorporate a density equation, which is solved by a high-resolution finite volume method ensuring mass conservation and max–min property based on a vertical boundary-fitted coordinate system; and (4) validate the developed model by using four tests including two internal seiche waves, lock-exchange flow, and internal solitary wave breaking. Comparisons of numerical results with analytical solutions or experimental data or other model results show reasonably good agreement, demonstrating the model’s capability to resolve internal waves relating to complex non-hydrostatic phenomena.

M2, overtides and compound tides generation in the Strait of Messina: the response of a non-hydrostatic, finite-element ocean model

ABSTRACT Quattrocchi, A.; Gorman, J.G., Piggott M.D. and Cucco, A., 2016. M2, overtides and compound tides generation in the Strait of Messina: the response of a non-hydrostatic, finite-element ocean model. In: Vila-Concejo, A.; Bruce, E.; Kennedy, D.M., and McCarroll, R.J. (eds.), Proceedings of the 14th International Coastal Symposium (Sydney, Australia). Journal of Coastal Research, Special Issue, No. 75, pp. 657–661. Coconut Creek (Florida), ISSN 0749-0208. The nonlinear interactions between the main tidal constituents and geographical constraints are investigated using Fluidity-ICOM, a non-hydrostatic finite-element ocean model that was implemented in the region of the Strait of Messina (Mediterranean Sea). As a first crucial approach, the barotropic tidal dynamics was investigated and the model solution was interpreted in terms of tidal maps and spectra obtained by harmonic analysis of the sea surface elevation. The model domain is represented by an unstructured mesh with variable spatial resolution...

GPU Acceleration of a Non-hydrostatic Ocean Model with a Multigrid Poisson/Helmholtz solver

To meet the demand for fast and detailed calculations in numerical ocean simulations, we implemented a non-hydrostatic ocean model on a graphics processing unit (GPU). We improved the model's Poisson/Helmholtz solver by optimizing the memory access, using instruction-level parallelism, and applying a mixed precision calculation to the preconditioning of the Poisson/Helmholtz solver. The GPU-implemented model was 4.7 times faster than a comparable central processing unit execution. The output errors due to this implementation will not significantly influence oceanic studies.

A Lagrange-Galerkin hp-Finite Element Method for a 3D Nonhydrostatic Ocean Model

We introduce in this paper a Lagrange-Galerkin hp-finite element method to calculate the numerical solution of a nonhydrostatic ocean model. The Lagrange-Galerkin method yields a Stokes-like problem the solution of which is computed by a second-order rotational splitting scheme that separates the calculation of the velocity and pressure, the latter is decomposed into hydrostatic and nonhydrostatic components. We have tested the method in flows where the nonhydrostatic effects are important. The results are very encouraging.

About interface conditions for coupling hydrostatic and nonhydrostatic Navier-Stokes flows

In this work we are interested in the search of interface conditions to couple hydrostatic and nonhydrostatic ocean models. To this aim, we consider simplified systems and use a time discretization to handle linear equations. We recall the links between the two models (with the particular role of the aspect ratio δ = H/L 1) and introduce an iterative method based on the Schwarz algorithm (widely used in domain decomposition methods). The convergence of this method depends strongly on the choice of interface conditions: this is why we look for exact absorbing conditions and their approximations in order to provide tractable and efficient coupling algorithms.

On the Conditional Frazil Ice Instability in Seawater

AbstractIt has been suggested that the presence of frazil ice can lead to a conditional instability in seawater. Any frazil forming in the water column reduces the bulk density of a parcel of frazil–seawater mixture, causing it to rise. As a result of the pressure decrease in the freezing point, this causes more frazil to form, causing the parcel to accelerate, and so on. This study uses linear stability analysis and a nonhydrostatic ocean model to study this instability. The authors find that frazil ice growth caused by the rising of supercooled water is indeed able to generate a buoyancy-driven instability. Even in a gravitationally stable water column, the frazil ice mechanism can still generate convection. The instability does not operate in the presence of strong density stratification, high thermal driving (warm water), a small initial perturbation, high background mixing, or the prevalence of large frazil ice crystals. In an unstable water column, the instability is not necessarily expressed in frazil ice at all times; an initial frazil perturbation may melt and refreeze. Given a large enough initial perturbation, this instability can allow significant ice growth. A model shows frazil ice growth in an Ice Shelf Water plume several kilometers from an ice shelf, under similar conditions to observations of frazil ice growth under sea ice. The presence of this instability could be a factor affecting the growth of sea ice near ice shelves, with implications for Antarctic Bottom Water formation.

The Effect of Meltwater Plumes on the Melting of a Vertical Glacier Face

AbstractFreshwater produced by the surface melting of ice sheets is commonly discharged into ocean fjords from the bottom of deep fjord-terminating glaciers. The discharge of the freshwater forms upwelling plumes in front of the glacier calving face. This study simulates the meltwater plumes emanated into an unstratified environment using a nonhydrostatic ocean model with an unstructured mesh and subgrid-scale mixing calibrated by comparison to established plume theory. The presence of an ice face reduces the entrainment of seawater into the meltwater plumes, so the plumes remain attached to the ice front, in contrast to previous simple models. Ice melting increases with height above the discharge, also in contrast to some simple models, and the authors speculate that this “overcutting” may contribute to the tendency of icebergs to topple inwards toward the ice face upon calving. The overall melt rate is found to increase with discharge flux only up to a critical value, which depends on the channel size. The melt rate is not a simple function of the subglacial discharge flux, as assumed by many previous studies. For a given discharge flux, the geometry of the plume source also significantly affects the melting, with higher melt rates obtained for a thinner, wider source. In a wider channel, two plumes are emanated near the source and these plumes eventually coalesce. Such merged meltwater plumes ascend faster and increase the maximum melt rate near the center of the channel. The melt rate per unit discharge decreases as the subglacial system becomes more channelized.

A nonhydrostatic, isopycnal-coordinate ocean model for internal waves

We present a nonhydrostatic ocean model with an isopycnal (density-following) vertical coordinate system. The primary motivation for the model is the proper treatment of nonhydrostatic dispersion and the formation of nonlinear internal solitary waves. The nonhydrostatic, isopycnal-coordinate formulation may be preferable to nonhydrostatic formulations in z- and σ-coordinates because it improves computational efficiency by reducing the number of vertical grid points and eliminates spurious diapycnal mixing and solitary-wave amplitude loss due to numerical diffusion of scalars. The model equations invoke a mild isopycnal-slope approximation to remove small metric terms associated with diffusion and nonhydrostatic pressure from the momentum equations and to reduce the pressure Poisson equation to a symmetric linear system. Avoiding this approximation requires a costlier inversion of a non-symmetric linear system. We demonstrate that the model is capable of simulating nonlinear internal solitary waves for simplified and physically-realistic ocean-scale problems with a reduced number of layers.

Physical vs. numerical dispersion in nonhydrostatic ocean modeling

Many large-scale simulations of internal waves are computed with ocean models solving the primitive (hydrostatic) equations. Under certain circumstances, however, internal waves can represent a dynamical balance between nonlinearity and nonhydrostasy (dispersion), and thus may require computationally expensive nonhydrostatic simulations to be well-resolved. Most discretizations of the primitive equations are second-order accurate, inducing numerical dispersion generated from odd-order terms in the truncation error (3rd-order derivatives and higher). This numerical dispersion mimics physical dispersion due to nonhydrostasy. In this paper, we determine the numerical dispersion coefficient associated with common discretizations of the primitive equations. We compare this coefficient to the physical dispersion coefficient from the Boussinesq equations or KdV equation. The results show that, to lowest order, the ratio of numerical to physical dispersion is Γ = Kλ 2, where K is an O(1) constant dependent on the discretization of the governing equations and λ is the grid leptic ratio, λ ≡ Δ x/ h 1, where Δ x is the horizontal grid spacing and h 1 is the depth of the internal interface. In addition to deriving this relationship, we verify that it indeed holds in a nonhydrostatic ocean model (SUNTANS). To ensure relative dominance of physical over numerical effects, simulations require Γ ≪ 1. Based on this condition, the horizontal grid spacing required for proper resolution of nonhydrostatic effects is λ < O(1) or Δ x < h 1. When this condition is not satisfied, numerical dispersion overwhelms physical dispersion, and modeled internal waves exist with a dynamical balance between nonlinearity and numerical dispersion. Satisfaction of this condition may be a significant additional resolution requirement beyond the current state-of-the-art in ocean modeling.

Solving the Poisson equation on small aspect ratio domains using unstructured meshes

We discuss the ill-conditioning of the matrix for the discretised Poisson equation in the small aspect ratio limit, and motivate this problem in the context of nonhydrostatic ocean modelling. Efficient iterative solvers for the Poisson equation in small aspect ratio domains are crucial for the successful development of nonhydrostatic ocean models on unstructured meshes. We introduce a new multigrid preconditioner for the Poisson problem which can be used with finite element discretisations on general unstructured meshes; this preconditioner is motivated by the fact that the Poisson problem has a condition number which is independent of aspect ratio when Dirichlet boundary conditions are imposed on the top surface of the domain. This leads to the first level in an algebraic multigrid solver (which can be extended by further conventional algebraic multigrid stages), and an additive smoother. We illustrate the method with numerical tests on unstructured meshes, which show that the preconditioner makes a dramatic improvement on a more standard multigrid preconditioning approach, and also show that the additive smoother produces better results than standard SOR smoothing. This new solver method makes it feasible to run nonhydrostatic unstructured mesh ocean models in small aspect ratio domains.

The surface boundary condition in nonhydrostatic ocean models

With increasing resolution in numerical ocean models, nonhydrostatic pressure effects have to be accounted for. In sigma-coordinate mode split ocean models, this pressure may be regarded as a pressure correction. An elliptic equation must be solved for the nonhydrostatic pressure, and the gradients are used to correct the provisional hydrostatic velocity components in each time step. The focus in the present work is on the surface boundary condition for the elliptic equation. In the literature, both Dirichlet and Neumann boundary conditions are suggested and applied. To investigate the sensitivity of the numerical results to the choice of boundary condition, three numerical experiments are performed. The first and second experiments are studies of the propagation and steepening of nonlinear internal waves. The first study is on tank scale and the second experiment is on ocean scale. In the tank-scale experiment, the density and the flow fields are very robust to the choice of boundary condition. In the ocean-scale experiment, the waves produced with a Dirichlet boundary condition become more damped than the waves produced with a Neumann boundary condition. The third study involves a surface buoyant jet. It is shown that well-known characteristics of the plume front are reproduced with a Neumann boundary condition, but the rotating turbulent core of this front is lost with a Dirichlet condition. It is accordingly argued that the appropriate surface boundary condition in mode split nonhydrostatic ocean models is the Neumann condition.

Effects on acoustics caused by ocean solitons, Part B: Acoustics

Large amplitude internal solitary waves in the ocean can interfere with underwater acoustic signals. For certain acoustic parameters (source depth, receiver depth, and frequency) a redistribution of acoustical energy to higher-order acoustic modes can occur as the acoustic signal propagates through the solitary wave train. Depending on the ocean bottom composition, this can result in a significant loss in acoustical signal intensity. In order to simulate and understand this phenomena, it is necessary to generate a solitary wave train that is realistic and contains the spectral components that will interact with the acoustic signal. This requires a primitive equation nonhydrostatic ocean model. Part A [A.C. Warn-Varnas, S.A. Chin-Bing, D. King, J. Hawkins, K. Lamb, Effects on acoustics caused by ocean solitons: Part A. Oceanography, Nonlinear Anal., in press ( doi:10.1016/j.na.2009.02.104)] of this two-part article addressed that issue. In this work we use acoustic model simulations to demonstrate the mode linkage between the solitary wave train and the acoustic signal that can result in a large signal loss. Understanding the linkage allows for pre-examination of the solitary wave characteristics to determine if they could affect propagating acoustic signals. The procedure has the advantage that only a few computer runs are required from the ocean model and no lengthy computer runs are required from acoustic propagation models. Only a quick calculation of the acoustic mode eigenvalues is needed.

Non-hydrostatic pressure in σ-coordinate ocean models

Numerical ocean modelling is computationally very demanding. Traditionally, the hydrostatic approximation has been applied to reduce the computational burden. This approximation is valid in large scale studies with coarse grid resolution. With faster computers and gradually smaller grid sizes, we may expect that more studies will be performed with non-hydrostatic ocean models. In recent papers several methods for including non-hydrostatic pressure in σ-coordinate models have been suggested. In this paper the sensitivity of the non-hydrostatic pressure field, the velocity fields, and the density fields to changes in the method for computing non-hydrostatic pressure in σ-coordinate ocean models is addressed. The first test case used involves the propagation and breaking of an internal wave at an incline in a tank. The other test case concerns tidally driven flow over a sill in a stratified fjord. The results from our numerical exercises suggest that the velocity and density fields are very robust to the model choices investigated here. The differences between the model results are of the same order as the uncertainty due to the internal pressure gradient error, and they are smaller than an estimate of the uncertainty due to subgrid scale closure.