Transport of shoreline-released tracer from the surfzone across the shelf can be affected by a variety of physical processes from wind-driven to submesoscale, with implications for shoreline contaminant dilution and larval dispersion. Here, a high-resolution wave–current coupled model that resolves the surfzone and receives realistic oceanic and atmospheric forcing is used to simulate dye representing shoreline-released untreated wastewater in the San Diego–Tijuana region. Surfzone and shelf alongshore dye transports are primarily driven by obliquely incident wave breaking and alongshore pressure gradients, respectively. At the midshelf to outer-shelf (MS–OS) boundary (25-m depth), defined as a mean streamline, along-boundary density gradients are persistent, dye is surface enhanced and time and alongshelf patchy. Using baroclinic and along-boundary perturbation dye transports, two cross-shore dye exchange velocities are estimated and related to physical processes. Barotropic and baroclinic tides cannot explain the modeled cross-shore transport. The baroclinic exchange velocity is consistent with the wind-driven Ekman transport. The perturbation exchange velocity is elevated for alongshore dye and cross-shore velocity length scales < 1 km (within the submesoscale) and stronger alongshore density gradient ∂ρ/∂y variability, indicating that alongfront geostrophic flows induce offshore transport. This elevated ∂ρ/∂y is linked to convergent northward surface along-shelf currents (likely due to regional bathymetry), suggesting deformation frontogenesis. Both surfzone and shelf processes influence offshore transport of shoreline-released tracers with key parameters of surfzone and shelf alongcoast currents and alongshelf winds.
Shoreline-released wastewater or runoff enters the surfzone (SZ, region of depth-limited wave breaking) delivering pathogens and contaminants to coastal regions, threatening the health and sustainability of coastal ecosystems (e.g., Ahn et al. 2005; Steele et al. 2018). For example, 35 million gallons per day (mgd) of untreated wastewater is released 10 km south of the Pacific U.S.–Mexico border at Pt. Bandera, Mexico (e.g., Orozco-Borbón et al. 2006). Shoreline tracer dilution occurs through exchange across the surfzone, inner shelf, and farther offshore. Similarly, the coastal connectivity of intertidal invertebrates (e.g., Becker et al. 2007; Shanks et al. 2010) also requires cross-shelf exchange. A fraction of shoreline-released runoff or small river input is transported alongshore in the SZ, dependent on the flow rate and waves (Wong et al. 2013; Rodriguez et al. 2018). Obliquely incident wave breaking vertically mixes tracers (e.g., Feddersen 2012) and drives SZ alongshore currents (Longuet-Higgins 1970; Feddersen et al. 1998), transporting SZ tracer alongcoast up to 10 km (Grant et al. 2005; Feddersen et al. 2016). SZ tracers are cross-shore transported (exchanged) to the inner shelf due to transient rip currents (Hally-Rosendahl et al. 2014, 2015; Suanda and Feddersen 2015) or bathymetric inhomogeneities (e.g., Castelle and Coco 2013; Brown et al. 2015). Inner-shelf tracer dilution occurs through transport to the midshelf (MS) and outer shelf (OS). Thus, understanding cross-shore transport pathways across these coastal regions is important for understanding the fate of shoreline-released tracer.
A variety of processes across a range of time scales can induce cross-shelf tracer transport, and relevant processes are reviewed here. On subtidal (>33 h) time scales, winds, waves, bathymetric variability, or regional-scale (10–100 km) alongshelf pressure gradients (APG) can drive offshore tracer transport. For an alongshelf uniform shelf, alongshelf wind driven Ekman layers induce cross-shore tracer transport for both stratified and unstratified conditions (e.g., Austin and Lentz 2002; Lentz and Fewings 2012); however, as the depth becomes shallower than the Ekman depth, Ekman transport shuts down. Cross-shore winds also induce transport albeit weaker than alongshore winds (e.g., Fewings et al. 2008; Horwitz and Lentz 2016; Wu et al. 2018). Bathymetric variations such as headlands, capes, and shoals can steer mean flow at a variety of length scales (e.g., Gan and Allen 2002; Castelao and Barth 2006; Radermacher et al. 2017), inducing cross-shore transport. A barotropic APG-driven flow induces a bottom Ekman transport with a compensating interior cross-shore geostrophic current (Lentz 2008; Marchesiello and Estrade 2010), yielding cross-shore exchange.
Semidiurnal and diurnal barotropic (BT, surface) and baroclinic (BC, internal) tides also play an important role in cross-shelf exchange (e.g., Pineda 1994; Walter et al. 2014). Cross-shore dye transport is induced by nonzero covariance between the cross-shore velocity and dye over a tidal cycle, analogous to tidal pumping inducing estuary–ocean exchange (e.g., Lerczak et al. 2006; Geyer and MacCready 2014). On the inner shelf, strong BT tides can induce residual flow via tidal rectification (e.g., Ganju et al. 2011). In the Southern California Bight (SCB), semidiurnal (e.g., Lerczak et al. 2003; Buijsman et al. 2012; Kumar et al. 2016; Sinnett et al. 2018) and diurnal BC tides are ubiquitous (Lerczak et al. 2001; Nam and Send 2013; Kumar et al. 2015) even though the diurnal frequency is subcritical at SCB latitudes. The semidiurnal BC tide can drive an onshore cold-water transport through tidal pumping (e.g., Walter and Phelan 2016), and can cross-shore export heat from the nearshore (Sinnett and Feddersen 2019). BC tides enhance model simulated horizontal and vertical tracer dispersion in 30–50-m water depth (Suanda et al. 2018), and enhance mixing, reducing stratification (Suanda et al. 2017) and inducing residual cross-shelf flow. The diurnal internal tide is primarily responsible for the diurnal offshore and onshore advection of a SCB beach-released dye (Grimes et al. 2020).
Submesoscale flows have O(1) km length scales, are often seen as fronts and filaments, and have O(1) or greater Rossby number Ro = ζ/f where ζ is the vertical relative vorticity and f is the local planetary vorticity (e.g., McWilliams 2016). Fronts and filaments can be a primary driver for offshore tracer transport hundreds of kilometers from shore (e.g., Nagai et al. 2015), and may be important in cross-shelf transport within 5 km of shore (e.g., Romero et al. 2016). Submesoscale flow variability is ubiquitous in coastal drifter observations (e.g., Ohlmann et al. 2017) and high-resolution coastal models (e.g., Dauhajre et al. 2017). Generation mechanisms of open ocean (deep water) submesoscale variability include mixed layer instability (e.g., Boccaletti et al. 2007), turbulent thermal wind balance (e.g., McWilliams et al. 2015), and deformation flow induced frontogenesis (e.g., Hoskins 1982). However, frontogenesis processes in shallow (relative to horizontal scales) and frictional coastal waters are less clear and understudied. In coastal environments, bathymetric variations (e.g., Pringle 2002) or wind forcing (e.g., Tilburg and Garvine 2003) can drive surface deformation flow that may induce frontogenesis.
As the effects of the Pt. Bandera shoreline-released untreated wastewater are unknown, we seek to understand the offshore transport pathways in the San Diego Bight using a realistic model that includes the surfzone. Regional studies include examination of San Diego Bay (SDB) tidal outflow (e.g., Chadwick and Largier 1999), Pt. Loma upwelling (Roughan et al. 2005), episodic small river plumes (e.g., Warrick et al. 2007), and the evolution of shoreline-released dye (Hally-Rosendahl et al. 2014, 2015; Grimes et al. 2020). HF radar estimated regional currents reveal rich variability on length scales of ≈10 km with O(1) Rossby number (Kim 2010). However, the HF radar cannot resolve length scales < 5 km. Dye observations within 1 km of shore in this region reveal spatial variability from 0.01 to 1 km (Hally-Rosendahl et al. 2015; Grimes et al. 2020). Thus, within 10 km of shore, significant submesoscale flow variability is likely present.
Dye and drifters in realistic shelf models have been extensively used to study Lagrangian transport and dispersion processes (e.g., Uchiyama et al. 2014; Romero et al. 2013; Giddings et al. 2014; Romero et al. 2016; Dauhajre and McWilliams 2019), with varying grid resolutions. Shoreline larval dispersion patterns for the San Diego Bight were simulated with a realistic model at 600-m grid resolution (Rasmussen et al. 2009), but without a surfzone. A Lagrangian transport model forced with HF radar currents was used to estimate Pt. Bandera exposure kernels (Kim et al. 2009). However, currents were unknown within 1–2 km of shore, surface trapped tracking, and objective mapping velocity smoothing lead to significant exposure uncertainty (Kim et al. 2010). Coupled wave and circulation models allow both the surfzone and shelf to be resolved (Kumar et al. 2012), which is crucial as shoreline-released tracer often is surfzone alongshore transported (Feddersen et al. 2016; Grimes et al. 2020).
Here, shoreline-released tracer is simulated from the SZ to the outer shelf in the San Diego Bight using a high-resolution wave–current coupled model with realistic forcing. The specific goal is to elucidate the principal mechanism(s) responsible for offshore tracer transport in the mid- to outer-shelf region using relatively simple metrics. The mechanisms considered include Ekman, barotropic and baroclinic tidal, and submesoscale flows. Analysis focuses on a 3-month (midsummer to fall) period characterized by weak to moderate winds, prevalent southerly incident surface gravity waves, strong alongshore pressure gradients, and active internal waves. The model configuration and methods are given in section 2. The spatiotemporal tracer variability and transport are presented in section 3. Alongshore transport mechanisms are examined in section 4. Mechanisms for cross-shore transport across a MS–OS boundary are diagnosed in section 5. The relationships between key parameters governing dye transport and the role of other processes are discussed in section 6. Section 7 provides a summary.
a. Model configuration
Surfzone and shelf circulation is simulated using the Coupled Ocean–Atmosphere–Wave–Sediment–Transport (COAWST) model system (Warner et al. 2010; Kumar et al. 2012). Three one-way nested parent models (from LV1 to LV2 to LV3) and one high-resolution child run (LV4, see Fig. 1a for the four grids) are run from midsummer to early winter 2015 spanning a seasonal transition from summertime mostly southerly incident surface waves to wintertime mostly northerly incident waves. The LV1–LV3 runs employ the Regional Ocean Modeling System (ROMS; Shchepetkin and McWilliams 2005), a three-dimensional, hydrostatic model using a stretched terrain-following vertical coordinate (Shchepetkin and McWilliams 2005). The LV4 run includes surface gravity waves by coupling ROMS with the Simulating Waves Nearshore model (SWAN; Booij et al. 1999), a phase-averaged wave model that solves the wave action balance equation. NAM surface flux fields (wind stress, heat, and precipitation) are used for all grids. Over the 5-month simulation period, the 21 days of NAM data gaps are filled with the Coupled Ocean–Atmosphere Mesoscale System (COAMPS) fields. The vertical viscosity and diffusivity are estimated using a k–ϵ scheme (Umlauf and Burchard 2003). A logarithmic bottom drag scheme is used with a bottom roughness z0 = 0.1 cm, following Kumar et al. (2015). The horizontal eddy viscosity and diffusivity are constant at 0.5 m2 s−1. For the LV4 model, SWAN and ROMS are two-way coupled at 10-min intervals allowing current effects on waves and wave effects on currents through surfzone wave breaking and vortex force.
1) Parent run grids and setup
The model grids of the three parent runs downscale from 2-km horizontal resolution for the SCB (LV1 with 253 × 390 horizontal grid cells), to 600-m resolution resolving the southern SCB (LV2 with 266 × 398 grid cells), and to 200-m resolution for the greater San Diego shelf region (LV3 with 251 × 413 grid cells) (see Fig. 1a). All three domains have 40 terrain-following vertical levels with enhanced resolution near the surface and bottom. Grid bathymetry is derived from the 3-arc-second NOAA/NGDC coastal relief dataset.
The outermost LV1 domain inherits the boundary and initial conditions from the California State Estimate (CASE) solution, an implementation of the z-level, primitive equation MIT general circulation model (MITgcm; Marshall et al. 1997). CASE assimilates a variety of remote and in situ observations, including satellite altimetry data, satellite measured sea surface temperature, temperature and salinity profiles from Argo and Spray glider, expendable bathythermograph (XBT) temperature transects, autonomous pinniped bathythermograph (APB) temperature profiles, and shipboard CTD profiles (Zaba et al. 2018). Daily averaged CASE solutions are linearly interpolated from z- to σ-level coordinates, and then horizontally interpolated onto the LV1 model grid and open boundaries. CASE does not include tides. Barotropic tidal elevation and velocities of 10 tidal constituents (M2, S2, N2, K2, O1, P1, Q1, K1, M4, and M6) are prescribed on the LV1 open boundaries with the amplitudes and phases from the ADCIRC tidal database (Westerink et al. 1993), allowing generation and propagation of internal waves within the model domain (e.g., Kumar et al. 2015; Suanda et al. 2017; Kumar et al. 2019).
The LV1 solutions provide initial and boundary conditions for LV2. Subsequently, the LV2 solutions are used for LV3, and the LV3 solutions are used for LV4. Chapman and Flather radiation boundary conditions are used for the sea level and the barotropic (depth independent) velocity (Flather 1976; Chapman 1985). For the baroclinic (depth dependent) flow and tracers, the Orlanski radiation condition is used together with nudging to constrain the interior solution to the parent (Marchesiello et al. 2001). In LV1–LV3 domains, the nudging time scale for outgoing baroclinic flow and tracers along open boundaries is 365 day−1, and the nudging time scale for the incoming baroclinic flow and tracers is 6 h−1. All solutions are saved at 1-h intervals. The LV1 model was initialized at 1200 UTC 1 July 2015, and the LV2 and LV3 models were initialized at 1200 UTC 4 July and 1200 UTC 7 July, respectively, allowing 3–4 days of spinup in each parent grid.
2) Lv4 run grid and setup
The LV4 grid (with 486 × 1142 grid cells, area 15 × 36 km2) spans the outer shelf to surfzone in the southern San Diego Bight (Fig. 1b), which includes the SDB and the headland Pt. Loma to the north. Southward of the SDB entrance, the shoreline first curves and then straightens, passing the Tijuana River Estuary (TJRE), the U.S.–Mexico border and Punta Bandera (PB) within Mexico. South of the curvature, the bathymetry is largely alongshore uniform, except for a broad shoal seaward of TJRE that extends offshore (Fig. 1b). The LV4 grid cross-shore resolution increases from 110 m at the western open boundary to 8 m along the coastline, and alongshore resolution varies from 110 m at the southern and northern open boundaries to 8 m near the TJRE mouth. The stretched vertical domain has 15 levels and the NOAA 1/3-arc-second coastal digital elevation is used for bathymetry.
The LV4 SWAN model has 25 frequencies between 0.04 and 0.29 Hz and 42 directional bands spanning from 145° to 355° (wave direction in nautical convention), covering all potential incidence angles. The shoreline normal direction south of 32.6°N is approximately 265°. CDIP wave model frequency-directional wave spectra are used for open boundary conditions (O’Reilly et al. 2016). The wave-breaking parameter γ = 0.5 is used following Kumar et al. (2015). Wind-wave generation is also included. Note that because SWAN is a wave-averaged model, the LV4 simulation has bathymetric rip currents but does not have transient rip currents which require a wave-resolving model (Feddersen 2014).
The LV4 ROMS component receives freshwater inputs from PB, TJRE, and the Sweetwater River within SDB (see Fig. 1b). The parent grids do not receive freshwater input. At PB, freshwater, representing untreated wastewater, is released onto the beach at a constant discharge rate Qr = 1.53 m3 s−1 (35 mgd). TJRE freshwater discharge is given by in situ measurements at International Boundary and Water Commission (IBWC) gauging station and discharge primarily occurs during rainfall events (Fig. 2a). Additional coastal runoff emanating from the Sweetwater River within SDB is also incorporated in LV4. The Sweetwater River discharge rate is approximately estimated by multiplying the observed flow rate at a nearby river (San Diego River) by the ratio of the drainage area (Archfield and Vogel 2010). As PB and Sweetwater River inflow temperatures are unknown, a 30-day low-pass-filtered in situ Tijuana River (Oneonta Slough) temperature measurement is applied for all three sources to remove weekly and higher-frequency variations that could be variable among sites. Passive tracer (dye) of constant concentration D = 1, representing untreated wastewater, is added to the PB freshwater discharge. The off-diagonal radiation stress tensor term Sxy (Longuet-Higgins 1970; Feddersen et al. 1998) is estimated as
where θw denotes the mean wave angle relative to the shoreline normal, Hs denotes the significant wave height, ρ0 = 1025 kg m−3 denotes a reference density, and cg and c denote the group and phase speeds, respectively. A positive Sxy corresponds to southerly incident waves, which drives northward alongshore surfzone currents. Subtidal filtering is performed with the PL64 filter (Limeburner et al. 1985) with 33-h cutoff.
The LV4 coupled model was initialized on 12 July 2015, allowing 5 days of parent LV3 model spinup, and integrated to 25 December 2015. Relevant time series are shown in Fig. 2. Modeled barotropic tidal amplitudes (Fig. 2b) and phases of M2, S2 and K1 compare well with in situ measurements (not shown here). NAM winds (Fig. 2c) are consistent with nearby buoy winds (not shown) and are frequently southward directed with low (|Uw| < 5 m s−1) to moderate 5–8 m s−1 speeds. Wave and alongshelf velocities are given at a 30-m-depth central location denoted SB (Fig. 1b). At SB, Hs varies between 0.5 and 1.5 m (Fig. 2d). The LV4 simulation has multiple periods of southerly incident waves (i.e., Sxy > 0, yellow shading in Fig. 2e). These wave statistics compare well with local buoy observations, consistent with previous Southern California results (O’Reilly et al. 2016). At SB, the subtidal depth-averaged alongshore flow VSB ranges from −0.1 to 0.3 m s−1 and is mostly northward (positive) (Fig. 2f). The modeled dye concentration and dye transport are analyzed during the analysis period from 22 July to 18 October 2015 (dashed line in Fig. 2). This allows a 10-day period of LV4 dye spinup (22 days from LV1 initialization) at which point LV4 domain averaged enstrophy has equilibrated. After 18 October, occurrences of northerly incident incoming waves end northward transport of PB dye.
b. Analysis methods
To facilitate the cross- and alongshore dye transport analysis, two regions are defined (Fig. 1b). The first region is the nearshore (NS) spanning from the shoreline with a cross-shore width of . The NS reaches 10-m water depth and includes both the SZ and the shallow portion of the inner shelf. The NS southern boundary is located 4.5 km north of PB allowing dye to adjust upon release and the alongshore length is . The northern extent is set to avoid the rapidly curving isobath farther north. The second region is the MS–OS boundary (Fig. 1b) with mean depth of 25 m (Fig. 1c) and alongshore length of . The MS–OS boundary curves offshore such that southern and northern ends are approximately 3.1 and 6.3 km from the shoreline, respectively (white dashed lines in Fig. 1b). As offshore transport can be bathymetrically induced, the MS–OS boundary was chosen such that the depth-averaged and time-averaged (over the analysis period) flow across all parts of the MS–OS boundary is zero (MS–OS boundary is a mean streamline). Note, at any time step, the depth and along-boundary averaged cross-boundary flow is not necessarily zero. Analysis of cross-shore dye transport will be performed on the MS–OS boundary. To facilitate alongshore analysis, a MS–OS transition zone is defined centered at the MS–OS boundary with a width of (dashed magenta, Fig. 1b).
Hereafter, the cross-shore (x, positive onshore) and alongshore (y, positive northward) directions are locally defined as the normal to and parallel to the MS–OS boundary. Time averages (over the analysis period) are denoted with ⟨⋅⟩. Standard deviations are represented by std(⋅). Because dye is positive definite and does not have a Gaussian distribution, the time-averaged dye ⟨D⟩ is based on a logarithmic average such that (e.g., Hally-Rosendahl et al. 2014). Temporal mean plus (minus) standard deviation of dye are defined as . Within NS the dye is volume averaged ().
Averaging operators are defined for analysis along the MS–OS boundary. For a variable ψ, an along-MS–OS boundary and depth-average is defined as
where the along-boundary (y) integral is over the MS–OS boundary length and the average MS–OS boundary total depth (the average still water depth is ). The along-boundary integral is only calculated over the wet portion of the MS–OS boundary, thus is an area average through the MS–OS y–z surface. In addition, an along-boundary average is defined as
where, because the MS–OS boundary has variable depth, L(z) is the along-boundary length of the wet portion of the MS–OS boundary at each z level. If the z level is always wet, then .
The cross-MS–OS-boundary velocity u is decomposed into three components,
representing the depth- and along-boundary averaged transport , the along-boundary averaged baroclinic (vertically varying) velocity , and the along-boundary perturbation velocity u′. These are defined as
where uL is the model Lagrangian horizontal velocity (Eulerian plus Stokes drift derived from the wave model) and n is the normal to the MS–OS boundary. The time average of the first component by definition as the MS–OS boundary is a streamline of the depth- and time-averaged velocity. Overall, fluctuates between −0.01 and 0.02 m s−1 (Fig. 2g), attributed to bathymetric induced transport due to the significant inverse correlation (r = −0.4) between and VSB, consistent with mass conservation. The along-boundary averaged baroclinic velocity [Eq. (5b)] has by definition . The third component u′ [Eq. (5c)] is associated with short-scale (relative to ) variability and by definition . This MS–OS boundary decomposition is also applied to the dye and density,
At the MS–OS boundary, the total cross-shore dye transport Qx is
where again the along-boundary integral is only calculated over the wet portion. The Qx is composed of three components defined as
Note that the terms , , , , , and are all zero by definition, which is also confirmed numerically. For the baroclinic and along-boundary fluctuating dye transports, corresponding MS–OS boundary exchange velocities are defined as
The is ascribed to baroclinic flow developed over a regional (<10 km) scale, while is attributed to shorter-scale processes with significant u′ and D′. Both and are only estimated when dye is present at the MS–OS boundary , which removes 7% of the data. Note, the exchange velocity is similarly defined but just equals .
Within NS and the MS–OS transition zone, an along-region averaged along-region dye transport velocity is similarly defined as
where (r) represents NS or MS/OS, x1 and x2 represent the offshore and onshore region locations [e.g., for NS, , and for MS–OS, x1 and x2 represent the dashed lines bounding 500 m on either side of the MS–OS boundary in Fig. 1b], and the Lagrangian alongshore velocity is υL. We note that dye mass balances within control volumes close, confirming the dye and dye transport estimates. To average over the tidal and higher-frequency variability, the quantities , Qx, , , and are subtidally filtered.
3. Results: Spatiotemporal dye variability
a. Example of an offshore dye transport event
Upon Pt. Bandera shoreline release at concentration D = 1, dye is advected at a range of spatiotemporal scales spanning surfzone and shelf. Four snapshots of dye D, density anomaly σt = ρ − 1000 kg m−3, and currents covering 18 h are presented to show the spatiotemporal evolution of an offshore dye transport event (Fig. 3). Event winds and wind stresses were weak (Fig. 2c) with average wind onshore at 2 m s−1 and a ±1 m s−1 diurnal rotating sea breeze. At the event start, southerly incident waves (yellow shading in Fig. 2e) have just arrived driving northward NS dye transport (Figs. 3a1–4). Dye concentration is low (D < 10−4) on the northern shelf and within SDB (Figs. 3a1–4). An alongshore density gradient is present and the shelf surface currents are primarily northward (Figs. 3b1–4).
At the first time step (1400 UTC 7 August, Fig. 3a1), surface D is advected northward. A high concentration (D > 10−4) and meandering dye patch extends seaward and obliquely crosses the MS–OS boundary with width (where D > 10−4) of 3.8 km and peak D = 2.7 × 10−3. At the dye offshore leading edge, the offshore current is 0.1 m s−1. The surface density perturbation has additional small-scale variability (Fig. 3b1). Along the cross-shore transect (green line in Fig. 3a1), dye is concentrated within upper 5-m layer above the thermocline T = 19°C (Fig. 3c1). The cross-shore current is relatively weak (magnitude < 0.05 m s−1, Fig. 3d1). Six hours later (2000 UTC 7 August, Fig. 3a2), the surface dye patch has been advected farther northward and elongated into a filament with MS–OS boundary width of 580 m at high peak concentration D = 2.4 × 10−4. At the dye offshore leading edge, the offshore current is 0.07 m s−1. The elongation is due to a positive ∂u/∂x associated with a negative ∂υ/∂y (not shown), which also enhances the alongshore density gradient at the MS–OS boundary (Fig. 3b2). At the cross-shore transect, the surface layer onshore current strengthens (reaching 0.2 m s−1, Fig. 3d2), leading to near-surface onshore dye transport. As a result, south of the filament, surface D is completely contained within the NS (Fig. 3a2). Subsurface, the T = 19°C isotherm and dye layer deepen to z = −14 m within 1 km from shore mostly within the NS (see the vertical dye contour D = 10−4 in Fig. 3c2).
Twelve hours later (0200 UTC 8 August, Fig. 3a3), the dye filament has been advected farther north, orienting more cross-shore. Seaward of the MS–OS boundary, the filament is advected offshore by the surface currents at 0.2 m s−1 (Fig. 3b3). At the MS–OS boundary, the dye filament has a width of 1.3 km with max D = 1.4 × 10−4. South of the filament, surface dye is almost completely contained within the NS (Fig. 3a3). At the cross-shore transect, the onshore (offshore) current within the surface (bottom) layer is well developed (reaching 0.1 m s−1 at surface, Fig. 3d3). Both the 20°C isotherm (Fig. 3c3) and 23.6 kg m−3 isopycnal (Fig. 3d3) are tilted down toward shore, and the near-bed dye is advected offshore to z = −18 m, but remains within 1.75 km of shore. Eighteen hours later (0800 UTC 8 August), the dye filament decays (Figs. 3a4,b4). Onshore of the MS–OS boundary and south of the strong alongshore density gradient, surface dye is advected offshore extending the D = 10−4 contour onto the MS–OS boundary, as surface and bottom cross-shore currents reverse (Figs. 3c4,d4). Subsurface, the isotherms and isopycnals flatten and the near bottom dye is advected back onshore. As wind forcing was weak, weak Ekman transport is expected. The short alongshore scales of this offshore dye transport event suggest that submesoscale flows are responsible for the offshore dye transport.
b. Dye and density statistics
The spatial variability of temperature T, salinity S, density anomaly σt, and dye statistics is investigated next (Fig. 4). For depth < 25 m, the time-mean temperature ⟨T⟩ is elevated on the northern shelf and SDB by ≈0.5°C relative to the southern shelf near PB (Fig. 4a1). This alongshore ⟨T⟩ signal is in the parent LV3 results (not shown here), indicating it is not due to the temperature of the LV4 freshwater discharge. The largest temperature standard deviation std(T) (>1.4°C) occurs in the nearshore (h < 10 m, Fig. 4b1), due to stronger diurnal warming and cooling in shallow water, without a clear alongshore gradient. Relatively low ⟨S⟩ and high std(S) (>0.3) occur near the three freshwater sources (Figs. 4a2,b2). The shelf ⟨S⟩ also has a north–south gradient, fresher to the north (Fig. 4a2), resulting in an alongshelf gradient of ⟨σt⟩ (Fig. 4a3) with northern ⟨σt⟩ lower by 0.2 kg m−3. The std(σt) (>0.3 kg m−3) is elevated near freshwater sources and in the nearshore (h < 10 m), implying combined contributions from S and T. Along the MS–OS boundary, the ⟨σt⟩ gradient is 5.8 × 10−6 kg m−4 but the std(σt) is largely alongshore uniform. The surface dye statistics ⟨D⟩− and ⟨D⟩+ (defined in section 2b, Figs. 4a4,b4) highlight the upper and lower ranges of typical dye, varying from 10−2 to 10−5. Both show northward dye transport and dilution away from PB, resembling a diffusive plume, with higher concentrations onshore, implying net offshore dye transport (Fig. 4b4). Onshore of the MS–OS boundary, ⟨D⟩+ is about 100 times that of ⟨D⟩−.
Large values of the surface density horizontal gradient magnitude |∇Hρ| and rms(ζ/f) > 1 indicate fronts and filaments and high Rossby numbers. Elevated values of both quantities are seen in the NS (Figs. 4a5,b5), attributed to surface wave breaking, bathymetric irregularities, freshwater inputs and surface heat fluxes. Within the NS rms(ζ/f) > 3 suggesting an active submesoscale (Fig. 4b5). Both statistics are also elevated near the SDB mouth. Onshore of the MS–OS boundary, the TJRE shoal has elevated rms(|∇Hρ|) and rms(ζ/f), while for the rest of the domain, the distribution of both is largely alongshore uniform. At the MS–OS boundary, rms(|∇Hρ|) = 1.3(±0.2) × 10−4 kg m−4 and rms(ζ/f) = 0.8(±0.1), consistent with other high-resolution coastal simulations (e.g., Dauhajre et al. 2017; Dauhajre and McWilliams 2019), and indicates that submesoscale dynamics may be important in this region. Although the PB (1.5 m3 s−1) and TJRE (Fig. 2a) freshwater input rates are typically small, salinity gradients contribute an average ≈30% to the rms(|∇Hρ|), and are elevated during time of TJRE freshwater input (Fig. 2a).
The temporal and alongshore variability of MS–OS boundary surface density perturbation ρ′ [Eq. (6b)] and dye D are examined in Fig. 5. The ρ′ has a persistent negative south to north gradient (Fig. 5a), consistent with Fig. 4a3. Significant variability is present with times of alongshore uniform (e.g., 6–10 October, Fig. 5a) or step function like ρ′ (9 August, near 32.56°N). MS–OS boundary surface dye also has very strong variability (Fig. 5b). The northern end almost always has D > 10−5, consistent with the mean dye fields (Figs. 4a4,b4). The temporal D variability is largely diurnal to subtidal. The along-boundary D variability is often patchy with a few km or less length scales, consistent with this offshore dye ejection event (Fig. 3). Dye can be present along much of the boundary (26 July), only in a limited region (9 August), or nearly absent (D < 10−5, 10 September).
The logarithmic scale in Fig. 5b obscures the strong along-boundary dye gradients, and the dye patchiness suggests that the MS–OS boundary dye has relatively short alongshore length scales. A MS–OS boundary alongshore surface dye length scale is defined as
which is evaluated at each time step when MS–OS boundary averaged dye is >10−6, and also subtidally filtered (Fig. 5c). The varies from 0.3 to 4 km, confirming that the MS–OS boundary D is patchy, and is also consistent with the example event 0.58–3.8-km D length scales (Fig. 3). The subtidal varies from 0.4 to 3 km. Length scales for density ρ and cross-boundary velocity u are similarly estimated along the MS–OS boundary (Lu and Lρ) which also vary between 0.4 and 4 km, jointly suggesting the importance of submesoscale dynamics (e.g., McWilliams 2016).
The vertical structure of the time mean and std of , , and , and the square root of time alongshore mean of squared ρ′, D′, and u′ are examined at the MS–OS boundary (Fig. 6). During the summer and early fall, the MS–OS boundary is strongly stratified (Fig. 6a1) with time mean (or N ≈ 0.02 s−1) that is two orders of magnitude or larger than the along-boundary stratification rms(|∇Hρ|) (Fig. 4a5). Associated with the vertical stratification, time mean baroclinic dye is surface intensified and decays exponentially downward with 6.3-m decay scale (Fig. 6a2). The temporal std of and is elevated at the surface and reaches a minima at 10-m water depth (Figs. 6a1,a2). The time-mean cross-boundary baroclinic velocity has magnitude usually <0.01 m s−1 and a three-layer profile with offshore directed flow at surface and bottom, and onshore flow in between (Fig. 6a3). The temporal std of the subtidally averaged varies between 0.01 and 0.03 m s−1, is elevated at surface (Fig. 6a3), and has a variability minimum near z = −10 m. The profile statistics are consistent with a combined Ekman and mode-1 baroclinic process. Note, the subtidal filter removes baroclinic tidal velocities, which will be discussed in section 5a.
In addition, significant temporal and alongshore variability is evident in ρ′, D′, and u′ (Figs. 6b1–b3) through their time and alongshore standard deviation (i.e., ). The is elevated in the upper 10 m near 0.07 kg m−3 (Fig. 6b1), about 10% of the top to bottom difference (Fig. 6a1), but is comparable to the temporal std, suggesting significant alongshore density fluctuations, consistent with Fig. 5a. Similarly, is near-surface elevated, decaying with a ≈10 m vertical scale. The is comparable to the std, consistent with Fig. 5b, and the O(1) km inferred (Fig. 5c). The increases from 0.02 m s−1 near bed to 0.04 m s−1 near surface (Fig. 6b3), slightly larger than the temporal std of (Fig. 6a3). This decomposition makes clear that MS–OS boundary offshore dye transport can occur due to both alongshore-uniform baroclinic processes and alongshore variable processes.
c. Dye and dye transport temporal evolution
During the analysis period, seven individual south swell events (i.e., Sxy > 0, Fig. 2e) occur (yellow shading in Figs. 7 and 2e) each with distinct subtidal NS-averaged dye peaks mostly >10−3 (Fig. 7a). In between south swell events, drops by a factor of 10–100. The MS–OS boundary mean [Eq. (6a)] is significantly weaker (≈10%) than (Fig. 7a). The peaks are qualitatively lagged relative to peaks, indicating a transport pathway. The depth std of largely follows and is usually slightly above (Fig. 7a), as dye is surface intensified. During peaks in offshore , the alongshore-depth std of D′ usually is just larger than both and indicating along-boundary dye patchiness. The 10–12 September time period (cyan rectangle in Fig. 7a) is notable because is strongly elevated (>10−3) for an extended duration yet the MS–OS boundary , , and are all weak (mostly <10−6), and will be discussed in more detail in section 6c.
The depth- and along-boundary mean cross-MS–OS-boundary dye transport [Eq. (8a)] switches sign frequently (Fig. 7b), has a time mean of −0.05 dye m3 s−1 (representing 19% of total transport) and std of 0.12 dye m3 s−1, and is attributed to bathymetrically driven flows. The baroclinic dye transport [Eq. (8b)] is the largest of the three dye transports and is primarily offshore (i.e., negative, Fig. 7d) with time mean of −0.15 dye m3 s−1 (57% of total transport) and std of 0.33 dye m3 s−1 over a number of quasi intermittent events. The along-boundary perturbation transport [Eq. (8c)] is also intermittent and fluctuates with a std of 0.17 dye m3 s−1 relatively large relative to the time mean of −0.06 dye m3 s−1 (24% of total transport). Overall is the second-largest term although occasionally it is bigger than (e.g., 10 August, Fig. 7b). Nearly all of the net (time averaged) offshore occurs prior to 1 September, and for this period time mean and are similar. For the post 1 September time period, is negligible but is the same. The baroclinic and alongshore perturbation cross-shore dye exchange velocities [Eqs. (9)] are largely negative (seaward directed, Fig. 7c), varying between 0.02 and −0.07 m s−1. The time mean and both are −0.01 m s−1 with slightly larger than . Thus, for the time period when and could be calculated (when ), both baroclinic processes and alongshore-perturbation processes have similar transport potential. Next, we explore the mechanisms driving alongshore and cross-MS–OS-boundary dye transports.
4. Alongshore dye transport mechanisms
The 500-m-wide NS region (Fig. 1b) typically has a 100-m-wide surfzone for the incident wave heights (Fig. 2d). Surfzone alongshore currents are driven by the breaking of obliquely incident waves (e.g., Longuet-Higgins 1970; Feddersen et al. 1998). The nearshore subtidal alongshore dye transport velocity [Eq. (10)] varies between −0.1 and 0.1 m s−1 (Fig. 8) corresponding to 9–17 km day−1. Consistent with surfzone-dominated transport, the subtidal is largely positive (northward directed) during southerly wave events (yellow shading in Fig. 8a) and is highly correlated with Sxy/ρ0 with squared correlation of r2 = 0.63 (p < 0.05), best fit slope of ≈1, and near-zero intercept. This best fit slope is factor of 2 consistent with a simple surfzone alongshore wave forcing and linear bottom friction balance assuming a 100-m surfzone width and linear drag coefficient of 3 × 10−3 m s−1 (e.g., Lentz et al. 1999). In contrast, was only weakly related (r2 = 0.14, p > 0.05) to the alongshore wind stress. This demonstrates the primary role of obliquely incident surface gravity waves in driving alongshore NS dye transport over long (tens of kilometers) distances consistent with previous observations and modeling (Grant et al. 2005; Hally-Rosendahl et al. 2015; Hally-Rosendahl and Feddersen 2016; Feddersen et al. 2016). Other mechanisms such as wind driven currents in the outer NS, tidal currents, and shear dispersion can play a secondary role.
The 1-km-wide MS–OS transition zone has a subtidal alongshore dye transport velocity [Eq. (10)] that is significantly stronger than within the nearshore (Fig. 8a). The subtidal is strongly correlated to the subtidal depth-averaged SB alongshelf velocity VSB (Fig. 2f) with r2 = 0.92 (p < 0.05) and at 75% the magnitude of VSB (Fig. 8a), attributed to a reduced current (i.e., current shear) onshore of SB. Note, alongshore transport is larger than cross-shore transport as is larger than the three Uex. The analysis period is characterized by weak to moderate alongshore wind forcing (Fig. 2c). Subtidal alongshore wind stress is uncorrelated (r2 = 0.05, p > 0.05) with VSB and has magnitude 4 times too weak (using a linear friction of 3 × 10−4 m s−1; Lentz and Winant 1986) to explain VSB, suggesting that other dynamics are driving the alongshelf current. Previous SCB studies (e.g., Hickey et al. 2003) have shown that during fall the barotropic APG is a significant driver of alongshelf flow even in 15 m depth (Lentz and Winant 1986). The subtidal barotropic APG is estimated from north to south in 15-m depth within the LV4 domain (S1 and S2 in Fig. 1b). The resulting barotropic APG largely varies between ±2 × 10−6 m s−2, is mostly northward directed as the alongshelf flow and largely varies on fortnightly time scales (Fig. 8b). The barotropic APG is reasonably correlated with VSB (r2 = 0.49, p < 0.05) and has the correct magnitude for a frictionally balanced flow (using 15-m depth and linear friction of 3 × 10−4 m s−1). Thus, the regional alongshore current which drives midshelf alongshore dye transport is primarily driven by the barotropic APG.
5. Cross-shore dye transport mechanisms at the MS–OS boundary
As baroclinic [Eq. (8b)] and the along-boundary perturbation [Eq. (8c)] were the largest transport terms, here, we separately examine potential mechanisms for and , including BT and BC tides, wind-driven Ekman transport, barotropic alongshore pressure gradients, and submesoscale flows. Recall that can be due to alongshore uniform but baroclinic flows such as Ekman transport or BC tides, while can be induced by alongshore variable flows over the 15-km-long MS–OS boundary. To elucidate the hydrodynamic processes, we compare the dye exchange velocity components [Eq. (9a)] and [Eq. (9b)] with relatively simple derived hydrodynamic velocities for each process.
a. Barotropic and baroclinic tidal cross-shore transport
One potential mechanism for the cross-shore dye transport across the MS–OS boundary could be barotropic (BT) and baroclinic (BC) tides. In the example offshore dye transport event (Fig. 3), the diurnal BC tide advects surface dye onshore and offshore and vertically, similar to the observed evolution of shoreline-released dye near 32.59°N during the analysis period (Grimes et al. 2020). If BT and BC tidal exchange mechanisms (see section 1) were a primary driver of the cross-shore exchange, then the subtidal BC exchange velocity would be less than and correlated with the slowing varying tidal velocity amplitude. Here, in both the semidiurnal (SD) and diurnal (DU) bands, alongshore mean (and std) MS–OS tidal surface velocity amplitudes are estimated for BT , and BC tides (see the appendix). Both BC and BT tidal variability are largely along-boundary coherent (Figs. A1b,c), thus the tidal velocity amplitudes are compared to the MS–OS boundary baroclinic exchange velocity. However, as some BT and BC tide alongshore phase variation is present, tidal velocity amplitudes are also compared to the alongshore perturbation exchange velocity.
The relationship of barotropic tidal velocity amplitudes and and are examined first (Figs. 9a,b). The alongshore averaged fluctuates fortnightly between 0.01 and 0.02 m s−1 with very small std along the boundary, indicating that the single BT velocity amplitude is representative. However, is much weaker than and uncorrelated (r2 = 0.005, p > 0.05, see Table 1) with (Figs. 9a,b). The alongshore mean is generally <0.01 m s−1 and poorly correlated (r2 = 0.02, p < 0.05) with . Similarly, is uncorrelated with , and not significantly correlated with (r2 = 0.18, p > 0.05), albeit is half the magnitude of . Thus, BT tides cannot be the main driver of the or cross-MS–OS-boundary dye transports. For BC tides, the surface DU BC tidal amplitude is generally larger than the surface SD amplitude (Fig. 9c), similar to previous observations (e.g., Kim et al. 2011; Johnston and Rudnick 2015) and modeling (e.g., Kumar et al. 2015) in the SCB. The is smaller than and uncorrelated (r2 = 0.03, p > 0.05) with (Figs. 9a,c). The alongshore averaged is comparable in magnitude (reaches 0.15 m s−1) to (Figs. 9a,c), but is uncorrelated (r2 = 0.007, p > 0.05). Similarly, the alongshore perturbation is uncorrelated with and (see Table 1). Thus, BC tides also cannot be the main driver of the or cross-MS–OS-boundary dye transports.
b. Wind-driven Ekman cross-shore transport
Wind-driven Ekman transport is a potential mechanism for offshore dye transport. Here, wind-driven Ekman surface velocity at the MS–OS boundary is estimated following Ekman (1905) and compared with the dye exchange velocity . Conceptually, Ekman transport is considered an alongshore uniform and baroclinic process. Thus, if the offshore dye transport were due to Ekman transport, surface Ekman velocities should be of similar magnitude and correlated to .
For a steady, alongshore-uniform shelf with no alongshelf pressure gradients, shelf currents have an Ekman component (balancing friction) and alongshelf geostrophic component (balancing the cross-shelf pressure gradient) such that the depth-averaged cross-shore velocity is zero (Ekman 1905). For constant depth h, steady conditions, a no slip seafloor boundary condition, no stratification, and a constant eddy viscosity Aν, Ekman’s analytic solution is used to estimate the (cross- and along-mean MS–OS boundary) velocity [Uek(z), Vek(z)] (following Estrade et al. 2008):
where , (τx, τy) are the (subtidal) wind stress components (estimated from SB), m = (1 + i)(f/2Aν)1/2, and Vg is the alongshelf geostrophic component obtained from the condition that the depth integrated cross-shore transport is zero. A constant Aν = 2.5 × 10−3 m2 s−1 is used based on the depth- and time-averaged modeled eddy viscosity at SB, a value consistent with Suanda et al. (2017). Neglected transient effects in the Ekman velocity (12) are weak as subtidal winds vary on synoptic time scales. The Ekman depth (2Aν/f)1/2 = 8 m is consistent with the z ≈ 10 m location of minimum and variability (Figs. 6a2,a3). A vertically averaged (over the upper 8 m) Ekman velocity is defined as for comparison to . With the predominant southward wind (Fig. 2c), surface is primarily negative (Fig. 9d). The correlation is maximum at zero time lag between and . A least squares fit between and yields a near-one slope of 0.97 with significant r2 = 0.40 (p < 0.05), indicating that Ekman transport is a principal driver for the baroclinic cross-MS–OS-boundary transport .
The time mean has a three-layer profile (Fig. 6a3); however, Uek(z) [Eq. (12)] has a two-layer profile. This difference together with the moderate r2 between and may result from the assumptions built into [Eq. (12)] or from other processes (e.g., advection). We note that the barotropic APG is mostly northward directed (Fig. 8b) which would induce near-bed offshore Ekman transport and onshore return flow in the remaining water column (e.g., Lentz 2008). Thus, APG driven near-surface flow is of the opposite sign to and cannot drive near-surface offshore tracer transport here.
c. Submesoscale-flow induced cross-shore transport
The simulation snapshots (Fig. 3) show offshore propagating cross-shore elongated dye structures with a width of 0.6–3.8 km. In general, the MS–OS-boundary dye is patchy with subtidal dye alongshelf length scales varying from 0.5 to 3 km (Figs. 5b,c), scales consistent with the coastal submesoscale (e.g., Dauhajre et al. 2017). Alongshelf surface density gradients are also consistently present at a variety of scales (Fig. 5a). The along-boundary perturbation dye transport is a significant component of the total dye transport (Fig. 7b), particularly prior to 1 September. Here, we examine the role of submesoscale flows in driving along-boundary perturbation via relationships between [Eq. (9b)], the surface alongshelf dye length scale LD [Eq. (11)], and the MS–OS boundary rms surface alongshelf density gradient rms(∂ρ/∂y)(MS,OS).
The subtidal exchange velocity is more negative for smaller (Fig. 10) with binned-mean squared correlation of r2 = 0.32 (p < 0.05). On average is 2 times larger for than for . The subtidal cross-boundary velocity length scale varies between 0.4 and 2 km and is more negative for smaller (not shown). Overall, the stronger offshore linked to smaller and suggests that the perturbation offshore dye transport is largely due to relatively short (<2 km) length scales associated with submesoscale flows.
An along-boundary density gradient could induce a cross-MS–OS-boundary (i.e., alongfront) flow, and thus , through geostrophic adjustment. As the density is time dependent (Fig. 5a), such cross-MS–OS-boundary flow must adjust and would be elevated and time-lagged for stronger rms(∂ρ/∂y)(MS,OS). The time-lagged squared correlation r2 between rms(∂ρ/∂y)(MS,OS) (Fig. 11a) has a maximum r2 = 0.19 (p > 0.05) when rms(∂ρ/∂y)(MS,OS) leads by 0.22 inertial periods (5 h). After adjusting for this time lag, is consistently more negative for larger rms(∂ρ/∂y)(MS,OS) (Fig. 11b), with binned mean r2 = 0.51 (p < 0.05), particularly for rms(∂ρ/∂y)(MS,OS) > 0.5 × 10−4 kg m−4. Using a scaled thermal wind relationship and assuming depth-uniform along-boundary density gradients (e.g., McWilliams 2016), a scaling for the rms geostrophic cross-boundary velocity rms(Ug) can be written as
where d is the average vertical scale of the density gradient. Here, the rms density gradients are largely coherent over the upper 10 m of the water column (not shown) setting an average vertical scale of d = 5 m in the scaling [Eq. (13)]. With this, a subtidal rms(∂ρ/∂y)(MS,OS) = 0.5 × 10−4 kg m−4 yields a rms(Ug)(MS,OS) = 0.03 m s−1, a factor 3 times larger than , consistent with this process driving exchange. Note, will be smaller than (Ug)(MS,OS) because frontal velocities are both onshore and offshore and recirculate dye. This indicates that when strong density gradients associated with O(1) km length scales are oriented alongshelf, cross-shore flows are induced with a time lag that transport dye offshore at these submesoscale length scales (Fig. 10). A model that does not adequately resolve these submesoscale flows will underestimate offshore tracer transport.
a. Potential mechanism for rms alongshore density gradient generation
Several mechanisms have been posited for the generation of open-ocean submesoscale variability such as mixed layer instability (e.g., Boccaletti et al. 2007), turbulent thermal wind balance (e.g., McWilliams et al. 2015), and deformation flow induced frontogenesis (e.g., Hoskins 1982). The presence of a shoreline boundary and shallow depths adds additional complexity. Here, we examine the role of an alongshore surface deformation flow in enhancing the existing mean alongshore density gradient potentially leading to the enhanced rms(∂ρ/∂y)(MS,OS). Detailed diagnosis of the submesoscale variability and dynamics will occur elsewhere.
Density and alongshelf current statistics are time-averaged and cross-shore averaged within the MS–OS transition zone. The time-averaged and cross-shore averaged surface density has a quasi-linear time-mean alongshelf density gradient (Fig. 12a), consistent with the mean density distribution (Fig. 4a3). The average alongshelf density gradient is weaker by a factor 5–30 times than rms(∂ρ/∂y)(MS,OS) (Fig. 11b). During the analysis period, the subtidal VSB is mostly northward (Fig. 2f). The surface alongshelf velocity is first cross-shore averaged within MS–OS transition zone and then time-averaged when VSB > 0 yielding , which has a consistent negative alongshelf gradient (Fig. 12b), indicating alongshelf convergence, likely due to the regional bathymetry. A linear best fit yields a , 2–5 times larger than the normal strain rate of mesoscale eddies (e.g., Chaigneau et al. 2008). If this northward convergent flow was responsible for strengthening rms(∂ρ/∂y)(MS,OS), then a stronger rms(∂ρ/∂y)(MS,OS) with increasing VSB would be expected. Indeed, the subtidal rms(∂ρ/∂y)(MS,OS) is consistently enhanced for stronger northward VSB (Fig. 13) with r2 = 0.49 (p < 0.05) and binned mean r2 = 0.88 (p < 0.05). This fit does not change if only times when could be calculated (dashed line in Fig. 7a) are used in the fit. This indicates a linkage between the two processes and supports the concept that an alongshelf convergent flow is promoting generation of rms alongshelf density gradients.
b. Role of other mechanisms in offshore dye transport
Although wind-driven Ekman transport and alongfront submesoscale flows were diagnosed as principal drivers of MS–OS boundary offshore tracer transport, the region is complex and other mechanisms may play a role. The stratification and circulation on subtidal time scales can be modified by BT (e.g., Ganju et al. 2011) and BC (Suanda et al. 2017) tides. Shelf BC tides enhanced 3D and 2D horizontal dispersion relative to simulations without BC tides (Suanda et al. 2018). Although DU BC tides did not principally drive subtidal offshore dye transport through a tidal exchange mechanism, BC tides could be similarly enhancing offshore transport relative to a no BC tide simulation (not performed). The TJRE shoal and large-scale bathymetry (coastline curvature, SDB entrance, and Pt. Loma) could also have a secondary effect on the MS–OS offshore dye transport. For example, enhanced rms(ζ/f) (well above 1) and |∇Hρ| variability at both the SDB entrance and the TJRE shoal (Figs. 4a5,b5) suggests strong local vorticity and buoyancy gradient generation. The surfzone also has strongly elevated vorticity and |∇Hρ| variability relative to the MS–OS transition zone (Fig. 4). Although bathymetric rip currents are present here, transient rip current (TRC) forcing is not included in this model. However, TRCs also drive further enhanced vorticity (Johnson and Pattiaratchi 2006; Suanda and Feddersen 2015) and buoyancy gradient (Kumar and Feddersen 2017a,b) variability on alongshore length scales of 50–500 m (Hally-Rosendahl et al. 2014, 2015; Kumar and Feddersen 2017b). As the NS and regions farther offshore are material transport linked (as evidenced by dye), offshore transport of NS or TJRE shoal vorticity and |∇Hρ| may seed submesoscale variability on the inner shelf and farther offshore, although of course neither is a conserved passive tracer. Although relatively weak, the freshwater sources at TJRE, Pt. Bandera, and within SDB could be an additional density gradient source.
Last, we examine the role of horizontal mixing in potentially inhibiting submesoscale variability generation. A horizontal eddy viscosity of νh = 0.5 m2 s−1, 1-km length scale, and 0.1 m s−1 velocity scale give a horizontal Reynolds number of 200, indicating a priori weak horizontal mixing. Surface momentum balance terms calculated at the MS–OS boundary have rms horizontal mixing at ≈17 times and ≈23 times weaker than the vertical mixing and nonlinear advective terms, respectively. Thus, the horizontal mixing has a minor effect on submesoscale variability. A similar conclusion was drawn for turbulent-thermal wind filament frontogenesis from simulations with a similar grid resolution (12.5 m versus a mean of 20–30 m within the MS–OS transition zone) using both νh = 0.5 m2 s−1 and νh = 0 m2 s−1 (McWilliams et al. 2015).
c. Relative roles of transport terms and key parameters with implications for coastal water quality
Cross-shore transport has direct implications on shoreline water quality by diluting shoreline concentrations. For nonzero cross-MS–OS-boundary dye transport, PB shoreline-released dye must be present at the MS–OS boundary, requiring first northward NS transport driven by waves incident from the south (Sxy > 0, Fig. 2e). The baroclinic cross-MS–OS-boundary transport is linked to the wind-driven Ekman transport. Although relatively small, near-surface Ekman velocities were directed offshore for extended durations (Figs. 6a3 and 9d). Together with the surface enhanced dye, this resulted in significant MS–OS-boundary dye transport with exchange velocity magnitude similar to the Ekman velocity. The along-boundary perturbation transport [Eq. (8c)] is linked to increased boundary rms(∂ρ/∂y)(MS,OS) through alongfront geostrophic flow. Although the submesoscale velocities u′ (Fig. 6b3) are as large as the velocities, because the u′ flows fluctuate onshore and offshore, significant recirculation is likely present and the resulting MS–OS-boundary is substantially weaker than the rms velocities [Eq. (13) and Fig. 11b]. For the regional geography south of Pt. Loma, the rms(∂ρ/∂y)(MS,OS) is enhanced during northward convergent mean flow, primarily APG driven (Fig. 8), allowing diagnosis of rms(∂ρ/∂y)(MS,OS) from VSB. Thus, for a shoreline-released tracer, the key parameters cross-MS–OS-boundary transport are Sxy, alongshelf winds, and VSB. Other headland bounded regions may be similar.
During the analysis period, Sxy/ρ0 is mostly positive (Fig. 2e) and VSB is largely northward with only a few short times of southward flow (Fig. 2f). The winds are mostly upwelling favorable, with only a few cases of sustained downwelling favorable conditions (10–12 September and 4–6 October, Fig. 2c). The 4–6 October downwelling winds occur during negative to weak positive Sxy, and thus relatively little dye is present in both the NS and at the MS–OS boundary (Fig. 7a). The 10–12 September downwelling wind conditions are interesting because Sxy was strongly positive leading to strongly elevated (Fig. 7a) and shelf flow was also strongly northward (VSB > 0.2 m s−1, Fig. 2f), but had the weakest MS–OS boundary dye, mostly <10−6 (Fig. 7a).
To explain the lack of MS–OS boundary dye, we examine this 10–12 September event with two model snapshots spanning 57 h (Fig. 14). At the event start (0100 UTC 10 September), shoreline-released dye streams north remaining mostly within the 500-m-wide NS, before the flow separates near 32.63°N, bends around Pt. Loma, and continues north (Fig. 14a). Toward the end of the event (1000 UTC 12 September), dye is still concentrated in the 500-m-wide NS, but there are ≈1 km alongshore scale plume structures extending 1–3 km offshore, and, farther north, the dye plume wraps around Pt. Loma (Fig. 14b). In both cases, the D = 10−4 contour is always at >2 km from the MS–OS boundary. For the entire event duration, MS–OS boundary never had D ≥ 10−5 (Fig. 5b, cyan rectangle). During this event, winds were only moderately northward (Fig. 2c), yet the MS–OS-boundary surface was onshore at roughly 0.05 m s−1 for many days. This is likely due to onshore surface Ekman transport but may also be enhanced by onshore flow balancing offshore bottom Ekman transport from the northward APG-driven flow. In conclusion, both SZ and shelf processes influence shoreline concentrations for shoreline-released tracers. Simultaneously occurring northward surfzone currents, downwelling winds, and northward forcing APG, (e.g., 10–12 September) is the worst-case scenario for regional beach water quality, as offshore transport from the NS is weak, and the shoreline is not diluted.
Here, we investigate the processes transporting a shoreline-released dye representing untreated wastewater in the San Diego–Tijuana region that has a curving shoreline, an estuary, a bay, and a headland. A high-resolution wave–current coupled model is used resolving the surfzone and receiving realistic air–sea forcing, tides, waves, and offshore boundary conditions inherited from a larger-scale data assimilated model. Model dye is shoreline released 10 km south of the U.S.–Mexico border representing untreated wastewater and analyzed from midsummer to fall with largely southerly incident waves, mostly northward alongshelf flow, strong shelf stratification, and moderate wind forcing. Analysis focuses primarily on the tracer transport across the midshelf (MS) to outer-shelf (OS) boundary (≈25-m depth) chosen as a streamline of the time mean and depth-averaged velocity.
Within 500 m of the shoreline, alongshore tracer transport is primarily driven by obliquely incident wave breaking. At the MS–OS boundary, alongshore density gradients are persistent and dye is surface enhanced and time and alongshore patchy with length scales from 0.3 to 4 km. Significant vertical (baroclinic) and along-boundary density and velocity variability is present. The cross-MS–OS-boundary dye transport has significant baroclinic and along-boundary perturbation components from which baroclinic and along-boundary perturbation dye exchange velocities are estimated. Barotropic tides and semidiurnal baroclinic tides cannot explain these two exchange velocity components. The baroclinic exchange velocity is significantly correlated and has similar magnitude to a simply Ekman transport velocity, indicating Ekman transport is driving the baroclinic cross-MS–OS-boundary tracer transport. The perturbation exchange velocity is elevated for smaller (<1 km) alongshore dye length scales and stronger root-mean-square (rms) alongshore density gradients ∂ρ/∂y, indicating geostrophic along-frontal submesoscale flows induce the along-MS–OS-boundary perturbation transport. During periods of northward flows, the surface alongshore current is convergent with a relatively strong mean deformation rate. Stronger northward flows are linked to elevated rms(∂ρ/∂y)(MS,OS), potentially generated by deformation frontogenesis. A model that does not adequately resolve these submesoscale flows will underestimate offshore tracer transport. Both surfzone and shelf processes influence offshore transport for shoreline-released tracers, and the key parameters governing cross-MS–OS-boundary dye transport are the incident Sxy, the alongshelf winds, and the APG-driven alongshelf current. When the co-occurrence of these parameters strongly inhibits offshore transport, shoreline concentrations are not effectively diluted leading to poor water quality.
This work was supported by the National Science Foundation (OCE-1459389) as part of the Cross-Surfzone/Inner-shelf Dye Exchange (CSIDE) experiment. Additional funding is through the Environmental Protection Agency through the North American Development Bank, however, it does not necessarily reflect the policies, actions or positions of the U.S. EPA or NADB. NK was supported by ONR Award N00014-17-1-2890 and NSF Award OCE-1735460. ONR supported GG via NOPP Award N000141512598. The numerical simulations were performed on the comet cluster at the San Diego Super Computer Center through the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation (ACI-1548562). NOAA provided the NAM atmospheric forcing fields and the bathymetries. SIO Coastal Data Information Program provided wave forcing. Chris Edwards from UC Santa Cruz provided COAMPS data generated by NRL. We also appreciate extra support from the Tijuana River National Estuarine Research Reserve, the Southern California Coastal Ocean Observing System. CASE solutions are available online (http://ecco.ucsd.edu/case.html) Geno Pawlak, Derek Grimes, Andre Amador, and Bruce Cornuelle provided useful feedback to this work. We appreciate insight from Jack McSweeney.
Barotropic and Baroclinic Tide Surface Velocity Amplitude Estimation
Here, the MS–OS boundary barotropic and surface baroclinic tidal velocities estimation method is described. The cross-shore velocity is first decomposed into barotropic (BT, depth averaged) and baroclinic (BC) components along the boundary. These two components are then bandpass filtered to obtain the semidiurnal (SD, 16−1 to 10−1 cph) and diurnal (DU, 33−1 to 16−1 cph) band components. Here, the SD-band analysis is described. First, we calculate the bulk BT SD tidal velocity amplitude (envelope) (USD), assumed narrow banded with form
where ωSD is the SD frequency, is the SD velocity amplitude that varies on longer time scale (denoted as εt), and θSD is a phase that will only vary slightly over the MS–OS boundary. As the BT tidal velocity is narrow banded, the amplitude is estimated via Hilbert transform at each alongshore location. The alongshore mean and std are presented in Fig. 9b.
BC velocities vary with depth and have much shorter length scales than the BT tidal velocities, leading to additional analysis. At each location yi along the MS–OS boundary, the SD baroclinic cross-shore current uSD(yi, z, t) is decomposed into a vertical EOF such that
where and are the EOF amplitude and vertical structure at yi, and N = 15 is the total number of vertical levels. At all MS–OS boundary locations, the first (n = 1) EOF accounts for >78% of the SD-band variance (>92% of the DU-band variance). For both bands, the first vertical EOF [ and ] is consistent with a first-mode baroclinic motions with midwater column sign change and is nearly alongshore uniform on the MS–OS boundary (Fig. A1a). Thus, surface BC tidal cross-shore velocities can be reconstructed with a single EOF at all alongshore locations, i.e., for the SD-band surface velocity
The alongshore coherent variability of reconstructed surface SD cross-shore velocity [and also DU ] is further examined with a complex (Hilbert) EOF (CEOF) (e.g., Horel 1984; Merrifield and Guza 1990). A complex time series is generated according to
where is the Hilbert transform of uSD and . The variability of is then CEOF decomposed into
where M is the total number of MS–OS boundary grid points, is the complex eigenvector, and is the complex amplitude.
The alongshore first CEOF of the SD and DU band explains 56% and 95% of the alongshore variability, respectively. The magnitude of the SD is maximum near the center of the MS–OS boundary, and is reduced about 50% at the northern and southern ends (Fig. A1b). The DU is largely alongshore uniform (Fig. A1b). The SD along boundary phase is estimated as
where and are the imaginary and real operators, respectively. For the SD band, the phase varies quasi-linearly by π/2 (90°) along the 15-km boundary (Fig. A1c), suggesting southward propagation of the SD BC tide. For the DU band, the phase is near zero for the southern half of the boundary, and varies π/6 (30°) over the northern half (Fig. A1c), also indicating mostly southward propagation, consistent with the regional observations of southward DU-band propagation in 12–15-m depth (Grimes et al. 2020). Overall, the phase variations indicate that the SD and DU BC tides alongshelf scale is substantially larger than the 15-km length of the MS–OS boundary. The amplitude of the first CEOF reconstructed SD and DU surface velocities are estimated as for the BT tidal velocity at each y location resulting in surface and . The alongshore mean and std are shown in Fig. 9c and described in section 5a.
Denotes content that is immediately available upon publication as open access.