Open Collections

UBC Faculty Research and Publications

On vertical advection truncation errors in terrain-following numerical models: Comparison to a laboratory.. Allen, Susan E.; Dinniman, M. S.; Klinck, J. M.; Gorby, D. D.; Hewett, A. J.; Hickey, B. M. 2003

You don't seem to have a PDF reader installed, try download the pdf

Item Metadata


Allen_AGU_2001JC000978.pdf [ 1.13MB ]
JSON: 1.0041900.json
JSON-LD: 1.0041900+ld.json
RDF/XML (Pretty): 1.0041900.xml
RDF/JSON: 1.0041900+rdf.json
Turtle: 1.0041900+rdf-turtle.txt
N-Triples: 1.0041900+rdf-ntriples.txt
Original Record: 1.0041900 +original-record.json
Full Text

Full Text

On vertical advection truncation errors in terrain-followingnumerical models: Comparison to a laboratory model forupwelling over submarine canyonsS. E. Allen,1M. S. Dinniman,2J. M. Klinck,2D. D. Gorby1A. J. Hewett,1and B. M. Hickey3Received 14 May 2001; revised 27 December 2001; accepted 10 April 2002; published 3 January 2003.[1] Submarine canyons which indent the continental shelf are frequently regions of steep(up to 45C176), three-dimensional topography. Recent observations have delineated the flowover several submarine canyons during 2–4 day long upwelling episodes. Thus upwellingepisodes over submarine canyons provide an excellent flow regime for evaluatingnumerical and physical models. Here we compare a physical and numerical modelsimulation of an upwelling event over a simplified submarine canyon. The numericalmodel being evaluated is a version of the S-Coordinate Rutgers University Model(SCRUM). Careful matching between the models is necessary for a stringent comparison.Results show a poor comparison for the homogeneous case due to nonhydrostatic effectsin the laboratory model. Results for the stratified case are better but show a systematicdifference between the numerical results and laboratory results. This difference is shownnot to be due to nonhydrostatic effects. Rather, the difference is due to truncation errors inthe calculation of the vertical advection of density in the numerical model. The calculationis inaccurate due to the terrain-following coordinates combined with a strong verticalgradient in density, vertical shear in the horizontal velocity and topography with strongcurvature. INDEX TERMS: 4255 Oceanography: General: Numerical modeling; 3337 Meteorology andAtmospheric Dynamics: Numerical modeling and data assimilation; 4520 Oceanography: Physical: Eddies andmesoscale processes; KEYWORDS: model, sigma, error, topography, upwelling, canyonCitation: Allen, S. E., M. S. Dinniman, J. M. Klinck, D. D. Gorby, A. J. Hewett, and B. M. Hickey, On vertical advection truncationerrors 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, 2003.1. Introduction[2] Both laboratory and numerical models have been usedtostudyprocessesincoastaloceancirculation.Recentstudieshave compared test simulations from different numericalmodels with the disquieting result that the solutions differnot only quantitatively but, in some cases, qualitatively[Haidvogel and Beckmann, 1998]. Given the wide varietyof processes and parameterizations in these models, thecause(s) of these differences have not been identified.[3] The traditional method for verifying numerical mod-els is to compare to an analytical solution. However, only alimited number of analytic solutions exist for stratified flowover steep topography. An alternative evaluation method isto compare results from different models. However dis-agreements between the models frequently tend to raisequestions instead of answering them.[4] We chose to compare results from numerical simu-lations to those from a physical (laboratory) model of bothhomogeneous and stratified fluids flowing over a submarinecanyon on a continental shelf. This problem activates anumber of physical processes (discussed below) making it agood test case [Haidvogel and Beckmann, 1998]. Thecirculation we consider starts impulsively; the adjustmentto steady state is compared between the two models.[5] Laboratory models are thought to be better analogs forgeophysical scale flow than numerical models since a widerrange of length scales are active. Regional circulation isproduced by forces acting on scales from 100 km to 1 mm(the turbulence dissipation scale), a range of 8 orders ofmagnitude. A laboratory model has scales from 1 m to 1 mm,a range of three orders of magnitude. Numerical modelsrepresent scales from the domain size to the grid spacing, theratio of which is typically 100 (the number of grid points inonedirection),arangeoftwoordersofmagnitude.Onemightargue that sub-grid-scale parameterizations represent scalessmaller than the grid, but these formulations range in qualityfrom approximately correct to merely plausible but unveri-fied [Gent et al., 1995]. These problems with formulationsJOURNAL OF GEOPHYSICAL RESEARCH, VOL. 108, NO. C1, 3003, doi:10.1029/2001JC000978, 20031Oceanography, Department of Earth and Ocean Sciences, University ofBritish Columbia, Vancouver, British Columbia, Canada.2Center for Coastal Physical Oceanography, Old Dominion University,Norfolk, Virginia, USA.3School of Oceanography, University of Washington, Seattle, Wa-shington, USA.Copyright 2003 by the American Geophysical Union.0148-0227/03/2001JC000978$09.003 - 1areaddedtoissuesoftruncation errorsinnumericalschemes,domain distortion (e.g., bottom-following coordinates), andother parameterizations (e.g., equation of state, estimates ofbottom or surface stress, etc.). Therefore, laboratory modelsshould represent good solutions (assuming appropriate labo-ratory technique) which can be used as a stringent test fornumerical models.[6] The extrapolation of laboratory analogs to a geo-physical scale system involves measures of similarity.Similarity theory applied to fluid dynamics has yieldedrich results over the centuries [Batchelor, 1967]. Flowswith the same values of nondimensional ratios, whichcompare sizes of physical processes, will be identical otherthan having different scales. One need not match all non-dimensional numbers in such a comparison, only the‘‘important’’ ones. In many cases, ‘‘important’’ is easy todetermine; in others it is open to interpretation. For advec-tion over a canyon the important nondimensional numbersare a Rossby Number, a Burger Number, a Froude Numberand a geometric number describing the shape of the canyon[Allen et al., 1999].[7] First we present an overview of the flow regime underconsideration. Next the two flow representations (laboratoryand numerical) are presented and their results are comparedfor both homogeneous fluid and stratified fluids. Thediscussion includes an estimate of the effect of nonhydro-static terms in the laboratory results and an estimate of theeffect of errors in density advection on the numericalresults.[8] Submarine canyons are common features of the con-tinental shelf. A typical canyon cuts from the continentalslope into the continental shelf about halfway to the coast.Data sets including off-axis information are available overAstoria Canyon [Hickey, 1997], Lydonia Canyon [Butman,1983], Barkley Canyon [Allen et al., 2001] and CarsonCanyon [Kinsella et al., 1987]. Submarine canyons areregions of enhanced upwelling during upwelling favourableconditions and thus are important in evaluating cross-shelfexchange in coastal regions. The problem used to comparethe two models is a short duration upwelling event. Previousnumerical studies have shown qualitative agreement withthe observations [Klinck, 1996; Allen, 1996; She and Klinck,2000].[9] Afteraccelerationoftheshelfbreakcurrent(inthefielddue to the wind or shelf waves) a flow pattern forms over thecanyon with strong variation in the vertical (Figure 1). Inparticular: level 1, the near surface current is not affected oris only weakly affected by the canyon and travels straightover the canyon; level 2, flow just above the canyon rimmoves over the canyon, descends into the canyon, turnstoward the canyon head and crosses the downstreamcanyon rim near the head; level 3, slightly deeper flowover the continental shelf bends into the canyon near themouth, traverses up the canyon and crosses the downstreamcanyon rim near the head; and level 4, deep flow within thecanyon turns cyclonically. Flows numbered 2 and 3 aboveform the active layer (Figure 1). Observations [Hickey,1997; Allen et al., 2001] show a trapped eddy over theFigure 1. A schematic of flow over a canyon. The near surface and deep layer are shown explicitly. Thedeep shelf flow and the upwelling current together form the active layer. The shelf-break depth isapproximately 200 m in the field and is 2.2 cm in the lab. The near surface flow extends down about 150 mfrom the surface in the field and 1.6 cm from the surface in the lab. The near surface flow-passes over thecanyon with little deflection. In the active layer the flow is strongly constrained by the topography. Flow inthe deep layer turns cyclonically. [after Allen et al., 2001]3 - 2 ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODELcanyon (in flow numbered 2 above). The presence orabsence of the eddy is due to the amount of stretching nearthe canyon rim [Allen et al., 2001].[10] The parameters governing the advective flow over asubmarine canyon during an upwelling event are [Allen etal., 1999] (1) a Rossby number: Ro = U/fr where U is theincident velocity at the half length of the canyon and at thedepth of the shelf break Hs, f is the Coriolis parameter and ris the radius of curvature of the shelf break isobath; (2) aFroude number: U/NHswhere N is the buoyancy frequencyat the depth of the shelf break; (3) a Burger number: (NHs)/fL. The length of the canyon L is measured, at shelf breakdepth, from the shelf break to the head (last isobathdeflected significantly by the canyon); (4) a Geometricparameter Ge = W/Wsbwhere Wsbis the width of the canyonat the shelf break. The width of the canyon W is measured athalf length. The geometric scales, strength of stratificationand rotation rate for the laboratory model were chosen togive values of the nondimensional parameters similar tothose for Astoria Canyon (Table 1).2. Model Details[11] The test problem is largely determined by the labo-ratory setup, in which the flow is started with a rapid changein rotation rate of the table. The numerical model (SCRUMVersion 3.1, a terrain-following coordinate model, Song andHaidvogel [1994]) is configured as closely as possible to thelaboratory setup with an initial, impulse-like body forcestarting the flow. The body force is applied over onerotation period to match the time over which the rotationrate of the table is changed.[12] Two sets of lab experiments are considered, one withhomogeneous fluid, and the other with stratified fluid. Thelaboratory and numerical geometries using a homogeneousfluid are identical. However, in the stratified case thehorizontal dimension was increased by a factor of 10 toreduce the bottom slope and allow the simulation toproceed. Appropriate adjustments were made to otherparameters (detailed below) to preserve the relevant non-dimensional ratios.2.1. Laboratory Scale Model[13] Physical or laboratory models have their specificlimitations. Generally, the scale of the model impliesmoderate Reynolds number (order 1000). Either the viscouseffects are minimized or the molecular viscosity in the tankis assumed to model a uniform eddy viscosity in the realworld. The second major laboratory limitation is the neces-sity of exaggerating the vertical scale. The total verticalscale must be large enough to reduce viscous effects; yet thetotal horizontal scale is limited by the size of the apparatus.In some cases (here, in the homogeneous case), verticalexaggeration leads to nonhydrostatic effects that would notbe expected in the real world. Although the eventual goal ofour research is to model the real world, here where we arecomparing two models, the numerical model is configuredto use uniform molecular viscosity.[14] The laboratory experiments were performed on a 1-mdiameter, computer-controlled rare earth motor-driven rotat-ing table at the University of British Columbia. The outsidetank is a cylindrical Plexiglas tank of radius 50 cm. It wasfitted with a Plexiglas insert to give a generic coastalbathymetry consisting of a flat abyssal plane, a steepcontinental slope (slope 45C176 in laboratory) and a continentalshelf (slope 5C176 in laboratory). The continental shelf andslope were fitted to the outside of the tank with the abyssalplane in the center. Without the canyon, the topography hascylindrical symmetry. A 22C176 slice was cut out of the genericbathymetry to allow insertion of a canyon shape (Figure 2).The canyon was made of layers of 1 mm thick sheets ofpolystyrene. The layers were joined using plastic cementand the canyon model was covered with automotive fillerand sanded to produce a smooth finish and so join smoothlywith the generic shelf topography. The canyon was paintedflat black to allow flow visualization. The exact bathymetryTable 1. Scales for Astoria Canyon [Hickey, 1997] and theLaboratory and Numerical Models Considered HereNumber Astoria Lab NumericalHs150 m 2.2 cm 2.2 cmU 0.2 m sC011.2 cm sC011.2 cm sC01f 1.05 C2 10C04sC010.52 sC010.052 sC01L 21.8 km 8 cm 80 cmW 8.9 km 2.4 cm 24 cmr 4.5 km 1.4 cm 14 cmN 7.5 C2 10C03sC012.2 sC012.2 sC01Wsb15.7 km 6.9 cm 69 cmRo 0.42 1.6 1.6Fr 0.18 0.25 0.25Bu 0.49 1.2 1.2Ge 0.57 0.35 0.35Figure 2. Schematic bathymetric map of the topographyin the tank. Contours are 1 cm apart. The vertical scale isexaggerated compared to the horizontal scale by a factor of10 (radius of tank is 50 cm). The ‘‘shelf’’slope is 5C176 and the‘‘slope’’ slope is 45C176. The numerical geometry is similarexcept: (1) for the stratified simulations the horizontal scaleis increased by a factor of 10 (so there is no verticalexaggeration in the numerical model) and (2) only 60% ofthe tank is simulated.ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODEL 3 - 3was determined by digitizing photographs. Horizontal lightsheets were used to illuminate the canyon surface. The lineof illumination marked an isobath with depth dependent onthe height of the horizontal light sheet. Each pictureproduced one bathymetric contour and the contours weretaken 1 cm apart. This bathymetry was used in the numer-ical model.[15] For the homogeneous experiments, the tank wasfilled with fresh water and then accelerated to its startingrotation speed ( f = 0.4 sC01). The table was allowed to rotateuntil the fluid reached solid body rotation and was sta-tionary with respect to the tank (about an hour).[16] For the stratified case the table was first acceleratedto its starting rotation speed ( f = 0.4 sC01). Then the tank wasfilled by gravity feed through a tube coming down from theceiling along the axis of rotation. The tank was filled fromthe bottom using the Osler method of stratification [Osler,1965]. In a tank with straight sidewalls the Osler methodresults in a linear stratification. Because of the geometry ofthe shelf slope topography, here the stratification is strongernear the top (Figure 3). The table was allowed to rotate untilthe fluid reached solid body rotation and was stationary withrespect to the tank (several hours). The tank stratificationwas replaced each day to reduce the formation of ‘‘mixed’’layers at the top and bottom of the tank.[17] The flow was forced by accelerating the table fromf = 0.4 sC01to f =0.52sC01linearly over 27.3 s. This periodis about one revolution of the tank and is comparable to awind event of one day in the real world. The change inrotation speed causes the bathymetry to accelerate relativeto the water, which due to inertia continues to rotate at theslower rate. Relative to the tank, the water moves withshallow water to its left. This forcing results in flow at theshelf break of about 1.2 cm sC01.[18] The flow in the vicinity of the canyon was visualizedusing neutrally buoyant particles (pliolite) which wereadded upstream of the canyon during the acceleration phaseof the experiment. The tank was lit by a slide projector witha bicolored slit. The projector was mounted above the tankattached to the table (i.e., it rotated with the table). The lightbounced off a 45C176 mirror set in the abyssal plain near theaxis of rotation. Two horizontal sheets of light of 1 cmthickness resulted. The green light illuminated depths of 12to 22 mm and the red light illuminated depths of 22 to 32mm. A Hi-8 video camera mounted above the tank (rotatingwith the table) recorded the movement of the particles.[19] The flow was analyzed for a 5 s interval, 2.7 s afterthe acceleration of the tank was complete. Using framecapture software, frames from the video were captured ontocomputer at 0.5 s intervals. Particles were traced manuallyfrom frame to frame. The experiment was repeated 5–15times. Composite photos of all the frames for the homoge-neous case and the stratified case were generated anddigitized versions of these composite figures are used formodel comparison.[20] Two major sources of error occur in the laboratory.First, motion in the initial flow; second, vertical motion ofthe pliolite with respect to the fluid. Errors are likely to bedifferent between different runs of the model. To estimateerror, two tracks very close in space but from different runswere compared (see section 4, Figure 7b.) The locationdifference after 2.5 s is 0.5 cm, with a resultant velocityerror of 0.2 cm sC01. Other minor sources of error include thedepth of the light sheets (estimated as plus or minus 0.1 cmtop and bottom), and errors in identifying and digitizing theparticles (estimated as plus or minus 0.1 cm).2.2. Numerical Model Configuration[21] The S-Coordinate Rutgers University Model Version3.1 [SCRUM, Song and Haidvogel, 1994] developed atRutgers University was used for the numerical modelingpart of this study. SCRUM is a hydrostatic primitiveequation ocean circulation model with a free surface usinga vertical s (terrain-following) coordinate well suited for usein simulations with variable bathymetry. The model alsoallows for a curvilinear horizontal coordinate system whichenabled us to represent a cylindrical model geometry.[22] The combination of steep topography and strongstratification can cause errors in numerical simulations[Haney, 1991]. The errors, usually in the form of meanalong isobath currents, are due to truncation errors in thecalculation of the horizontal pressure gradient term interrain-following coordinates. One current method to assesserrors is to initialize the model with horizontal isopycnalsand no forcing and evaluate the size of the observed errorflows [e.g., Petruncio, 1996]. The resulting flows are due topressure gradient errors and the boundary condition ondensity at the bottom which requires isopycnals to haveno normal gradient at the bottom. The net effect of pressuregradient errors on strongly stratified flows over steep top-ography have been successfully assessed using laboratoryexperiments [Kliem and Pietrzak, 1999]. Steep topographyand stratification can lead to other errors in s coordinatenumerical models (as we show here). Such errors cannot beidentified without forcing and have to date received littleattention.[23] The large slopes of the tank canyon were too steepfor the model to accurately compute horizontal pressuregradients. Therefore, for the stratification case only, thelength scale was increased to ten times that of the tank. ToFigure 3. Density profile in the tank and numericalsimulation. The initial density profile (solid) is calculatedassuming no extra mixing in the filling process. Theexpected diffusion after 1 hour of spin-up is shown by thedashed line. The numerical simulation used the solid profileas its initial condition.3 - 4 ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODELuse the same bathymetry, stratification and flow velocity asthe tank, the Coriolis parameter was reduced (from 0.52 sC01to 0.052 sC01) and the vertical viscosity was reduced (fromthe molecular value of 10C06m2sC01to 10C07m2sC01)byafactor of 10 in order to maintain the same Rossby, Burgerand Ekman numbers. The duration of the forcing wasincreased by a factor of 10 so that the model, like the tank,was accelerated over approximately one pendulum day.[24] The model domain included 242 points in the azi-muthal direction, 81 points in the radial, and 16 verticallevels. The grid spacing in the vicinity of the canyon wasabout 6.2 cm (Figure 4). The domain did not cover theentire 2p radians in the azimuthal direction, but rather 0.6 C22p (in order to save computation time) with a periodiccontinuation in that direction. The duration of the simula-tions was much less than the time required for disturbancesto wrap around the domain. In the radial direction, themodel domain was based on a tank with an inner wall 1 mfrom the center, an outer wall 5 m from the center (10 C2 theactual tank dimensions) and a no-slip condition for momen-tum at both lateral walls. The gridded model bathymetrywas derived from the digitized tank depth contours using anobjective analysis scheme.[25] The Laplacian horizontal mixing was along geo-potential surfaces for both momentum and density with acoefficient of 10C06m2sC01. Vertical mixing of momentumand density used constant coefficients of 10C07m2sC01and1.5 C2 10C010m2sC01(one tenth of the value for moleculardiffusivity of salt in water [Weast, 1979]), respectively. Asimplified linear equation of state dr = C0dT means thattemperature was used as an inverted proxy for density. Forthe stratified simulation, the initial vertical density profilewas the same as for the tank (Figure 3) and the startingdensity was horizontally uniform.[26] At the start of the simulation the flow was zero andthe free surface was level. The water was forced to moverelative to the canyon by a body force over the entire watercolumn. The forcing was zero at the start, increased to amaximum value over six seconds, held constant until 273seconds and then set to zero. The forcing was in the angulardirection only and increased linearly with radial distance inorder to simulate the situation in the tank where the water isinitially in rigid body rotation before the tank is acceleratedwith respect to the rotating water. Linear (Ekman) bottomstress was applied as a body force over the bottom layer(drag coefficient of 7.2 C2 10C05msC01).3. Homogeneous Comparison[27] To compare the numerical and laboratory results,tracer particles were introduced into the numerical simula-tion at the time (0.063 rotation periods after accelerationwas complete) and place tracer particles were observed tostart in the laboratory experiment. Similarity and differencesbetween the flows were determined by comparing the tracksof these particles (Figure 5). A value of 1 cm was chosen todelineate ‘‘good’’ and ‘‘poor’’ agreement. This value wasdetermined using the estimate of the error in the laboratoryof 0.2 cm sC01and the maximum track durations of 5 s. Formost tracks the final position of the particle according to thenumerical model is within 1 cm of the final positionaccording to the laboratory model (Figure 5a). The tracksshow good agreement upstream of the canyon and fairagreement both across the mouth of the canyon and nearthe head of the canyon.[28] All tracks which finish more than 1 cm apart aredownstream of the canyon inshore of the shelf break. Inthese cases the numerical tracks show stronger upwelling(stronger flow toward the coast) than the laboratory tracks(Figure 5b). The discrepancy suggests that greater cyclonicvorticity is generated within the canyon in the numericalmodel than in the laboratory model. All numerical tracks arebent strongly inshore over the canyon by this vorticity.[29] One explanation for the smaller vorticity generationin the laboratory comes from careful inspection of thevertical particle movement. The laboratory particles arecolor coded based on their depth. Upstream of the canyon,green particles become red particles as the flow drops downinto the canyon (Figure 5b). On the downstream side, redparticles become green particles as the flow is stronglyupwelled. These numerous changes in color illustrate thestrong vertical velocities occurring in the homogeneouscase. We suspect that because the laboratory model isstrongly vertically exaggerated nonhydrostatic flows occurin the laboratory for homogeneous flows.[30] To examine more carefully whether nonhydrostaticeffects are indeed occurring in the laboratory model avisualization of the vertical shear downstream of the canyonwas performed. In a homogeneous flow the pressure gra-dient is uniform with depth (barotropic) and vertical shearcan only be generated by differential body forcing (notpresent) or in boundary layers. However, besides a clear andFigure 4. The numerical grid for a cross-section throughthe canyon. The vertical scale is in cm. The horizontal scaleis in grid points with C1x C25 6.2 cm for this cross-section.The deepest cell at grid point 37 is used for the errorcalculation in section 5.ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODEL 3 - 5expected bottom Ekman layer, a middepth jet was observedwith a strong up-canyon component below a surface layerflowing more slowly downstream (not shown). Thus mid-depth vertical shear was strong in the laboratory consistentwith the presence of nonhydrostatic effects.[31] The numerical model is hydrostatic and thus con-strained to have no vertical shear above the bottom Ekmanlayer. Fluid columns passing over the canyon are stretched tothe full depth of the canyon, producing the observed strongcyclonic vorticity. In the laboratory model, columns passingover the canyon are only more moderately stretched.4. Stratified Comparison[32] To compare the numerical and laboratory results,tracer particles were placed in the numerical simulation atthe time and position tracer particles were observed to startin the laboratory experiment. In the stratified case, verticalshear is significant (Figure 6). In the laboratory, particledepths were determined, by color, only to within 1 cmdepth. This resolution proved inadequate because of thestrength of the vertical shear. Therefore in the numericalsimulation, particles were seeded throughout the 1 cm depthinterval and the track that finished closest to the laboratorytrack was chosen for comparison.[33] As discussed in the overview (section 1.1), flow inthe canyon can be divided into three layers (Figure 1). Near-surface flow is not affected or only weakly affected by thecanyon. Flow deep within the canyon turns cyclonically.Active upwelling onto the shelf occurs between these twolayers. The upwelling layer includes water from both theshelf and slope which turns cyclonically in the vicinity ofthe canyon, flows up the canyon and exits the canyon nearthe head. The discussion of the stratified results will be inreference to these layers.[34] The deep flow (below 2.15 cm) shows the expectedtrapped cyclonic eddy within the canyon (Figure 7d). Thesize of the eddy appears to be different between thenumerical and laboratory realizations: the eddy is widerand longer in the lab, and shorter and narrower in thenumerical simulations. Also in the region of the mouth ofthe canyon, the laboratory flow appears to slow more thanthe numerical flow.[35] Agreement between the numerical simulation andthe laboratory experiment is generally very good in theactively upwelling layer (1.50–2.05 cm depth, Figures 7band 7c). The flow shows upwelling flow over the down-stream rim and near the head of the canyon. Two system-atic differences are apparent. Upstream of the canyon overthe slope the laboratory tracks show weak downwellingflow (e.g., Figure 7b, track starting at x = C01, y = 2). Thisflow is not observed in the numerical simulation where thetracks follow the curvature of the tank. The difference inthe tracks is about 0.5 cm and occurs over about 2 s givinga velocity error of 0.2 cm sC01. Downstream of the canyonthe opposite situation occurs with weak downwelling flowin the numerical simulation but not in the laboratory flow(e.g., Figure 7b, track starting at x = C04, y = 2). The tracksof the numerical floats are consistent with the numericalvelocity so these differences are not due to the Lagrangiantracking scheme. Flow above the canyon (between 1.25 cmto 1.35 cm depth, Figure 7a) shows the same systematicdifference upstream of the canyon as seen in the activelayer.Figure 5. A plan view of the particle tracks in the laboratory experiment (unfilled red and greensymbols) and the numerical simulation (filled blue symbols) for the homogeneous case. The depthinterval of the laboratory tracers was determined by the color-sheet illuminating them. The laboratoryparticles here are color-coded by depth, green: 1.2–2.2 cm depth, red: 2.2–3.2 cm depth. Thebathymetric contours are shown in blue with depth in cm. The axes are in cm in the tank frame. The lastposition is shown as a larger box. (a) Tracks that finish within 1 cm of each other. (b) Tracks that finishmore than 1 cm apart. Note that these tracks are all downstream of the canyon. The weaker upwellingobserved in the laboratory case is due to non-hydrostatic effects (see text).3 - 6 ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODEL[36] To quantitatively measure the fit between the numer-ical results and the laboratory results, a normalized errorwas calculated for each particle at each point:Error ¼Position numericalðÞC0Position laboratoryðÞjjLabFluctuationVelocity C3 timeð1Þwhere the position is the (x, y) position in lab coordinates ofthe point, the laboratory fluctuation velocity is the velocityerror between runs estimated above at 0.2 cm sC01and timeis the time since the particle started (Table 2).[37] The deep point near the mouth of the canyon that iswithin the deep eddy in the lab and not in the numericalsimulation has an error at the end of 6.05. This point isvery close to a flow separation point, dividing flow travel-ing into the canyon (as the laboratory tracer does) and flowtraveling along the slope (as the numerical tracer does).Removing this point, the remaining points have an averageendpoint error of 0.7 and an average error of 0.8. If theerror estimate for the laboratory measurement was accurateand the model is reasonably accurate, a value near 1 wouldbe expected.5. Discussion[38] The utility of a numerical laboratory model compar-ison lies in resolving differences between them. In thefollowing, we demonstrate the limitation of the numericalmodel simulating stratified shear flows over steep topog-raphy. The largest errors here are not due to the pressureterms as investigated by previous authors [e.g., Haney,1991; Kliem and Pietrzak, 1999] because for our gridgeometry these errors are small (see below) but are due toerrors in vertical advection of the density field.[39] Although the average error in the stratified experi-ments is less than the estimated laboratory error, the resultsshow systematic differences. These differences are mostobvious near the upstream rim of the canyon where thelaboratory flow is downwelling and the numerical flow isnot. Downwelling flow in this region is observed in the field[Hickey, 1997; Allen et al., 2001]. A detailed search for theorigin of the differences was made. Running the model withno forcing showed that pressure gradient errors are 50 timestoo small to explain the error. The timescale of the simu-lation is so short that errors due to diffusion are also toosmall to explain the observed differences. The details of thetracer advection scheme (fourth-order centered or upstream-biased third-order (Gamma scheme)) make little difference.Below we present an estimate of the nonhydrostatic effectsin the laboratory and an estimate of the error in the verticaladvection of density in the numerical simulation. Nonhy-drostatic effects are shown to be too small and of the wrongsign to be the origin of the error.5.1. Estimate of Nonhydrostatic Effects[40] The homogeneous case presented in section 3 illus-trates that nonhydrostatic effects can sometimes explain thedifference between laboratory and numerical models. In thestratified case however the strong stratification should beenough to prevent nonhydrostatic flows (conditions calcu-lated by Boyer et al. [2000] are met). Here we willcarefully evaluate the possible strength of the nonhydro-static terms in the momentum equation for the stratifiedcase.Figure 6. A cross-section of the cross-canyon horizontal velocity at about mid-canyon for the stratifiedcase. Positive velocities are to the right. View is toward the head of the canyon. Note the large verticalshear over the canyon. The vertical scale is 5 cm from top to bottom. Units are cm sC01.ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODEL 3 - 7[41] The vertical momentum equation isro@w@tþ u@w@xþ v@w@yþ w@w@zC18C19¼C0@p@zC0 r0g þ m@2w@x2þ@2w@y2þ@2w@z2C18C19ð2Þwhere rois a constant reference density, (u, v, w)isthevelocity, p is the pressure, r0is the fluctuation density, g isthe gravitational acceleration and m is the viscosity.[42] Scale analysis shows that the third largest term(after the pressure and gravitational terms) is the horizontaladvection of the vertical momentum. After taking thehorizontal derivative of equation (2) with respect to xand retaining only the three largest terms, equation (2)becomesC01ro@2p@z@x¼gro@r0@xþ u@2w@x2ð3ÞFigure 7. Plan views of the particle tracks in the laboratory experiment (unfilled symbols) and thenumerical simulation (filled symbols) for the stratified case. The bathymetric contours are shown withdepth in cm. The axes are in cm in the tank frame. The ends of the tracks are shown as larger boxes. Notethe generally good agreement. However, there are systematic errors upstream and downstream of themouth of the canyon. In particular in the active layer and near surface layer at x = 1 cm, y =1–2cmthelaboratory tracks show downwelling flow compared to the numerical tracks. (a) Tracks between 1.25 and1.35 cm depth. (b) Tracks between 1.50 cm and 1.60 cm depth. (c) Tracks between 1.85 cm and 2.15 cmdepth. (d) Tracks between 2.15 (just above rim depth) and 3.15 cm depth.3 - 8 ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODELwhere the first term on the right is due to the hydrostaticpressure and the second term is the nonhydrostatic effect.[43] The largest terms in the cross-shore momentumequation are the pressure gradient and Coriolis terms (thatis, the flow is approximately geostrophic). Substitutingequation (3) into the z derivative of the cross-shore geo-strophic balance gives@v@z¼C0gf ro@r0@xC0uf@2w@x2ð4ÞSo the nonhydrostatic effect on the velocity shear isC0uf@2w@x2ð5Þ[44] To estimate the magnitude of the nonhydrostaticeffects in the laboratory where w was not measured weneed a magnitude for the vertical velocity. The numericalvertical velocity (Figure 8) is in error (see below). A localminimum in vertical velocity is expected over the upstreamrim of the canyon but the maximum should not occur untilthe downstream rim of the canyon based on observedisopycnals over Astoria Canyon [Hickey, 1997]. Estimatingthe vertical velocity as the radial velocity (u =1cmsC01)multiplied by the slope of the canyon edge in the tank givesan upper estimate for the magnitude of the vertical velocity(w C25C01cmsC01). Using the width of the canyon as anappropriate scale for x implies equation (5) is of orderC00.04 sC01. Because v is near zero near the surface, thenonhydrostatic effect gives positive (upwelling) flows, theTable 2. Initial Positions, Final Depth Position (Numerical) andErrors Between the Numerical and Laboratory SimulationsaPointNumber X Start Y Start Z Start Z Finish Aver Error End Error1 C03.8 5.3 3.2 3.2 0.99 0.8113 4.6 0.2 3.2 3.2 1.26 1.248 C00.7 2.5 3.2 3.1 2.63 2.1411 0.8 0.6 3.2 3.1 0.82 1.0427 C00.7 C04.4 3.2 3.0 0.79 0.2812 2.9 0.4 3.0 3.0 1.24 1.0323 4.4 C00.1 3.0 3.0 0.32 0.4718 1.8 0.1 2.5 2.4 0.78 0.5525 0.0 C02.0 2.2 2.2 5.13 6.0530 1.6 C05.0 2.2 2.2 0.42 0.726 4.5 5.5 2.2 2.1 0.38 0.197 C01.3 2.2 2.2 2.1 0.91 1.0526 1.6 C01.1 2.2 2.1 0.66 0.6017 1.1 C00.2 2.0 2.1 1.08 0.412 0.6 5.2 2.0 2.0 0.48 0.4928 0.5 C03.6 2.0 2.0 0.71 0.5133 3.5 0.8 1.9 1.9 0.17 0.1421 3.6 C00.9 1.9 1.8 0.29 0.1322 4.5 C00.2 1.9 1.8 0.53 0.8619 2.1 0.4 1.6 1.6 0.39 0.2824 5.0 C01.0 1.6 1.5 0.28 0.0629 1.3 C05.1 1.5 1.6 0.69 0.084 2.2 4.2 1.5 1.5 0.98 0.9220 3.5 C01.3 1.5 1.5 0.30 0.1215 6.1 1.0 1.5 1.4 0.63 0.243 1.2 6.3 1.4 1.4 0.26 0.3214 5.2 1.4 1.3 1.2 1.36 1.135 2.3 5.4 1.2 1.3 0.50 0.5031 2.7 C06.2 1.2 1.3 3.09 2.779 3.6 3.9 1.2 1.2 0.73 0.6510 C03.9 1.1 1.2 1.2 0.69 0.4516 C04.2 C00.7 1.2 1.2 0.13 0.14aData are sorted in-depth order.Figure 8. A cross-section of numerical model true vertical velocity at about mid-canyon for the stratifiedcase. Positive velocities are upward. The vertical scale is 5 cm from top to bottom. Units are cm sC01in thenumerical model. Note the small scale of the variations compared to the width of the canyon.ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODEL 3 - 9opposite of the observed effect. Assuming 1 cm is anappropriate depth for the strong downward flow gives anup-canyon, nonhydrostatic flow of 0.04 cm sC01. This erroris one order of magnitude too small (the observed differenceis 0.2 cm sC01). Since nonhydrostatic effects are the onlyreason the laboratory results could be ‘‘wrong’’, we nowturn to the numerical model to find the origin of the lab-numerical difference.5.2. Vertical Velocity and Vertical Advection[45] The numerical vertical velocity (Figure 8) showsunexpected spatial small-scale features. In particular theupward vertical velocity, just inside the rim of the canyon,is counter to our physical understanding based on observa-tions of isopycnals over Astoria Canyon [Hickey, 1997].The downward velocities at the downstream rim are alsoquestionable. We expect similar error near both rims buthere will concentrate on the upstream rim.[46] Just below the upstream rim the numerical horizon-tal cross-canyon velocity is about 0.2 cm sC01(Figure 6) theslope is about 1 in 20, the up-canyon velocity is about 0.05cm sC01and the slope is about 1 in 60. Thus from thebottom boundary condition (no flow through the canyonwall) we would expect a negative vertical velocity of about0.01 cm sC01, not the positive velocities observed.[47] The vertical velocity is not a dynamic variable withinthe model; it is calculated as a diagnostic variable. Temporalfluctuations can thus be large but here are less than 1%. Thevertically important dynamic within the model is the verticaladvection of density (Figure 9). The calculated isopycnalsreflect the calculated vertical velocity; they tilt gently downbefore the canyon and then abruptly upward over thecanyon. The total deflection of the isopycnals in the numer-ical model results is very small (about one tenth the magni-tude) compared to observations [Hickey, 1997]. If we assumethe isopycnals are tilted in the laboratory as they are in thefield we can estimate the generated down-canyon flow usingthe thermal wind equation (N2Di2/fDx) where the isopycnaldrop Diobserved in the field is about 25 m which wouldcorrespond to about 0.2 cm in the lab, and Dx is a lengthscalefor the drop taken as one quarter of the canyon width. Theresulting flow is 0.2 cm sC01in agreement with the observedflow error.[48] The calculation of the advection of temperature inthe numerical model uses the Shchepetkin and McWilliams[1998] Gamma advection scheme to reduce dispersive error.This scheme uses a large number of terms. The truncationerror in the vertical velocity can be calculated for thisscheme but the algebra is dense. To more easily illustratethe error, we will use the diagnostic calculation for thevertical velocity instead. The derivation for the Gammaadvection scheme is in the Appendix. In the following, (1)the link between density advection and vertical velocity willbe shown mathematically. Then (2) the truncation error inthe vertical velocity will be calculated, followed by (3) theresults (only) of the Gamma advection scheme.[49] The calculations will use the numerical values for a10 m diameter tank. However, unlike the nonhydrostaticcalculation above, the calculation of the advection error cansimply be scaled to give values for a field simulation.Factors (based on Astoria Canyon) are 23 C2 103,6.8C2103, 17, 5.0, for horizontal distance, vertical distance,Figure 9. A cross-section of the density at about mid-canyon. The vertical scale is 5 cm from top tobottom. Units are st. Very small isopycnal deflections are seen compared to field observations [Hickey,1997]. In particular, over the canyon near the upstream rim (yellow band) isopycnals tilt upwards to theright in a region observed to have strongly descending isopycnals.3 - 10 ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODELhorizontal velocity and vertical velocity, respectively. Per-centage errors are the same, regardless of scale.5.2.1. Step One: Density Advection and VerticalVelocity[50] Consider a case in which the density (temperaturehere) only varies vertically. Then the time variation oftemperature T is given by@T@t¼C0w@T@z: ð6ÞThe advection of temperature is calculated within the modelas@T@t¼C0u@T@xC0 v@T@yC0C10@T@sð7Þwhere s is the terrain-following vertical variable and C10 =Ds/Dt. Equation (7) can be rewritten as@T@t¼C0u@z@x@T@zC0C10@z@s@T@zð8Þwhere variation in y has been neglected and temperature hasbeen assumed to only vary vertically. The problem withequation (8) is that the two terms on the right-hand sidealmost cancel, so that the calculation of @T/@t is based onthe small difference between two large terms. Combiningterms we can relate the density advection in the model to thevertical velocity in the following way:@T@t¼C0 u@z@xþC10@z@sC18C19@T@zð9Þwhere the term in brackets is the model vertical velocity asshown in equation (6). The numerical model does not usedT/dz (only @T/@x and @T/@s) and so equation (9) is notcomputed. However, this equation illustrates how the modelcalculated vertical velocity is related to the actual advectionof a vertical temperature gradient. The vertical velocity isfound from the small difference between the two large termsu@z/@x and C0C10@z/@s. The actual calculation of theadvection of temperature in the model is complicated bythe chosen discretization. However, the error identified hereis due to the estimation of the vertical velocity. The smalldifferences due to the discretization will be discussed instep Step Two: Error in the Vertical Velocity[51] The second step is to estimate the truncation error inthe model vertical velocity given by equation (9). The flowover the canyon is about 0.5 cm sC01with a vertical shear ofabout 1.6 sC01. Only the shear contributes to the error in thevertical velocity, so we will set the background flow to zero.This choice reduces the magnitude of the terms in equation(9) so the fact the problem is the small difference betweentwo large terms will be partially obscured. For reference,with the observed background flow u@z/@x C25C00.02 cm sC01,from the model w =C10@z/@sC250.01 cm sC01and so the verticalvelocity is aboutC00.01 cm sC01. We will assume that the flowat the bottom is parallel to the topography and that the flowin the up-canyon (y) direction is zero. Variations in y will beneglected.[52] Consider a single grid cell (Figure 10) near theupstream rim of the canyon (Figure 4). The velocity at thelower depth is taken as zero, and the velocity at the upperdepth is U whereU ¼C0@u@z@z@xC1x ð10ÞThe diagram shows that @z/@x = C0tan C22a.[53] The vertical sigma velocity is w = C10Hzwhere C10 =Ds/Dt and Hz= dz/ds. The value of w at the top of the cell isfound fromw ¼C0ZrC1 uHzðÞds ð11Þwhich, using the expansion formula [Hedstro¨m, 1997], andassuming v = 0 giveswjs¼ wjsC0C1sðÞC0udzðÞjxþC1x=2ðÞC0 udzðÞjxC0C1x=2ðÞC16C17C1x¼ UC1ZC1xð12Þ[54] The vertical velocity is calculated at the centre pointasw ¼wjsþC1s=2ðÞþwjsC0C1s=2ðÞC16C172þu@z@xC0C1C12C12xþC1x=2ðÞþ u@z@xC0C1C12C12xC0C1x=2ðÞC16C172¼UC1Z2C1xþU2@z@xC12C12C12C12xC0C1x=2ðÞð13ÞThese two terms are of opposite sign near the upstreamcanyon rim. The first, C1Z/C1x is positive whereas thesecond, @z/@x is negative. Substituting for U gives theexpression for the numerical wwnum¼12@u@ztan C22a C1Z C0C1xtana‘ðÞð14Þwhere a‘= 1/2 (a1+ a2).[55] Now w at the bottom (here taken on the straight linebetween the two u grid points) is C0u tan abwhere ab=1/2(a2+ a4) and as there is no horizontal change in u, w isconstant with depth. Sowreal¼12@u@zC1Z0C0C1xtanabðÞtanabð15ÞThe correct velocity is also found from the difference of twoterms. There is truncation error in both numerical termswhich becomes even more significant because the verticalvelocity is the difference between them.[56] The fourangles for thegridcell at theupstream cornerare: a1= 0.0094, a2= 0.012, a3= 0.086, a4= 0.105. Thesize of this grid cell is C1 x = 6.2 cm by C1Z = 0.31 cm.The vertical shear was estimated from the numerical results(Figure6)as@u/@z=1.6sC01givingwnum=0.0103cmsC01.Thelength C1Z0= C1Z + C1x/2 [tan(a2) + tan(a4) C0 tan(a1) C0tan(a3)] is 0.373 cm and wreal= 0.00047 cm sC01. This differ-ence results in an upward velocity error of 0.010 cm sC01ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODEL 3 - 11in the same direction as the observed error. The error in thevertical velocity is of the same size as the actual velocity.This positive vertical velocity error at the upstream canyonedge gives a positive density anomaly and therefore positive(upwelling) flow error upstream of the canyon.[57] The grid cell at the upstream corner and the corre-sponding grid cell at the downstream corner are expected togive the largest errors in the domain because they havestrong topography curvature (angles on the right are differ-ent from the angles on the left) and are in regions of strongvertical shear.[58] Equations (14) and (15) can be combined and assum-ing C1x and C1Z are sufficiently small, the error can beshown to be of second order in these quantities. However, itis important to note that for the case considered here C1x andC1Z are too large for the error to be dominated by the leadingterms.5.2.3. Step Three: Error in Advection[59] We have shown that (1) the vertical velocity and thevertical advection of temperature are linked and (2) thatthere are large truncation errors in the calculated verticalvelocity calculation. However, the vertical velocity calcu-lated above is not directly used by the model in theadvection of temperature. One can, however, similarlycalculate the equivalent vertical velocity used in the Gammascheme for advection of a vertical gradient in temperature(see Appendix A). This givesC0@T@t¼@T@z0:007 cm sC01C0C1ð16ÞFigure 10. Sketch of a single grid cell with a sharp change in orientation. This grid cell sits at theupstream rim of the canyon (Figure 4). Four angles (a1–a4) and the vertical (C1Z) and horizontal (C1x)grid sizes characterize its geometry. Horizontal flow is sheared with U at the higher left-hand grid-edgeand 0 at the lower right-hand grid-edge. The vertical velocity is w and the cross-grid ‘‘vertical’’ velocityis w.3 - 12 ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODELThe error is smaller (by 30%) but still very significant atabout 70% of the expected real vertical velocity (includingthe background flow).[60] To illustrate the effects of the truncation error onvertical velocity the model was run with two finer grids, onewhich does not reduce the size of the changes in bottomslope (Figure 11a) and one which does (Figure 11c). Thefirst finer grid has a horizontal grid size one half that of theoriginal grid (Figure 4) but the same bottom topography.The smaller horizontal grid size leads to smaller scales inthe vertical velocity (Figure 11b) and no decrease in themagnitude of the positive velocities over the upstream rim.Close examination of equations (14) and (15) reveals thatthe difference is primarily dependent on the change ofbottom angle. The first fine grid still contains a cell withthe same angle change and thus the same large error invertical velocity.[61] The second finer grid slightly smoothes the bottomtopography (Figure 11c). The maximum change in angle isreduced by 25%. The vertical velocity (Figure 11d) showsa similar reduction in error. However, the vertical velocityerror is still very large. The density field and particletrajectories are not significantly altered. For the currentconfiguration, reducing the vertical grid size wouldFigure 11. Grids and vertical velocities near the upstream rim for two grids with smaller grid spacing inthe horizontal. (a) Grid with no smoothing. (b) Vertical velocity for grid with no smoothing. (c) Grid withsmoothing of bottom topography. (d) Vertical velocity for grid with smoothing. Vertical velocity units arecm sC01(positive upward) in the numerical model. Note the smaller vertical velocities in the smoothercase.ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODEL 3 - 13actually increase the error. Although the error is related togrid size and the numerical scheme does converge as thegrid size approaches zero, reducing the error by reducingthe grid size is onerous because both the vertical andhorizontal grids must be reduced. Other possibilities forreducing the error are smoothing the topography or usingan advection scheme that accounts for strong verticalshear. Calculation of the vertical advection errors for thesmooth broad canyons used by Perenne et al. [2001] showthat their truncation errors are small, consistent with thegood match they found between laboratory and numericalresults. However, canyons in the ocean have steep, sharptopography.6. Conclusions[62] Upwelling over submarine canyons has been inves-tigated using both a numerical and a physical model.Although qualitative results are similar, the results showsystematic differences between the models. In the homoge-neous case, the differences can be accounted for by non-hydrostatic effects in the laboratory model.[63] In the stratified case, nonhydrostatic effects are toosmall and have the wrong sign to account for the observeddifferences. However, the numerical truncation error in thecalculation of the vertical advection was shown to be largeenough and to have the correct sign to account for thedifferences. The presence of this error requires (1) strongvertical gradients in density, (2) vertical shear in thehorizontal velocity, (3) topography with strong verticalcurvature, and (4) the use of terrain-following coordinates.As the first three conditions are typically found in coastalenvironments, a solution to this error will be important forsimulations that use realistic topography. In the modelconfiguration used here, the error in the vertical advectionis calculated to be 70%. As the canyon shape used is similarto Astoria Canyon the error for field-scale numerical sim-ulations is the same, 70% in vertical advection. A similarerror has been found by the authors in realistic simulationsof Astoria Canyon. Because the error is dependent on boththe vertical and horizontal grid size, reduction in grid size isan expensive approach to the problem. Higher order esti-mates of the vertical velocity that account for the verticalshear are likely required.Figure A1. Expanded view of grid in the vicinity of the upstream rim of the canyon. The center of thecell used for error estimation is marked by T. The geometry of the grid is needed farther upwards andupstream for the vertical advection calculation. Values of the angles and grid spacing are given. Note thatthe vertical scale is exaggerated.3 - 14 ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODELAppendix A: Advection of a Vertical Gradient inTemperature[64] In section 5.2, an error in the vertical velocity wascalculated. Here we will calculate the dynamically impor-tant error, the error in the vertical density advection. Con-sider the same grid cell. To calculate vertical advection, wemust consider flow further upstream and further up. Con-sider the extended mesh shown in Figure A1.[65] The change in temperature at the center point(marked T) should be@T@t¼C0wreal@T@zassuming that the temperature only varies vertically andwrealis given by equation (15). The spatial discretization forthe calculation of the density advection is:C0@@tHzTðÞ¼1C1xuxþC1x2;sC18C19C1Hzx þC1x;sðÞþHzx;sðÞ2C1TxþC1x;sðÞþTx;sðÞ2C01C1xuxC0C1x2;sC18C19C1Hzx;sðÞþHzx C0C1x;sðÞ2C1Tx;sðÞþTxC0C1x;sðÞ2þ1C1sC10 x;s þC1s2C18C19C1 Hzx;s þC1s2C18C19Tx;s þC1sðÞþTx;sðÞ2C01C1sC10 x;s C0C1s2C18C19C1 Hzx;s C0C1s2C18C19Tx;sðÞþTx;s C0C1sðÞ2ðA1Þneglecting the additional ‘‘Gamma’’ term in the Shchepet-kin-McWilliams advection scheme. It will be discussedbelow. To calculate the value for the center of the grid cellwe will ignore the uniform horizontal velocity and includeonly the vertical shear. Thus without loss in generality weassume u(x + C1x/2, s) = 0 and u(x C0 C1x/2, s)=U.Asw(x, sC0 C1s/2) is on the topography it is zero and as before C10Hz=w. We note that 0.5 (Hz(x C0 C1x, s)+Hz(x, s)) = C1Z/C1s,Hz(x, s)=0.5(Z + Z0)/C1s, and that equation (12) gives w(x,s + C1s/2) = UC1Z/C1x.[66] Then equation (A1) becomesC0C1Z þC1Z02C1sC18C19@T@t¼C0UC1Z2C1xC1sTx;sðÞðþTxC0C1x;sðÞÞþUC1Z2C1sC1xTx;s þC1sðÞðþTx;sðÞÞ ðA2Þ[67] Using the geometry shown in Figure A1, we canexpandU ¼@u@ztan C22aC1xTx;sþC1sðÞ¼Tx;sðÞþ12dTdzC18C1ZþC1Z1þC1x2x tana2C0 tana0ðÞC19TxC0C1x;sðÞ¼Tx;sðÞþdTdzC1x2tana1þ tana2ðÞwhere a0and C1Z1are defined on Figure A1. ThusC0@T@t¼14@u@ztan C22a@T@zC1 C1Z þC1Z1ðC0C1x2tana2þ 2tana1þ tana0ðÞC192C1ZC1Z þC1Z0C18C19This expression is similar to wnumdT/dz with wnumgiven byequation (14). There are three differences. The expressionfor the vertical velocity wnumuses a2rather than theweighted average of a0, a1and a2and uses C1Z rather thanthe average of C1Z and C1Z1. The last factor does not appearin the calculation of wnum. The first change is a detriment,the second and third are improvements.[68] In particular, for the coarse grid:C0@T@t¼@T@z0:009cmsC01C0C1;which is an improvement of 10% reducing the error to 90%of the expected total vertical velocity including thatbackground horizontal velocity as well as the shear.[69] The Shchepetkin-McWilliams advection scheme hasan additional term (C0d1) in the compution of a tracer fluxthrough a grid cell surface where C0 = C01/4 and d1= T(x) C02T(x C0 C1x)+T(x C0 2C1x)forU positive. This term onlyaffects the horizontal advection. Also, for the verticaltemperature advection a four-grid point calculation is used.Using these expressions instead of the simple second ordercenter difference scheme, equation (A2) for our chosen boxbecomesC0C1Z þC1Z02C1ZC1sC18C19@T@t¼UC1x916Tx;s þC1sðÞC20C0116Tx;s þ 2C1sðÞC034TxC0C1x;sðÞþ18Tx;sðÞþ18TxC0 2C1x;sðÞÞC21where we note that the code uses the no-heat flux bottomboundary condition. To evaluate this expression it isnecessary to consider the structure of the grid furtherupwards and farther to the right. The geometric parametersare given on Figure A1 and the result isC0@T@t¼@T@z0:007cmsC01C0C1[70] This scheme is an improvement over the directcalculation of w and reduces the error to about 70% ofthe estimated vertical velocity. Thus the vertical advectionof the velocity field (which should be strongly downward inthis region) is reduced to 70%.[71] Acknowledgments. This work was supported by a NationalScience and Engineering Research Council of Canada Strategic Grant toS. Allen and by National Science Foundation Grants OCE 96-18239 to J.Klinck and OCE 96-18186 to B. Hickey. Kurt Nielson was instrumental inthe design and building of the laboratory apparatus. The authors thank ananonymous reviewer for his/her detailed and helpful comments on aprevious version of the manuscript.ReferencesAllen, S. E., Topographically generated, subinertial flows within a finitelength canyon, J. Phys. Oceanogr., 26, 1608–1632, 1996.ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODEL 3 - 15Allen, S. E., M. S. Dinniman, J. M. Klinck, and B. M. Hickey, Dynamics ofadvection-driven upwelling over a submarine canyon, Eos Trans. AGU,80, 221, 1999.Allen, S. E., C. Vindeirinho, R. E. Thomson, M. G. G. Foreman, and D. L.Mackas, Physical and Biological processes over a submarine canyonduring an upwelling event, Can. J. Fish. Aquat. Sci., 58, 671–684,2001.Batchelor, G. K., An Introduction to Fluid Dynamics, 615 pp., CambridgeUniv. Press, New York, 1967.Boyer, D. L., X. Zhang, and N. Pe´renne, Laboratory observations of rotat-ing, stratified flow in the vicinity of a submarine canyon, Dyn. Atmos.Oceans, 31, 47–72, 2000.Butman, B., Long term current measurements in Lydonia and Oceanogra-pher Canyons, Eos Trans. AGU, 64, 1050–1051, 1983.Gent, P. R., J. Willibrand, T. J. McDougal, and J. C. McWilliams, Para-meterizing eddy-induced tracer transports in ocean circulation models,J. Phys. Oceanogr., 25, 463–474, 1995.Haidvogel, D. B., and A. Beckmann, Numerical models of the coastalocean, in The Sea, edited by K. H. Brink, and A. R. Robinson, vol.10, chap. 17, pp. 457–482, John Wiley, New York, 1998.Haney, R. L., On the pressure gradient force over steep topography in sigmacoordinate ocean models, J. Phys. Oceanogr., 21, 610–619, 1991.Hedstro¨m, K. S., User’s manual for an S-coordinate primitive equationocean circulation model (SCRUM) version 3.0, Contrib. 97-10, Inst. ofMar. and Coastal Sci., Rutgers Univ., New Brunswick, N. J., 1997.Hickey, B. M., The response of a steep-sided narrow canyon to strong windforcing, J. Phys. Oceanogr., 27, 697–726, 1997.Kinsella, E. D., A. E. Hay, and W. W. Denner, Wind and topographic effectson the Labrador current at Carson canyon, J. Geophys. Res., 92, 10,853–10,869, 1987.Kliem, N., and J. D. Pietrzak, On the pressure gradient error in sigmaco-ordinate ocean models: Comparison with a laboratory experiment,J. Geophys. Res., 104, 29,781–29,799, 1999.Klinck, J. M., Circulation near submarine canyons: A modeling study,J. Geophys. Res., 101, 1211–1223, 1996.Osler, G., Density gradients, Sci. Am., 213, 70–76, 1965.Perenne, N., D. B. Haidvogel, and D. L. Boyer, Laboratory-numericalmodel comparisons of flow over a coastal canyon, J. Atmos. OceanicTechnol., 18, 235–255, 2001.Petruncio, E. T., Observations and modeling of the internal tide in a sub-marine canyon, Ph.D. thesis, 176 pp., Nav. Postgrad. School, Monterey,Calif., 1996.Shchepetkin, A. F., and J. C. McWilliams, Quasi-Monotone advectionschemes based on explicit locally adaptive dissipation, Mon. WeatherRev., 126, 1541–1580, 1998.She, J., and J. M. Klinck, Flow near submarine canyons driven by constantwinds, J. Geophys. Res., 105, 28,671–28,694, 2000.Song, Y. H., and D. B. Haidvogel, A semi-implicit ocean circulation modelusing a generalized topography-following coordinate system, J. Comput.Phys., 115, 228–244, 1994.Weast, R. C., editor, Handbook of Physics and Chemistry, Chem. Rubber.Publ. Co., 59th edition, 1979.C0C0C0C0C0C0C0C0C0C0C0C0C0C0C0C0C0C0C0C0C0C0S. E. Allen, D. D. Gorby, and A. J. Hewett, Oceanography, Department ofEarth and Ocean Sciences, University of British Columbia, 6270 UniversityBoulevard, Vancouver, B. C., Canada, V6T 1Z4. ( S. Dinniman and J. M. Klinck, Crittenton Hall, Center for CoastalPhysical Oceanography, Old Dominion University, Norfolk, VA 23529,USA.B. M. Hickey, School of Oceanography, University of Washington, Box357940, Seattle, WA 98195, USA.3 - 16 ALLEN ET AL.: COMPARISON OF A NUMERICAL AND A PHYSICAL MODEL


Citation Scheme:


Usage Statistics

Country Views Downloads
United States 20 9
China 10 15
Canada 8 0
Japan 3 0
Germany 3 2
Russia 2 0
Ukraine 2 0
France 1 0
Sweden 1 0
City Views Downloads
Unknown 7 2
San Francisco 6 0
Beijing 6 0
Ashburn 5 2
Kingston 5 0
Mountain View 5 2
Tokyo 3 0
Shenzhen 2 14
North York 2 0
Redmond 1 0
Atlanta 1 0
Novosibirsk 1 0
Seattle 1 1

{[{ mDataHeader[type] }]} {[{ month[type] }]} {[{ tData[type] }]}
Download Stats



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