A Preliminary Look at Methods to Compute Baroclinic Pressure Gradients in FE Models

Abstract

In shelf regions with steep bathymetry in the presence of density gradients, the computation of the baroclinic pressure gradient term in 3D shallow water models may either become unstable or physically unrealistic. This manuscript examines four algorithms to compute the baroclinic pressure gradient term in finite element models. Two common systems for discretizing the vertical are sigma or zlevel coordinates. In turn, permutations of these two coordinate systems serve as the basis for the four different algorithms examined herein. All are implemented in the context of the finite element hydrodynamic model, ADCIRC. Several density gradients that vary horizontally and vertically, with the pycnocline occurring at different depths, are used to evaluate the methods in a three-dimensional box grid with variable bathymetry. Initial testing, the subject of this work, focuses on model behavior for simplified problem with known analytic solutions. Long term goals for this study (the subject of subsequent papers) are three-fold: 1) to determine the vertical node placement algorithm that produces the most stable and physically realistic results; 2) to determine the interplay of vertical and horizontal resolution (and bathymetry and density profile), while considering simulation time; and 3) to produce accurate 3D flow fields and baroclinic pressure gradients. Outcomes from the study will be used to direct on-going modeling research in the Mississippi Sound. Introduction Background An accurate prediction of tides and circulation patterns using hydrodynamic models is useful in a wide variety of civilian and military applications that range 1. School of Civil Engineering and Environmental Science, 202 W. Boyd, Room 334, University of Oklahoma, Norman, OK, 73019, USA, dresback@ou.edu, kolar@ou.edu 2. Oceanography Division, Code 7322, Naval Research Laboratory, Stennis Space Center, MS, 39529, USA, blain@nrlssc.navy.mil 2 Dresback, Blain and Kolar from forecasting storm surges from extreme events, such as tsunamis and hurricanes, in order to plan development in coastal areas to predicting circulation patterns in order to guide fleet operations. Moreover, the accuracy of water quality models depends heavily on an accurate, robust hydrodynamic model. One class of models that finds frequent use for these and other applications are the so-called shallow water models. Shallow water equations in 2D are obtained by vertically averaging the mass and momentum balances over the depth of the fluid. For our applications, we solve the equations using algorithms amenable to irregular triangular meshes (e.g., finite element, finite volume). Such a meshing technique allows grids to be easily refined in shallow coastal areas or in areas with changing topography where high resolution is needed to accurately resolve tidal constituents and circulation patterns. Also irregular coastal boundaries are more accurately described with triangular elements, and flux-type boundary conditions enter the problem naturally. Models in this class include ADCIRC [16], QUODDY [22], and UTBEST [5]. Finite difference based models, such as SIMSYS2D[15], CH3D [4], and POM [30] that use a staggered grid on an orthogonal curvilinear mesh, are stable and mass conservative, but their grids are not as inherently flexible as those based on irregular triangular meshes [15, 29]. Early finite element solutions of the shallow water equations were often plagued by spurious oscillation modes. Various methods were introduced to eliminate the oscillations, but most include some type of artificial or numerical damping [27, 31]. Lynch and Gray [21] present the wave continuity equation as a means to eliminate spurious oscillations without resorting to artificial damping. It is currently the base algorithm in ADCIRC (ADvanced 3D CIRCulation model [16]), the model that is the subject of this research. Since the inception of the wave continuity formulation in 1979, the original algorithm has been modified in a number of substantial ways: Kinnmark introduced a numerical parameter to provide a general means of weighting the continuity equation [11]; viscous dissipation terms were incorporated [12, 20]; and 3D barotropic simulations were realized by resolving the velocity profile in the vertical [16, 25]. Most recently, the model has incorporated the following features: a wetting/drying algorithm for near-shore elements [17]; spherical coordinates for large domains via a map projection technique [13]; a parallel domain decomposition algorithm for the 2D portion of the code [26]; a fully implicit time marching algorithm [6]; options for hydraulic structures (e.g., levees, breakwaters); and a diagnostic baroclinic algorithm [18]. These advances are all present in the current version of the ADCIRC model, a model that represents over 20 years of research and development by scientists and engineers at six different institutions. ADCIRC has been extensively tested using analytical solutions and field data, and it is being used by a diverse group of researchers and practitioners to model the hydrodynamic behavior of many coastal and oceanic regions [2, 7, 8, 10, 23, 34].