USGS -- SMIG --
Surface-water quality and flow Modeling Interest Group

Modeling Surface Trapped River Plumes: A Sensitivity Study

by Jason Hyatt1 and Richard P. Signell1

1 U.S. Geological Survey, Woods Hole, MA, U.S.A.

Please direct correspondence to:
  Jason Hyatt
  USGS, Woods Hole Field Center
  384 Woods Hole Road
  Woods Hole, MA 02543-1598
  Internet: jhyatt@usgs.gov
  Phone: (508) 457-2330
  FAX: (508) 457-2309

Editor's note:
This article was written for the 6th International Conference on Estuarine and Coastal Modeling, held November 3-5, 1999 in New Orleans, LA. The document available here is based on the final draft provided to the editors. Minor discrepancies between this document and the published version, therefore, may exist.

This version of the article has all of the figures converted to thumbnails with links to the larger images. A version with all of the figures inline is also available; the download time may be longer, but the inline figures may be more convenient for viewing and printing.

This paper has links to animations in AVI format. If you don't have a movie player that allows fine control over frame advancement and reverse, consider obtaining one. For the Mac and Windows environments, try the free QuickTime viewer from Apple. For Unix systems, try xanim (also free).

Citation:
Hyatt, J. and Signell, R.P., 2000, Modeling surface trapped river plumes: A sensitivity study, in Estuarine and Coastal Modeling, 6th Int. Conf., ASCE, New Orleans, LA, November 3-5, 1999. Editors: Malcolm L. Spaulding and Alan F. Blumberg.

Contents

Abstract

To better understand the requirements for realistic regional simulation of river plumes in the Gulf of Maine, we test the sensitivity of the Blumberg-Mellor hydrodynamic model to choice of advection scheme, grid resolution, and wind, using idealized geometry and forcing. The test case discharges 1500 m3/s of fresh water into a uniform 32 psu ocean along a straight shelf at 43o north. The water depth is 15 m at the coast and increases linearly to 190 m at a distance 100 km offshore. Constant discharge runs are conducted in the presence of ambient alongshore current and with and without periodic alongshore wind forcing. Advection methods tested are CENTRAL, UPWIND, the standard Smolarkiewicz MPDATA and a recursive MPDATA scheme. For the no-wind runs, the UPWIND advection scheme performs poorly for grid resolutions typically used in regional simulations (grid spacing of 1-2 km, comparable to or slightly less than the internal Rossby radius, and vertical resolution of 10% of the water column), damping out much of the plume structure. The CENTRAL difference scheme also has problems when wind forcing is neglected, and generates too much structure, shedding eddies of numerical origin. When a weak 5 cm/s ambient current is present in the no-wind case, both the CENTRAL and standard MPDATA schemes produce a false fresh- and dense-water source just upstream of the river inflow due to a standing two-grid length oscillation in the salinity field. The recursive MPDATA scheme completely eliminates the false dense water source, and produces results closest to the grid-converged solution. The results are shown to be very sensitive to vertical grid resolution, and the presence of wind forcing dramatically changes the nature of the plume simulations. The implication of these idealized tests for realistic simulations is discussed, as well as ramifications on previous studies of idealized plume models.

Introduction

Surface trapped river plumes are important features around the world, carrying freshwater, nutrients, and pollutants into the coastal ocean. In the western Gulf of Maine, river plumes have been linked with the along-coast delivery of red tide organisms (Franks and Anderson, 1992), making prediction of plume behavior a major regional goal. Realistic simulation of surface trapped river plumes presents a particular challenge, because plumes respond rapidly to wind and baroclinic forcing, and often exhibit strong vertical and horizontal density gradients. Numerical calculations of river plume dynamics are therefore particularly sensitive to the advection scheme, grid resolution, and wind forcing used.

Many previous investigations have studied idealized plume responses. Chao and Boicourt (1986) studied the onset of the plume on a flat-bottomed coastal ocean without ambient current, observing a bulge of anticyclonic surface flow at the estuary mouth and a bore intrusion along the shelf. Chao (1988) added a sloping bottom and classified the plumes into four regimes, supercritical, subcritical, diffusive-supercritical, and diffusive sub-critical. The plume in this study is diffusive-subcritical. Oey and Mellor (1992) used a flat bottom domain to examine the unsteady aspects of the plume and front system, noting pulsating detached pools. Kourafalou et al. (1996) also observed the bulge and bore and included wind in the their calculations over a flat bottom domain. Table 1 summarizes the geometry of these studies. Most of these studies resolved the internal Rossby radius with just a few grid cells in the horizontal and the thickness of the river plume with just a few grid cells in the vertical dimension. It therefore seemed a useful contribution would be to examine the role of advection scheme as well as vertical and horizontal resolution on idealized plume cases under both no-wind and wind-forced conditions.

Table 1. Summary of some previous idealized plume studies geometry and resolution

table 1

Modeling Framework

We used the Estuary, Coastal, and Ocean Model (ECOM-3D), the HydroQual, Inc. version of the Blumberg and Mellor sigma-coordinate model (1987). This version allows for easy switching between four advection schemes, including CENTRAL, UPWIND, and the two forms of multi-dimensional positive definite advective transport algorithm (MPDATA) (Smolarkiewicz, 1984, Smolarkiewicz and Clark, 1986). Each of these advection schemes includes advantages and disadvantages, especially in the presence of sharp gradients in scalar variables. The CENTRAL scheme has a long history of use in a variety of domains and is second order accurate with no numerical diffusion (Roach, 1976). However, it introduces numerical dispersion resulting in ripples and subsequent under- and over-shoot of initial values. The UPWIND scheme has no dispersion and is positive definite, but introduces strong numerical diffusion as a result of the first-order truncation. MPDATA is an iterative scheme that attempts to correct the diffusion errors in the UPWIND scheme. After an initial calculation of the scalar field using the UPWIND scheme, an "anti-diffusion" velocity field based on the local first-order truncation error is generated. The UPWIND scheme is then implemented with this velocity field to "bring back" the scalar field that has spread too far during the time step due to the numerical diffusion. This restorative step may be exercised a number of times, which depends upon the users discretion and computational resources. Here SMOLAR_2 refers to the scheme that makes exactly two corrective sweeps for numerical diffusion.

A significant advance was made by Margolin and Smolarkiewicz (1989) who realized that a recursion relation for the antidiffusion velocities could be determined analytically, allowing the benefits of many iterations to be realized in a single correction step (albeit with a more complex determination of the antidiffusion velocities). Here we term this scheme SMOLAR_R, implemented in ECOM-3D by Gomez-Reyes and Blumberg (1995).

The scale of the model domain was chosen to be representative of the Androscoggin-Kennebec river system in the western Gulf of Maine. The base case was run on a grid with a horizontal resolution of 50 x 70 and 13 sigma layers (Figure 1). The grid cells at the coast measured 1.5 x 3 km and linearly increased in size to 3 x 3 km at the eastern open boundary. The bottom sloped linearly seaward from 15 m at the coast out to a depth of 190 m. A freshwater source of 1500 m3/s flows into a uniform coastal ocean of 32 psu with an ambient southward current of 5 cm/s specified at the northern boundary. The southern outflow boundary had an Orlanski radiation condition on elevation and a no-gradient condition on salinity and temperature.

fig. 1
Figure 1. The Base Case has 13 Sigma Layers, 50 by 70 grid cells, and a sloping bottom of 15 meters at the shore to a depth of 200 meters at the seaward bound. An ambient current of 5 cm/s flows southward.

All cases had a bottom roughness of 0.003 meters. Smagorinsky mixing is employed in the horizontal with a coefficient of 0.05. In the vertical, the Mellor-Yamada level 2.5 scheme was used, with the minimum vertical background mixing held at 5 x 10-6 m2/s. The CENTRAL scheme was always used for advecting momentum.

Results

Our study was strongly motivated by our initial finding that the base case run with the four different advection schemes led to markedly different results (Figure 2). This led to an analysis of the sensitivities of the run to horizontal and vertical resolution and wind forcing. In each test, the base case was modified and the four advection schemes tested.

fig. 2 movie of fig. 2
Figure 2. The original problem of four runs with identical forcing, but different advection schemes yielding very different results. (This figure is also available as a movie.)

Sensitivity to Horizontal Resolution

The initial approach to determining the "correct" solution was to increase the horizontal resolution until all schemes converged on the same answer. The hope was that we could identify which of the original resolution simulations best approximated the convergent solution.

Doubling the horizontal grid resolution, however, did not cause the four schemes to converge (Figure 3). Another doubling of the horizontal resolution also had surprisingly little effect.

fig. 3 movie of fig. 3
Figure 3. Same calculations as in Figure 1, but using a grid with doubled horizontal resolution (Figure 2). (This figure is also available as a movie.)

Sensitivity to Vertical Resolution

Doubling the vertical resolution had a much more dramatic effect on improving convergence of the solutions (Figure 4). It revealed quite clearly that the initial distribution of vertical levels was seriously underresolving the surface trapped plume. This became especially true as the surface sigma layers expanded in the vertical in the deeper offshore waters.

fig. 4 movie of fig. 4
Figure 4. Increased vertical resolution: With 28 layers, the plume structures converge towards the solution most similar to the Base Case SMOLAR_R solution. (This figure is also available as a movie.)

In order to further explore the effect of vertical resolution, an exponential function for describing the distribution of vertical layers was defined as

equation

where y are the sigma levels, gamma is a shape parameter, and n is the number of vertical layers. Runs of varying vertical resolution were then carried out using the CENTRAL and SMOLAR_R difference schemes. For the runs described herein, gamma was held at 3. This shape function is useful in that it provides for a bottom layer of equal thickness for all tests (Figure 7) and allows for a simple direct comparison of runs with any integer number of vertical layers with packing of levels near the surface.

The increase in vertical resolution results in a gradual convergence in the model solutions (Figures 5) for the CENTRAL scheme. The base case 13 sigma level spacing was not adequate for the CENTRAL scheme to resolve the plume, whereas the SMOLAR_R at this resolution produced results qualitatively similar to the converged solution. Increasing the vertical resolution results in a convergence of the solution, with the CENTRAL runs showing little change from the 28 layer to the 62 layer run.

fig. 5 movie of fig. 5
Figure 5. The eddies were observed to be reduced with increased vertical resolution. Note that this 13 layer run has different vertical segmentation than the original 13 layer run (see Figure 7). (This figure is also available as a movie.)

Figure 6 shows the problem of vertically underresolving a surface trapped river plume. In order to separate sigma-coordinate effects from vertical resolution effects, plume simulations were run at various vertical resolutions over a flat bottom and a similar solution convergence with higher vertical resolution was observed. This was confirmed over a flat bottom of various depths.

fig. 6
Figure 6. For the CENTRAL scheme, adequate vertical resolution of the plume is necessary to avoid eddies of numerical origin.

Figure 7 compares the vertical resolutions tested here with those used in past studies (Chao and Boicourt, 1986; Chao 1988; Oey and Mellor, 1993). While the different depths, slopes, plume structures and numerical models prohibit direct comparisons, the previous studies appear to be in the range where the choice of advection scheme and vertical resolution negatively affect the solutions.

fig. 7
Figure 7. Comparison of vertical segmentation. "Orig" refers to the base case. "CB 86" refers to Chao and Boicourt, 1986. "Chao 88" refers to Chao, 1988. "Oey 93" refers to Oey and Mellor, 1993.

False Dense/Fresh Water Source

Numerical dispersion resulting from advection calculation of scalar variables can result in ripples and over- and under-shoot. In the simulation of a surface trapped river plume, this can represent a significant source of false dense and fresh water where the ambient current encounters the plume and both higher than ambient salinity water and bogus fresh water form in the near-stationary upstream ripples. The dense plume sinks to the bottom and flows along the slope and offshore. This occurs in both the dispersive CENTRAL and SMOLAR_2 schemes, but not in the UPWIND or SMOLAR_R schemes (Figure 8.) Although this dense-water source did not play a dynamically significant role in our study of the surface trapped river plumes, studies interested in bottom boundary layer phenomena should be aware that this seemingly interesting dynamical feature is, in fact, a purely numerical artifact. Also, the additional fresh water can result in artificial freshening of the plume. This effect can easily go unnoticed by total freshwater budgets, a common method of validating a plume model, because the spurious freshwater is exactly offset by the spurious dense water.

fig. 8
Figure 8. The bottom sigma layer shows the false dense water flowing downstream and downslope. SMOLAR_R does not produce the false dense water.

It is interesting to note that neither increased vertical nor horizontal resolution eliminates this artifact, as dispersive schemes will generate 2 dx ripples upstream of the salinity front regardless of resolution. However, the average bottom layer salinity over the entire domain does come closer to the initial background salinity of 32 with increased vertical resolution, but not horizontal resolution (Figure 9). This comparison is valid because the bottom layer of each of these runs is equally thick (Figure 7).

fig. 9
Figure 9. The average bottom layer salinity, a proxy for the amount of false dense water produced, as a function of vertical and horizontal resolution at time = 37.25 days. The solid line shows the dependence on vertical resolution and the dashed line shows the dependence on horizontal resolution.

The observation in Figure 9 also shows an association of the false dense- and fresh-water production with the numerical eddy shedding. It is with the 13 Layer CENTRAL run that these eddies disappear (Figure 5) and the average bottom layer salinity approaches the background salinity. In addition, the SMOLAR_R and UPWIND schemes produce neither the scalar over- and under-shoots nor the eddies. This association is the subject of future work.

Sensitivity to wind forcing

While idealized no-wind simulations over month time-scales are interesting, in most regions wind varies significantly in the 2 to 10 day weather band. In order to determine whether the presence of wind forcing affects our conclusions, a sinusoidal north-south wind with a 4 day period and a 5 m/s maximum was imposed. The domain responds to the alongshore wind forcing with periodic upwelling and downwelling, driving the plume on and offshore consistent with Ekman dynamics.

This run demonstrates the dramatic effect of wind on the behavior of surface trapped river plumes (Figure 10). Strong mixing occurs during upwelling wind conditions (Fong, 1997). This mixing can overshadow some of the subtleties caused by changes in advection scheme and resolution.

fig. 10 movie of fig. 10
Figure 10. The effect of wind forcing on plume behavior: wind effects can mask some of the differences caused by advection scheme and resolution (compared to Figure 2). (This figure is also available as a movie.)

Implications to the Realistic Simulation

A model of the Western Gulf of Maine using realistic geometry and forcing was examined in light of the sensitivity analysis. Due to the dominant role of wind mixing on the plume, the differences due to resolution and advection scheme are much more subtle than in the simulations with no wind forcing.

We examined these subtle differences between advection schemes and resolution in the presence of realistic wind forcing. First, we find the elimination of the false dense water source by using the SMOLAR_R advection scheme highly desirable. The modeled biology of red tide includes the process of cell cysts, or seeds, in the sediments germinating and swimming upward in response to environmental factors, including salinity and temperature. The presence of a false dense water source has the definite potential to influence this source of cells.

Next, we wished to test whether or not our optimal advection scheme and resolution combination, SMOLAR_R with normal vertical resolution, provided a significant improvement in the agreement with the data (Figure 11). The plumes appear qualitatively similar, with the SMOLAR_R scheme and CENTRAL scheme with high vertical resolution in good agreement, as expected (two center panels). The CENTRAL scheme with the normal vertical resolution (far left), however, also appears qualitatively similar, again suggesting that advective effects are being overshadowed by wind-driven mixing and other processes. The fourth panel shows dramatically the effect of wind mixing by revealing the simulated salinity field without wind. It is temping to favor the CENTRAL scheme with normal vertical resolution, as this would be the lowest-cost solution. Because our biological growth functions are based on temperature and salinity, however, over- and under-shoots and a dense water source are highly undesirable. We therefore are using SMOLAR_R in our ongoing coupled physical-biological simulations in realistic Gulf of Maine domains.

fig. 11 movie of fig. 11
Figure 11. Surface salinity in the Gulf of Maine runs with identical forcing, except for the no-wind case on the right. The no-wind case uses the CENTRAL Scheme. (This figure is also available as a movie.)

Computational Expense

Both higher vertical resolution and higher order advection schemes add computational expense to a model simulation. In the case of the idealized domain, using SMOLAR_R at a lower vertical resolution proved much cheaper than using CENTRAL at a higher vertical resolution (Figure 12). The extent of the savings provided by a higher-order scheme depends on the number of vertical levels and the degree of time splitting between internal and external modes.

fig. 12
Figure 12. Computational Expense for the idealized river plume domain.

Conclusions

These studies demonstrate that for horizontal and vertical resolutions commonly used in river plume modeling, advection schemes can have a large effect, particularly under low wind conditions. Higher order advection schemes such as the recursive Smolarkiewicz scheme can provide results that can only be obtained with much higher grid resolution and at higher cost. The SMOLAR_R scheme also had the desirable characteristic of eliminating the formation of spurious fresh- and dense-water plume upstream of the river source.

While the advection scheme made a great difference in the no-wind cases, the addition of typical wind forcing greatly modified the plume characteristics, and reduced the sensitivity to the choice of advection scheme. In our realistic simulations with western Gulf of Maine topography and wind forcing, the choice of advection scheme made little qualitative difference in the salinity fields, thus suggesting that the inexpensive CENTRAL scheme might be the optimal choice if false dense water source is not a major concern. In our case, however, consideration of the biological ramifications of over- and under-shoots caused us to choose the more expensive SMOLAR_R scheme. This highlights the fact that numerical models should be designed to allow for a choice of advection schemes with different degrees of accuracy and cost, so that simulations where advection plays a stronger or weaker role can be optimized accordingly.

References

Blumberg, A.F. and Mellor, G.L., 1987, A description of a three-dimensional coastal model, in Three-Dimensional Coastal Ocean Models, N. Heaps [ed], Coastal and Estuarine Sciences, 4, p. 1-16.

Fong, D.A., Geyer, W.R., and Signell, R.P., 1997, The wind-forced response on a buoyant coastal current: Observations of the western Gulf of Maine plume, J. Marine Sys., 12, 69-81.

Franks, P.J.S. and Anderson, D.M., 1992, Alongshore transport of a toxic phyoplankton bloom in a buoyancy current: Alexandrium tamarense in the Gulf of Maine, Marine Biology 112, 153-164.

Gomez-Reyes, E. and Blumberg A.F., Pollutant Transport in Coastal Water Bodies, Computer Modelling of Seas and Coastal Regions II, Computational Mechanics Publications, Southhampton, U.K. 1995, pp. 87-94.

Margolin, L.G. and Smolarkiewicz, P.K, 1989, Antidiffusive velocities for multipass donor cell advection, Lawrence Livermore National Laboratory, Report W-7405-Eng-48.

Roach, P.J., 1976, Computational Fluid Dynamics, Hermosa Publishers, Albequerque, New Mexico.

Smolarkiewicz, P.K., 1983, A simple positive definite advection transport scheme with small implicit diffusion, Monthly Weather Review, 111, 479-486.

Smolarkiewicz, P.K., 1984, A fully multidimensional positive definite advection transport algorithm with small implicit diffusion, J. Comptutional Physics, 54, 325-362.

Smolarkiewicz, P.K., and T.L. Clark, 1986, The Multidimensional Positive Definite Advection Transport Algorithm: Further Development and Applications, J. Comptutional Physics, 67, 396-438.

Back to the SMIG Features Page

button bar

SMIG Home Mailing List Features Conferences Classes Reading Room Model Archives Feedback Home | Mailing List | Features | Conferences | Classes | Reading | Model Archives | Feedback

Stewart Rounds, SMIG coordinator <sarounds@usgs.gov>
U.S. Geological Survey
http://smig.usgs.gov/SMIG/features_0300/plumes.html
Last modified Thursday, 10-Mar-2005 17:35:49 EST
Privacy Statement · Disclaimer · FOIA · Accessibility