UBC Faculty Research and Publications

Solution convergence of flow over steep topography in a numerical model of canyon upwelling. Dawe, Jordan T.; Allen, Susan E. May 31, 2010

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata


52383-Allen_AGU_2009JC005597.pdf [ 907.13kB ]
JSON: 52383-1.0041906.json
JSON-LD: 52383-1.0041906-ld.json
RDF/XML (Pretty): 52383-1.0041906-rdf.xml
RDF/JSON: 52383-1.0041906-rdf.json
Turtle: 52383-1.0041906-turtle.txt
N-Triples: 52383-1.0041906-rdf-ntriples.txt
Original Record: 52383-1.0041906-source.json
Full Text

Full Text

Click Here for Full Article Solution convergence of flow over steep topography in a numerical model of canyon upwelling Jordan T. Dawe1 and Susan E. Allen1 Received 23 June 2009; revised 10 December 2009; accepted 28 December 2009; published 11 May 2010. [1] The convergence as resolution is increased is examined for a numerical model simulation of upwelling over a continental shelf canyon. A series of idealized one‐ and two‐dimensional models are used to plot the dependence of bottom Ekman layer velocity structure and transport on model resolution. Using these results as a guide, a bottom Ekman layer‐resolving three‐dimensional numerical model is constructed. This model is used to simulate a laboratory model of canyon upwelling. Alongshore velocity, volume of fluid upwelled through the canyon onto the shelf, and change in depth of the canyon water are examined at three grid resolutions. Numerical model particle trajectories are compared with a laboratory model of canyon upwelling to validate the model. The results suggest that current computational cluster power is sufficient to accurately simulate all aspects of upwelling in a steep canyon with the exception of flow separation from the canyon wall. Citation: Dawe, J. T., and S. E. Allen (2010), Solution convergence of flow over steep topography in a numerical model of canyon upwelling, J. Geophys. Res., 115, C05008, doi:10.1029/2009JC005597. 1. Introduction [2] The interaction of flows with steep topography is important to a wide range of ocean processes, including continental shelf fluid exchange [Allen, 2004], tidal energy conversion [Laurent et al., 2003], and deep ocean mixing [Kunze and Smith, 2004]. Canyons on the edge of the con- tinental shelf are a good test case for the study of these kinds of flow. First, canyons enhance shelf exchanges during upwelling or downwelling events, which has impacts on shelf biology and deep ocean water properties [Allen et al., 2001; Allen, 2004;Bosley et al., 2004; Ladd et al., 2005;Williams et al., 2006], making the problem of interest. Second, the importance of canyons in shelf exchange processes means that there have been several field studies done on a variety of coastal canyon systems which give a reasonable picture of the dynamics of the system [Hickey, 1997; Sobarzo et al., 2001; Wang et al., 2003; Palanques et al., 2005]. Third, coastal canyon systems are easily represented in a laboratory tank model [Pérenne et al., 2001a; Allen et al., 2003; Boyer et al., 2006], giving us an idealized theoretical test bed to compare with model simulations. [3] The elements needed to represent a canyon system are a continental shelf cut by a canyon, density stratification, and an alongshore geostrophic current. Coastal upwelling or downwelling events in the ocean are triggered when the wind and Coriolis forces conspire to push fluid toward or away from the land. This sets up a barotropic pressure gra- dient that drives a geostrophic current that runs parallel to the coastline. Friction then creates a bottom Ekman layer that transports fluid onshore or offshore, relaxing the barotropic pressure gradient [Thorpe, 1987]. However, fluid that is transported up or down the steep continental shelf slope advects the vertical density gradient and creates a baroclinic pressure gradient which eventually arrests the bottom Ekman flow [MacCready and Rhines, 1993]. [4] Canyons serve to enhance these events in two ways: first, they act like a ramp and advect fluid less vertically than the shelf, delaying the shutdown of the boundary layer; second, the walls of the canyon prevent alongshore flow within the canyon, breaking geostrophic balance and expos- ing the interior of the canyon to the barotropic pressure gra- dient [Freeland and Denman, 1982; Klinck, 1996]. Canyons can further be divided into “long” canyons, in which the breaking of geostrophy within the canyon is largely accom- plished by convergence and divergence of the topographic contours, and “short” canyons, in which it is the result of advection [Allen, 2000; Waterhouse et al., 2009]. [5] Laboratory models that include all these effects can be created using a cylindrical rotating tank containing a cylin- drically symmetrical topography representing the continental shelf and canyon and filled with a vertically stratified fluid [Pérenne et al., 1997]. Alongshore currents can be created by varying the tank rotation rate [Pérenne et al., 2001a]. Several studies have been performed using this type of setup and have been compared with numerical models [Boyer et al., 2000; Pérenne et al., 2001b; Allen et al., 2003; Boyer et al., 2004; Jaramillo, 2005; Haidvogel, 2005]. The results, in general, show good agreement between tank and numerical models for canyons with shallow slopes and poor agreement for canyons with steep slopes. 1Department of Earth and Ocean Sciences, University of British Columbia, Vancouver, British Columbia, Canada. Copyright 2010 by the American Geophysical Union. 0148‐0227/10/2009JC005597 JOURNAL OF GEOPHYSICAL RESEARCH, VOL. 115, C05008, doi:10.1029/2009JC005597, 2010 C05008 1 of 13 [6] These continuing problems in accurately representing the flow field in these models are the motivation for the following work. Here we concentrate on a z level model, the Massachusetts Institute of Technology general circulation model (MITgcm) [Marshall et al., 1997], using it to simu- late a laboratory tank experiment of an upwelling event over a continental shelf canyon. By running the model in a variety of configurations we attempt to estimate the grid spacing necessary to resolve the numerical solution and identify which dynamical features are hardest to resolve. Then we use data from the laboratory experiment to examine how well our numerical solution represents reality. Our results suggest that the accuracy of numerical models of flow over steep topog- raphy may be limited by the impact of small‐scale nonlinear processes such as flow separation. 2. Two Models [7] In this study we discuss two models: a laboratory model of an upwelling continental shelf canyon system and a numerical model simulation of the laboratory model. The laboratory model data we refer to are the result of previous work by Allen et al. [2003], Mirshak and Allen [2005], Bowie [2006], Waterhouse [2006], and Waterhouse et al. [2009]. 2.1. Laboratory Model [8] The laboratory model consists of a 1 m diameter cylindrical Plexiglas tank placed upon a rotating table. The tank contains a radially varying topography intended to simulate a continental shelf (Figure 1). A 22° section of this shelf is designed to be removable to accommodate a variety of canyon inserts. [9] A vertical density gradient is created in the tank by varying the salinity of the water used to fill the tank using the two‐bucket method of Oster [1965] [see also Hill, 2002]. This method uses two buckets, one filled with fresh water and the other filled with salt water, and the two are connected by a small tube with a valve that is initially closed. Fluid from the freshwater bucket is slowly drained into the laboratory tank while the valve between the buckets is opened, allowing salt water into the freshwater bucket where a mechanical stirrer ensures that the bucket is well mixed; thus, the fluid entering the tank increases in salinity linearly with time. The initial salinity of the saltwater tank was chosen so that the buoyancy frequency N at the model shelf break would be 2.2 s−1. [10] The tank is filled while rotating with Coriolis parameter f = 0.4 s−1 and once filled is left to spin up for 3 h to allow the fluid to equilibrate into solid body rotation. Upwelling events are simulated in the tank by varying the rotation rate of the table to create an alongshore flow. During an experiment the tank’s rotation rate is increased from f = 0.4 s−1 to f = 0.52 s−1 with a constant angular acceleration over a period of 27.3 s, a duration chosen to simulate a day‐long upwelling event. Data were collected immediately after the rotational acceleration ended. 2.2. Numerical Model [11] We use the MITgcm, a primitive equation z coordi- nate model, to construct a numerical model of the laboratory tank. Since the flow in the tank is laminar, we set the model viscosity and diffusion coefficients to the molecular values for water: n = 10−6 m2 s−1 for viscosity and  = 1.5 × 10−9 m2 s−1 for salinity diffusion. Temperature is set to be constant, and the model is initialized with a vertical salinity gradient similar to that which would be produced via the Oster [1965] two‐bucket method, modified by 3 h of vertical diffusion (Figure 2). A linear equation of state is used for density calculations. [12] The model is configured to run in nonhydrostatic mode in order to capture internal wave dynamics. However, comparison of the nonhydrostatic runs presented here with equivalent hydrostatic runs shows essentially identical results in all model fields, save for the vertical velocity. Maxima in the hydrostatic vertical velocity field tend to lie directly above or below minima, while the nonhydrostatic field has a more complex spatial pattern. However, the vertical velocity dif- ferences occur in regions where the vertical velocity is small and thus have little effect on the simulation. Nevertheless, we use the nonhydrostatic dynamics to ensure that our simulation is as realistic as possible. [13] The model is run in cylindrical coordinates to take advantage of the cylindrical symmetry of the tank. Mea- surements of the laboratory tank are used to create the numerical model topography, the accuracy of which is Figure 1. Model topography viewed (top) in cross section along the canyon axis and (bottom) from above. The tank is 1 m in diameter and filled to depths of 9.1 cm over the abyss, 2.1 cm at the shelf break, and 3 mm at the outside wall. The base of the shelf slope has a radius of 20.5 cm, and the shelf break has a radius of 27.5 cm. The central cyl- inder added to the numerical model to ensure numerical sta- bility is indicated by dashed lines and is 10 cm in radius. The dotted line indicates the height of the canyon walls. DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 2 of 13 enhanced through the use of MITgcm’s partial grid cell representation of topography [Adcroft et al., 1997]. The tank topography is modified in the numerical model by the addition of a 10 cm radius cylinder in the center of the tank, in order to avoid Courant‐Friedrichs‐Lewy criterion viola- tions where the cylindrical grid spacing is small. Friction with the wall of this cylinder is disabled to minimize its influence on the model flow. A variety of different grid spacings were used in this study to examine the convergence of the numerical simulation. [14] To simulate the change in rotation rate in the tank experiments that drives the tank flow, an azimuthal body force is applied to the tank for 27.3 s. The Coriolis param- eter is maintained at f = 0.52 s−1 throughout the model run; analytical calculations indicate that this constant Coriolis parameter results in negligible errors. 3. Resolving the Ekman Layer [15] If one ignores nonlinear effects, such as flow sepa- ration within the canyon, the smallest scale that must be resolved in the tank is the viscous bottom Ekman layer thickness. This thickness is set by the Ekman depth,  ¼ ffiffiffiffiffi 2 f s ; ð1Þ which equals 1.96 mm for our values of n and f. Thus, one might expect that a vertical resolution of the order of 2 mm would be necessary to resolve the velocities in the Ekman layer. [16] However, the velocity structure of the Ekman layer is of secondary importance compared to the amount of fluid transported. Ekman transport is the primary mechanism by which fluid moves from areas of high pressure to low pressure, reducing surface height gradients and spinning down the geostrophic flow. If the Ekman transport is incorrect, the spin‐down rate of the geostrophic flow will be also be incorrect, and this will affect the strength of the pressure gradients which drive upwelling on the shelf and in the canyon. Furthermore, the Ekman transport also affects the advection of density gradients up the continental shelf slope. This in turn will alter the rate at which the Ekman transport is shut down by the baroclinic pressure gradient created as dense water advects up the slope [MacCready and Rhines, 1993]. In sections 3.1 and 3.2 we use a set of simple one‐ and two‐dimensional models to estimate the model resolution needed to resolve bottom Ekman transports. 3.1. One‐Dimensional Model [17] Consider an infinite f plane with a steady geostrophic flow U in the x direction, with a no‐slip condition at the bottom. The equations governing the evolution of this sys- tem are @u @t ¼ fv 1 0 @P @x þ  @ 2u @z2 ð2Þ @v @t ¼ fu 1  @P @y þ  @ 2v @z2 ; ð3Þ where u and v are the velocities in the x and y directions, respectively; P is the pressure; and r is the density of the fluid. The initial conditions for this system are then u =U, v = 0, and P = P0 − rfUy, where P0 is a constant dependent on the depth of the fluid. Bottom friction will quickly result in the creation of an Ekman velocity spiral, reducing the bottom u velocity and driving a cross‐flow v velocity within the fric- tional boundary layer. [18] Discretization of these equations in the z direction and integration forward in time with a fourth‐order Runge‐ Kutta scheme allow us to examine the effect of vertical resolution on the representation of this Ekman layer. A 2 cm deep model was configured with the velocity at the top of the domain fixed to u = U and the viscosity, density, and f parameters set to appropriate values for the laboratory tank. The model was then run with a variety of vertical resolu- tions, from Dz = 2 cm to Dz = 0.02 mm, for 300 s of model time (about 12 inertial periods). This was sufficient to let the solution converge to a steady state where velocities were changing <0.01% s−1. [19] Comparison of the theoretical solution of the bottom Ekman boundary layer [Kundu, 1990] with a low‐resolution (Dz = 2 mm) run shows that while the velocity structure of the run is reasonable, the cross‐geostrophic flow Ekman transport is overestimated (Figure 3a). Calculating the mod- el’s Ekman transport as a function of resolution (Figure 3b) indicates that in order for the numerical transport to be within 1% of the theoretical value a vertical resolution of ∼0.2d is required; for the transport to be within 2%, ∼0.3d is required. 3.2. Two‐Dimensional Model [20] The experimental tank topography is far more com- plicated than an infinite flat plane, with several changes in slope from the central abyss to the shallow continental shelf. The continental shelf’s 4.5° slope allows for more compli- cated dynamics such as upslope advection of the vertical density gradient and flow interactions with the edges of the Figure 2. Model initial density profile. DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 3 of 13 shelf and the tank. To simulate these effects, we turn to the MITgcm. [21] For now, we ignore the effect of the canyon and simply look at the Ekman layer in a tank with an unbroken shelf. This allows us to run theMITgcm in a two‐dimensional mode under the assumption that the flow is cylindrically symmetric. This simplification lets us perform many model runs and compare the Ekman layer representation as the model resolution changes.Wemake a set of runs with a radial resolution varying from 1.25 to 6.25 mm and a vertical resolution varying from 0.4 to 2.0 mm, with an additional ultrahigh‐resolution run done at 0.625 mm horizontal and 0.2 mm vertical resolution in order to judge the convergence of the lower‐resolution models. [22] First, we integrate the radial velocity at the end of the ultrahigh‐resolution run (t = 30 s) from the bottom of the tank to the surface to generate the r‐z stream function. The r‐z stream function shows two distinct overturning cells, one in the abyss and one on the shelf, separated by the steep continental slope (Figure 4a). No bottom layer flow is apparent between these cells as the steepness of the 45° slope quickly results in baroclinic pressure gradients that shut down the Ekman transport. [23] Comparison of the 6.25 mm horizontal × 2.0 mm vertical resolution model stream function with the ultrahigh‐ resolution stream function (Figure 4b) shows several differ- ences. In agreement with the one‐dimensional model results, the 2 mm vertical resolution model overestimates the Ekman transport over the flat abyssal plane relative to the higher‐ resolution model. However, over the sloped continental shelf the situation is reversed, with the low‐resolution model having less Ekman transport than the high‐resolution model. One obvious reason for this is the staircase‐like representa- tion of the smooth continental shelf in the lower‐resolution model, which obstructs flow within the Ekman layer. The effects of these sudden depth transitions are mitigated in MITgcm via the inclusion of a partial grid cell representation of bottom topography, but this is not enough to completely alleviate the problem. [24] A clearer view of how the Ekman transport is affected by resolution and topographic slope can be found by exam- ining how the net Ekman transport varies over the radius of the tank (Figure 5). To calculate the cross‐geostrophic Ekman transport, we assume that all flow directed radially outward within 1 cm of the bottom is in the Ekman layer and sum all positive radial transports over the entire model run. We do this for five of the two‐dimensional model runs: one with low vertical and low horizontal resolution (6.25 mm horizontal × 2 mm vertical), one with high vertical and high horizontal Figure 4. Two‐dimensional radial‐depth stream functions for models run at (a) 0.5 mm horizontal × 0.2 mm vertical and (b) 6 mm × 3 mm grid spacings. Figure 3. One‐dimensional Ekman layer model. (a) Steady state solution of one‐dimensional Ekman layer model cross‐ mean flow velocity in mm s−1 with Dz = 2 mm (bars) com- pared to the analytical solution (line). Transport is represented by the area of the bars or the area inside the curve. (b) Model cross‐mean flow transport at a variety of vertical resolutions, normalized by the theoretical transport. Vertical resolution is given in units of the Ekman depth. DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 4 of 13 resolution (1.25mmhorizontal × 0.4mm vertical), twomixed resolution runs (6.25 mm × 0.4 mm and 1.25 mm × 2 mm), and the ultrahigh‐resolution run (0.625 mm × 0.2 mm). [25] The resulting transports show two maxima, one in the abyss and one over the shelf, and a minima over the steep continental slope. In the abyss the transport shows little sensitivity to horizontal resolution but shows the same overestimation of transport at low vertical resolution that the one‐dimensional model did. At 0.4 mm vertical resolution (0.2d) the transport is almost identical to the ultrahigh‐ resolution run, and the model representation appears to have converged. Over the slope all but the very lowest resolution run appear to have converged, indicating that steep topo- graphic slopes need little resolution to be well resolved, despite the presence of a sizable transport. Even in the low‐ resolution run the transports only differ from the ultrahigh‐ resolution run near the shelf break and appear to be more related to the shelf break Ekman transports which pull water up from the abyss. [26] The shelf transports show the strongest sensitivity to resolution, with the low‐resolution runs underestimating the transport by ∼15%. Increasing either the horizontal or the vertical resolution improves the transports roughly equally. However, the 2 mm vertical resolution models both display a small radial oscillation in the transport strength due to tiny vertical recirculation cells that occur where the flow meets a step in the ill‐resolved topography. This suggests that the vertical resolution is still slightly more important to achieving a realistic representation of the transport. The high vertical and high horizontal resolution run shows good overall agreement with the ultrahigh‐resolution transport, although the transport is still slightly underestimated. [27] Clearly, the shelf Ekman layer is the most difficult part of the flow to resolve, and so we focus on this to help pick a suitable resolution for a full three‐dimensional model. We run the two‐dimensional model at a variety of horizontal and vertical resolutions and measure the total Ekman trans- port over the shelf at r = 0.3625 m, a point chosen to be near the middle of the shelf and to be at the minimum of the low‐ resolution transport oscillations, in order to maximize the transport errors. Comparing these values as a percentage of the ultrahigh‐resolution (0.625 mm × 0.2 mm) transport allows us to plot out a two‐dimensional space representing how shelf transport varies with resolution (Figure 6). From this it is clear that both horizontal and vertical resolution are equally important in resolving the shelf Ekman layer flow. 4. Representation of Upwelling in a Continental Shelf Canyon [28] Having examined the grid spacing needed to resolve two‐dimensional radial Ekman transport over a model con- tinental shelf, we move to a full three‐dimensional simulation Figure 6. Net midshelf bottom Ekman layer transport as a function of vertical and horizontal grid resolution. Trans- ports are normalized by the model transport at a resolution of 0.625 mm horizontal × 0.2 mm vertical. Figure 5. Bottom Ekman layer transport (m2 s−1, given in units of flux per centimeter of azimuthal dis- tance) as a function of radius for a variety of model resolutions. DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 5 of 13 of a continental shelf canyon upwelling system. To minimize our computational requirements, we employ a stretched grid in the azimuthal direction in order to concentrate model resolution over the canyon (Figure 7). Grid resolution is unevenly distributed, with more grid cells downstream of the canyon to better represent the canyon’s wake. The spacing of the stretched grid was chosen so that grid cells would be approximately square over the canyon. Three main runs are performed with grid sizes of 64 azimuthal × 81 radial × 32 depth (5mmhorizontal grid spacing × 2.8mmvertical grid spacing), 128 × 134 × 91 (3 mm × 1 mm), and 256 × 268 × 182 (1.5 mm × 0.5 mm). Referring to Figure 6, these resolutions resulted in Ekman transports of ∼85%, ∼92.5%, and ∼97% of the ultrahigh‐resolution run. We choose a variety of physical quantities that are meaningful in the simulation of canyon upwelling and examine their represen- tation in each model run for evidence of convergence of the solution. [29] First, we look at the alongshore velocity on the opposite side of the domain from the canyon. Far away from the area of flow disturbed by the canyon, the alongshore velocity structure is fairly similar at all three model resolu- tions with a radially increasing velocity in the abyss, a velocity maximum in a small jet just offshore of the shelf break, and a broad region of roughly constant velocity over the shelf (Figure 8). Viscous boundary layers are apparent close to the topography. As the model resolution increases, the maximum speed within the jet slightly increases while the velocity over the shelf decreases. These changes are consis- tent with the two‐dimensional model results that showEkman transport increasing with resolution. More upslope Ekman transport over the shelf results in greater reduction of surface pressure gradients, which in turn reduces the alongshore barotropic geostrophic velocity. At the same time, the downslope return flow above the Ekman layer results in a depression of isopycnals at the edge of the shelf break, driving a baroclinic pressure gradient that creates the off- shore jet; increased Ekman transport thus results in a stronger offshore jet. [30] Between the 5 × 2.8 mm model (Figure 8a) and the 3 × 1 mm run (Figure 8b), the strength of the offshore jet increases by 2%, while the velocity over the shelf decreases by 3%–10%. Between the 3 × 1 mm (Figure 8b) and 1.5 × 0.5 mm (Figure 8c) runs, the changes are much smaller, with a <1% change in the jet velocity and a reduction of 1%–4% over the shelf. Furthermore, the velocity structure over the shelf in the 5 × 2.8 mm model appears to be significantly influenced by the low‐resolution step changes in the shelf topography. This effect is much smaller at 3 × 1 mm reso- Figure 7. Stretched grid used in the 3‐D MITgcm model runs showing the concentration of resolution over the can- yon region. The grid shown here only has a resolution of 12.5 mm, for clarity. The grid cells were evenly subdivided to produce the higher resolutions used in the model runs. Figure 8. Azimuthal velocity in cm s−1 on the opposite side of the tank from the canyon at (a) 5 mm horizontal × 2.8 mm vertical, (b) 3 mm × 1 mm, and (c) 1.5 mm × 0.5 mm grid resolutions. DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 6 of 13 lution and is not apparent at all in the 1.5 × 0.5 mm run. Overall, it appears that the velocity structure over the abyss and continental slope regions is well resolved in all the models and is likely within 5% of being resolved over the shelf, consistent with the results from the two‐dimensional model. This is unsurprising as far away from the canyon the low‐resolution laminar flow should behave much like the two‐dimensional model results. [31] Next, we turn to measures of the upwelling within the canyon and examine the depth change of water upwelled within the canyon over the 30 s of the model run. Assuming that diffusion of density is negligible implies that density changes over the model run are solely due to advection of the initial model density structure. Since the initial density is stably stratified and varies only with depth, every density is associated with a single initial depth. By subtracting the depth of each point by the initial depth associated with the density at that point, we get the depth change of the water at that model point. [32] All models show a general pattern of strong upwell- ing within the canyon with a maximum at the canyon head, smaller upwelling in the Ekman layer over the shelf, and a slight downwelling elsewhere (Figure 9). In the 5 × 2.8 mm model (Figure 9a), the maximum depth change at the can- yon head is ∼9 mm, while in the higher‐resolution runs (Figures 9b and 9c) the maximum upwelling is ∼11 mm. In addition, the higher‐resolution runs show a larger region of downwelling over the shelf than the 5 × 2.8 mm model. In contrast to these differences, the 3 × 1 mm and 1.5 × 0.5 mm resolution model upwelling fields appear nearly identical, suggesting that the canyon upwelling solution is fully converged. [33] For a second measure of the canyon upwelling, we look at the total volume of water upwelled through the top of the canyon over the model run. We calculate this by summing the vertical and radial volume fluxes through the top of the canyon, defined as the surface of the topography if there was no canyon. The fluxes are summed only for radii >27.9 cm (the shelf break radius) so as to only include upwelling out the top of the canyon. Fluxes are sampled once per second and are assumed to be constant over that second. We perform this calculation on the three main model runs and additionally include the results from a run performed at 2 × 0.67 mm resolution (Table 1). The total fluid volume upwelled increases by ∼19% between the 5 × 2.8 mm and 3 × 1 mm runs but in contrast is nearly identical between the three highest‐resolution runs, with a variation of <1% in the volume of fluid upwelled. Again, this suggests that the 1.5 × 0.5 mm model’s canyon upwelling solution is fully converged. [34] The structure of the vertical velocity across the can- yon shows the structure of standing internal waves gener- ated by the flow over the canyon and serves as a test of the convergence of the MITgcm’s nonhydrostatic dynamics (Figure 10). A series of areas of alternating upward and downward motion are apparent as the flow moves over the canyon. The 5 × 2.8 mm run tends to overestimate the mag- nitude of these velocities (Figure 10a) relative to the higher‐ resolution runs. The magnitudes and overall patterns of the vertical velocity in the 3 × 1 mm and 1.5 × 0.5 mm runs are very similar, but many small‐scale differences are still apparent, especially where the flow enters the canyon on the upstream side. Nevertheless, the overall changes between the 3 × 1mm and 1.5 × 0.5mm runs are small and suggest that the 3 × 1 mm model solution has nearly converged. [35] Next, we examine the density field at the end of the model runs. The model density is interesting for a few rea- Figure 9. Vertical displacement in millimeters of fluid between the start and end of each model run along the canyon axis at (a) 5 mm horizontal × 2.8 mm vertical, (b) 3 mm × 1 mm, and (c) 1.5 mm × 0.5 mm grid resolutions. Table 1. Net Fluid Volume Upwelled Onto the Shelf From Within the Canyon During Each MITgcm Model Run Resolution (mm) Upwelling (cm3) 5 × 2.8 16.63 3 × 1 20.43 2 × 0.67 20.57 1.5 × 0.5 20.48 DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 7 of 13 sons. First, it is the result of the integrated three‐dimensional density advection over the model run and thus is sensitive to the details of the time evolution of the system. Second, it lets us compare the density structure from this z coordinate model with the density from a previous s coordinate model exam- ined by Allen et al. [2003]. Third, Hickey [1997] measured temperature structure over Astoria Canyon, and these mea- surements allow us to compare the model density structure to ocean data. [36] During upwelling events in Astoria Canyon, Hickey [1997] found that the density field within the canyon formed a wedge‐shaped structure (Figure 11a), with fluid layers near the top of the canyon thinning on the upstream side and growing thicker on the downstream side as the flow pushed fluid downstream. A similar effect occurs in the model simulations (Figures 11b–11e), where we consider the “wedge” to be the fluid bounded by the 21 s and 24 s surfaces. [37] Models of s coordinates are generally considered more suitable for simulating flows over topographic features than z coordinate models. Allen et al. [2003] performed a study similar to the present one using the s coordinate Rutgers University model (SCRUM) [Song and Haidvogel, 1994] and found systematic problems with SCRUM’s rep- resentation of the laboratory model. To compare our current results with this previous research, we performed a simu- lation with SCRUM using our current model’s topography and forcing. All other aspects of the model were configured as given by Allen et al. [2003]. The SCRUM model run was configured with 16 s levels, which results in a vertical res- olution of ∼4 mm within the canyon, making the SCRUM resolution slightly coarser than the resolution of the 5 × 2.8 mm MITgcm simulation. [38] In both the SCRUM simulation (Figure 11b) and the 5 × 2.8 mm MITgcm simulation (Figure 11c) the density wedge is not well defined, with neither model doing a particularly good job representing the density wedge. This suggests that terrain coordinate models are not necessarily better than z coordinate models in representing this type of flow over steep topography. As we increase the model resolution (Figures 11d and 11e), the thickness difference between the upstream and downstream sides of the wedge consistently increases, and a well‐defined wedge structure emerges. The change in the downstream wedge thickness between the 3 × 1 mm and 1.5 × 0.5 mm runs is small, but the upstream thickness changes significantly between the runs. This thickness, however, is likely limited by the ver- tical resolution of the two models. [39] Finally, we compare the representation of potential vorticity (PV) in the model as resolution changes. PV con- servation represents a strong constraint on geostrophic flows and furthermore depends upon gradients in both velocity and density. For these reasons, PV should be more sensitive to changes in the flow structure than the velocity or density fields themselves. Because the model topography is much steeper than topography in the ocean, the horizontal vor- ticity in the model is of comparable magnitude to the ver- tical vorticity, and thus, we must use the full Ertel potential vorticity in our comparison: PV ¼ 2Wþr uð Þ   r: ð4Þ [40] To examine the representation of PV in the model, we calculate the difference between the initial and final PV in the model and then subtract from this the PV change due to the models’ azimuthal body force. This results in the PV changes due to advection and molecular diffusion (Figure 12). A three‐layer structure is apparent in the background PV changes, with increased PV in the canyon, decreased PV just over the shelf, and increased PV near the surface. The initial PV structure of the model is horizontally uniform and nega- tive, consisting of near‐zero PV near the surface and bottom Figure 10. Vertical velocity in mm s−1 across the canyon mouth at a radius of 30 cm at the end of the model run at (a) 5 mm horizontal × 2.8 mm vertical, (b) 3 mm × 1 mm, and (c) 1.5 mm × 0.5 mm grid resolutions. Figures 10a–10c are plotted looking toward the head of the canyon at 30 cm radius; mean alongshore flow is from left to right. DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 8 of 13 (due to the diffusion of the surface and bottom density gradients during the experiment setup) with a PV minimum near the shelf break. The PV changes in the canyon and near the surface can be explained as vertical advection of this background gradient as higher PV is advected upward in the canyon and downward near the surface. The middle band of decreased PV is the result of friction with the shelf topography. [41] Significant small‐scale changes exist between the PV representations at all model resolutions. These changes mostly represent frictional effects and so can be thought of as a map of the influence of the boundary layer on the flow. At 5 × 2.8 mm resolution the PV change has little small‐ scale structure, though there are some isolated hot spots of PV increase due to viscous effects near the walls of the canyon. At 3 × 1 mm resolution the PV change shows more viscous forcing, with regions of increased PV developing on the upstream and downstream flanks of the upper canyon. Additionally, a dipole structure of PV decrease lying over PV increase develops at the upstream edge of the canyon. Figure 11. Density structure (a) across the canyon mouth of Astoria Canyon during an upwelling event (based on Figure 3 from Hickey [1997]) and simulated by the (b) SCRUM model run and the MITgcm model runs at (c) 5 mm horizontal × 2.8 mm vertical, (d) 3 mm × 1 mm, and (e) 1.5 mm × 0.5 mm grid resolutions. The 26.6 st and 26.8 st surfaces in Figure 11a and the 21 s and 24 s surfaces in Figures 11b– 11e are highlighted to emphasize the “density wedge” structure that forms in an upwelling canyon. Figures 11a–11e are plotted looking toward the head of the canyon at 30 cm radius; mean alongshore flow is from left to right. DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 9 of 13 This dipole is the PV signature of the density wedge structure and results from viscosity acting on vorticity forced by den- sity surface squashing and stretching in the density wedge. At 1.5 × 0.5 mm resolution the overall pattern remains the same, but regions with fine PV structure, such as near the canyon walls, show more definition and stronger gradients. The PV dipole at the upstream edge of the canyon is larger in mag- nitude and spatial extent, nearly reaching the downstream canyon wall. Clearly, the model’s simulation of PV has not converged at 3 × 1 mm resolution; it may have converged by 1.5 × 0.5 mm resolution, but it is impossible to tell this from the model runs we have performed. Thus, small‐scale details of the simulation of advection and viscosity still may not be adequately represented. 5. Comparison With Lab Experiment [42] Over sections 3 and 4 of this paper we have examined the convergence of the model solution; however, we have not examined the accuracy of the solution that the model is converging to. To do this, we compare the horizontal flow field of the model with the flow field generated by a labo- ratory tank experiment performed by Bowie [2006]. [43] Bowie utilized neutral density wax particles to track the flow field of an upwelling flow in the laboratory tank. The particles were illuminated by horizontal sheets of col- ored light that allowed estimation of the pellets’ vertical locations: blue illuminated particles were between 14.5 and 17 mm depth, green particles were 17–18.6 mm, yellow particles were 18.6–19.5 mm, white particles were 19.5– 22.8 mm, and red particles were 22.8–25.5 mm. A video camera mounted above the tank was then used to record the trajectories of particles moving in and around the canyon. The density gradient was slightly weaker in Bowie’s labo- ratory experiments compared to our previous model runs (N = 2.0 s−1 at the shelf break), but all other physical parameters were the same. We performed a model run at 1.5 × 0.5 mm resolution (the highest resolution in section 4) that duplicated Bowie’s setup and recorded the model velocity field at 0.5 s increments. Numerical particles were then placed in the saved model velocity fields, and their trajectories were tracked from 28 to 30 s after the start of the model run, using a fourth‐order Runge‐Kutta algorithm. [44] Comparing the particle tracks from the laboratory experiment with the model experiment (Figure 13) shows that the model does a reasonable job of reproducing the laboratory tank’s flow field. In the abyssal area away from the shelf, particle tracks are smooth circular arcs indicating solid body rotation. Flow entering the canyon from the shelf turns left and flows up the canyon. Near the head of the canyon, the flow is directly up the canyon axis onto the shelf, with stronger flow on the downstream side. Some flow from the canyon onto the downstream shelf is also apparent. [45] However, there are two features of the model particle tracks that do not agree with the laboratory experiment. First, the model tracks show many more midlevel particles than are apparent in the laboratory tracks. This is due to oversaturation of the camera CCD during the laboratory experiment, resulting in many particles being misidentified as white [Bowie, 2006]. Second, the laboratory experiment shows particles on the upstream flank of the canyon mouth moving offshore into the abyss across isobaths before turning back onshore in a wide arc and flowing up the canyon. This offshore flow is a robust feature of the laboratory experiments and also appears in field data from Astoria Canyon [Hickey, 1997]. In contrast, the model particles tend to flow directly off the upstream flank of the shelf and then make a hard turn when they reach the downstream canyon wall. Thus, while the model solution is identical to the laboratory experiment in Figure 12. Ertel potential vorticity change between the start and end of each model run at (a) 5 mm horizontal × 2.8 mm vertical, (b) 3 mm × 1 mm, and (c) 1.5 mm × 0.5 mm grid resolutions. The effect of the model forcing on the poten- tial vorticity has been removed. Figures 12a–12c are plotted looking toward the head of the canyon at 30 cm radius; mean alongshore flow is from left to right. DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 10 of 13 most respects, it fails to simulate the flow on the upstream canyon flank realistically. 6. Discussion [46] We believe the failure of the model to simulate the flow on the upstream flank of the canyon is due to inadequate representation of flow separation. Temperature data from Astoria Canyon [Hickey, 1997] suggest that water drops slightly into the canyon as it leaves the shelf (Figure 11a, shown by the depth of the 26.6st isopycnal). This drop results in stretching of the water column as it leaves the shelf, which would tend to induce a counterclockwise, positive relative vorticity. Positive relative vorticity would create an offshore flow near the canyon wall and would increase the onshore flow within the canyon. In contrast, the model flow hardly drops at all as it crosses the canyon (Figure 11e, the 21 s isopycnal). This suggests that the model flow is separating from the upstream flank of the canyon too early. [47] We can estimate the amount of stretching needed to produce the observed flow by estimating the vorticity of the flow from Figure 13a. Particles move ∼5 mm offshore over the 2 s displayed in Figure 13 as they leave the upstream flank of the canyon and move a similar distance back onshore as they near the canyon axis; thus, we take the radial velocity change in this flow to be 0.25 − (−0.25) = 0.5 cm s−1. The distance from the upstream canyon edge to the canyon axis is ∼3 cm. Thus, the vorticity z of the flow is  ¼ @v @x ¼ 0:17 s1: ð5Þ [48] Using the potential vorticity of a homogeneous layer of thickness h, (z + f)/h, we can estimate the amount of vortex stretching needed to generate this much vorticity. Before leaving the upstream flank of the canyon the flow has a vorticity of 0, so f hi ¼  þ f hf ; ð6Þ Figure 13. Particle tracks from 28 to 30 s from (a, b) laboratory experiments and (c, d) MITgcm model output. Figures 13a and 13c show particles below the shelf break depth, and Figures 13b and 13d show particles above the shelf break depth (1.7 cm). Contours represent model topography, with the thicker contour indicating model shelf break depth. Particles are shown every 0.3 s and become slightly larger with time. DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 11 of 13 where hi is the initial thickness of the layer and hf is the thickness of the layer after stretching. Then we write hf hi ¼  f þ 1 ¼ 1:32: ð7Þ Examining Figure 11a, a 30% stretching is consistent with the stretching observed in Astoria Canyon. [49] One significant dynamical difference between our laboratory canyon and an ocean canyon is turbulence. True coastal canyon systems are far from laminar, with large turbulent flows that result in eddy diffusivities that vary in space and time and are of comparable magnitude for both viscosity and material properties of the fluid. Laboratory vertical and horizontal viscosities, once scaled to the size of the real ocean, are the same order of magnitude as measured ocean values, but salinity diffusivity will be 3 orders of magnitude higher. Bottom boundary layer depths are ∼20 m in the coastal ocean and thus are the same relative size as the Ekman layer in the laboratory tank, but the bottom turbulence will rise and fall in intensity with the currents and the tides, resulting in a time‐varying bottom layer depth. Despite this deficiency, our lab results compare favorably with field ob- servations, suggesting that these turbulent effects are small. 7. Conclusions [50] This paper examined the convergence of a numerical model simulation of an idealized laboratory model of an upwelling continental shelf canyon system with increasing resolution. A series of increasingly complex models were used to examine the conditions necessary to resolve flow features. A simple 1‐D model of bottom Ekman layer trans- port indicates that achieving a transport within 1% of the true transport value requires a vertical resolution of ∼0.2d, where d is the Ekman layer thickness. [51] A cylindrically symmetric, 2‐D, z level model shows that both horizontal and vertical resolution are important to resolution of the Ekman layer over sloped topography, while only vertical resolution is important over flat topography. This indicates that two concerns are important: achieving sufficient vertical resolution to resolve the velocity structure of the Ekman layer and sufficient horizontal resolution to achieve a smooth representation of the shelf slope, without abrupt steps in the topography. The partial cell parametri- zation present in the MITgcm helps the representation of topography but does not allow for improvements in the Ekman layer velocity structure and thus does not resolve problems with simulating Ekman layer transports. [52] Three full 3‐D MITgcm runs were performed with horizontal and vertical resolutions of 5 × 2.8 mm, 3 × 1 mm, and 1.5 × 0.5 mm. Significant changes in the model velocity, total canyon upwelling, and density advection occur between the 5 × 2.8 mm and 3 × 1 mm resolution runs, but there is little change in these fields between the 3 × 1 mm and 1.5 × 0.5 mm runs. Representation of potential vorticity, on the other hand, shows significant differences between the 3 × 1 mm and 1.5 × 0.5 mm runs and does not appear to have converged. Comparison of the 1.5 × 0.5 mm resolution model with an identically configured laboratory experiment shows that we are able to simulate most features of the laboratory tank flow field. However, the model un- derestimates the amount of offshore flow on the upstream flank of the canyon. This problem is likely due to the flow separating from the canyon wall too early. For problems that are not dependent on the details of the flow field on the upstream flank of the canyon, such as advection over the shelf or total upwelling through the canyon, modern com- putational clusters appear to have sufficient power to ade- quately represent an upwelling continental shelf canyon. [53] Acknowledgments. Funding for this research was provided to the second author via a NSERC Discovery grant. We would like to thank Mike Dinniman and John Klinck from Old Dominion University for pro- viding SCRUM model output and the two anonymous peer reviewers whose comments significantly improved aspects of this manuscript. This research has been enabled by the use of WestGrid computing resources, which are funded in part by the Canada Foundation for Innovation, Alberta Innovation and Science, BC Advanced Education, and the participating research institutions. WestGrid equipment is provided by IBM, Hewlett Packard, and SGI. Figures 1–13 were generated using the matplotlib library in the Python programming language. References Adcroft, A., C. Hill, and J. Marshall (1997), Representation of topography by shaved cells in a height coordinate ocean model, Mon. Weather Rev., 125, 2293–2315. Allen, S. E. (2000), On subinertial flow in submarine canyons: Effect of geometry, J. Geophys. Res., 105, 1285–1297. Allen, S. E. (2004), Restrictions on deep flow across the shelf‐break and the role of submarine canyons in facilitating such flow, Surv. Geophys., 25, 221–247. Allen, S. E., C. Vindeirinho, R. E. Thomson, M. G. G. Foreman, and D. L. Mackas (2001), Physical and biological processes over a submarine can- yon during an upwelling event, Can. J. Fish. Aquat. Sci., 58, 671–684, doi:10.1139/cjfas-58-4-671. Allen, S. E., M. S. Dinniman, J. M. Klinck, D. D. Gorby, A. J. Hewett, and B. M. Hickey (2003), On vertical advection truncation errors in terrain‐ following numerical models: Comparison to a laboratory model for upwelling over submarine canyons, J. Geophys. Res., 108(C1), 3003, doi:10.1029/2001JC000978. Bosley, K. L., J.W. Lavelle, R. D. Brodeur,W.W.Wakefield, R. L. Emmett, E. T. Baker, and K. M. Rehmke (2004), Biological and physical processes in and around Astoria submarine canyon, Oregon, USA, J. Mar. Syst., 50, 21–37. Bowie, A. W. (2006), Upwelling in short submarine canyons: A physical study, M.S. thesis, Univ. of B. C., Vancouver, B. C., Canada. Boyer, D. L., X. Zhang, and N. Pérenne (2000), Laboratory observations of rotating, stratified flow in the vicinity of a submarine canyon,Dyn. Atmos. Oceans, 31, 47–72. Boyer, D. L., D. B. Haidvogel, and N. Pérenne (2004), Laboratory‐numerical model comparisons of canyon flows: A parameter study, J. Phys. Oceanogr., 34, 1588–1609. Boyer, D. L., J. Sommeria, A. S. Mitrovic, V. K. C. Pakala, S. A. Smirnov, and D. Etling (2006), The effects of boundary turbulence on canyon flows forced by periodic along‐shelf currents, J. Phys. Oceanogr., 26, 813–826. Freeland, H., andK.Denman (1982), A topographically controlled upwelling center off Vancouver Island, J. Mar. Res., 40, 1069–1093. Haidvogel, D. B. (2005), Cross‐shelf exchange driven by oscillatory baro- tropic currents at an idealized coastal canyon, J. Phys. Oceanogr., 35, 1055–1067. Hickey, B. M. (1997), The response of a steep‐sided, narrow canyon to time‐variable wind forcing, J. Phys. Oceanogr., 27, 697–726. Hill, D. F. (2002), General density gradients in general domains: The “two‐ tank” method revisited, Exp. Fluids, 32, 434–440. Jaramillo, S. (2005), Numerical simulation of flow in a laboratory tank using a z‐coordinate numerical model, M.S. thesis, Univ. of B. C., Vancouver, B. C., Canada. Klinck, J. M. (1996), Circulation near submarine canyons: A modeling study, J. Geophys. Res., 101, 1211–1223. Kundu, P. K. (1990), Fluid Mechanics, 638 pp., Academic, San Diego, Calif. Kunze, E., and S. G. L. Smith (2004), The role of small‐scale topography in turbulent mixing of the global ocean, Oceanography, 17(1), 55–64. DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 12 of 13 Ladd, C., P. Stabeno, and E. D. Cokelet (2005), A note on cross‐shelf exchange in the northern Gulf of Alaska, Deep Sea Res., Part II, 52, 667–679. Laurent, L. S., S. Stringer, C. Garrett, and D. Perrault‐Joncas (2003), The generation of internal tides at abrupt topography, Deep Sea Res., Part I, 50, 987–1003. MacCready, P., and P. B. Rhines (1993), Slippery bottom boundary layers on a slope, J. Phys. Oceanogr., 23, 5–22. Marshall, J., A. Adcroft, C. Hill, L. Perelman, and C. Heisey (1997), A finite‐volume, incompressible Navier Stokes model for studies of the ocean on parallel computers, J. Geophys. Res., 102, 5753–5766. Mirshak, R., and S. E. Allen (2005), Spin‐up and the effect of a submarine canyon: Applications to upwelling in Astoria Canyon, J. Geophys. Res., 110, C02013, doi:10.1029/2004JC002578. Oster, G. (1965), Density gradients, Sci. Am., 213(2), 70–76. Palanques, A., et al. (2005), General patterns of circulation, sediment fluxes and ecology of the Palamós (La Fonera) submarine canyon, northwestern Mediterranean, Prog. Oceanogr., 66, 89–119. Pérenne, N., J. Verron, D. Renouard, D. L. Boyer, and X. Zhang (1997), Rectified barotropic flow over a submarine canyon, J. Phys. Oceanogr., 27, 1868–1893. Pérenne, N., D. B. Haidvogel, and D. L. Boyer (2001a), Laboratory‐numerical model comparisons of flow over a coastal canyon, J. Atmos. Oceanic Technol., 18, 235–255. Pérenne, N., J. W. Lavelle, D. C. Smith IV, and D. L. Boyer (2001b), Impulsively started flow in a submarine canyon: Comparison of results from laboratory and numerical models, J. Atmos. Oceanic Technol., 18, 1698–1718. Sobarzo, M., M. Figueroa, and L. Djurfeldt (2001), Upwelling of subsur- face water into the rim of the Biobío submarine canyon as a response to surface winds, Cont. Shelf Res., 21, 279–299. Song, Y., and D. B. Haidvogel (1994), A semi‐implicit ocean circulation model using a generalized topography‐following coordinate system, J. Comput. Phys., 115, 228–244. Thorpe, S. A. (1987), Current and temperature variability on the continental slope, Philos. Trans. R. Soc. London, Ser. A, 323(1574), 471–517. Wang, D.‐P., L.‐Y. Oey, T. Ezer, and P. Hamilton (2003), Near‐surface currents in DeSoto Canyon (1997–99): Comparison of current meters, satellite observation, and model simulation, J. Phys. Oceanogr., 33, 313–326. Waterhouse, A. F. (2006), A physical study of upwelling flow dynamics in long canyons, M.S. thesis, Univ. of B. C., Vancouver, B. C., Canada. Waterhouse, A. F., S. E. Allen, and A. W. Bowie (2009), Upwelling flow dynamics in long canyons at low Rossby number, J. Geophys. Res., 114, C05004, doi:10.1029/2008JC004956. Williams, W. J., E. C. Carmack, K. Shimada, H. Melling, K. Aagaard, R. W. Macdonald, and R. G. Ingram (2006), Joint effects of wind and ice motion in forcing upwelling in Mackenzie Trough, Beaufort Sea, Cont. Shelf Res., 26, 2352–2366. S. E. Allen and J. T. Dawe, Department of Earth and Ocean Sciences, University of British Columbia, 6339 Stores Rd., Vancouver, BC V6T 1Z4, Canada. (sallen@eos.ubc.ca; jdawe@eos.ubc.ca) DAWE AND ALLEN: SOLUTION CONVERGENCE OF FLOW OVER CANYON C05008C05008 13 of 13


Citation Scheme:


Citations by CSL (citeproc-js)

Usage Statistics



Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            async >
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:


Related Items