Science of the Total Environment 854 (2023) 158685 Contents lists available at ScienceDirect Science of the Total Environment journal homepage: www.elsevier.com/locate/scitotenv Spatial heterogeneity in sediment and carbon accretion rates within a seagrass meadow correlated with the hydrodynamic intensity ⁎ Jiarui Lei a,b, Rachel Schaefer a, , Phil Colarusso c, Alyssa Novak d, Juliet C. Simpson e, Pere Masqué f,g,h, Heidi Nepf a a Department of Civil and Environmental Engineering, Massachusetts Institute of Technology, 77 Massachusetts Avenue, Room 1-290, Cambridge, MA 02139, United States of America b Department of Civil and Environmental Engineering, National University of Singapore, Block E1A, #07-03, No. 1 Engineering Drive 2, 117576, Singapore c U.S. Environmental Protection Agency, 5 Post Office Square, Suite 100, Boston, MA 02109, United States of America d Department of Earth & Environment, Boston University, 685 Commonwealth Avenue, Room 130, Boston, MA 02215, United States of America e MIT Sea Grant College Program, Massachusetts Institute of Technology, 77 Massachusetts Avenue, Room NW98-151, Cambridge, MA 02139, United States of America f School of Science and Centre for Marine Ecosystems Research, Edith Cowan University, 270 Joondalup Drive, Joondalup 6027, Australia g International Atomic Energy Agency Environment Laboratories, 4 Quai Antoine 1er, B.P. 800, 98012 Monaco Cedex, Principality of Monaco h Department of Physics, Institute of Environmental Science and Technology, Universitat Autònoma de Barcelona, Edifici C, 08193, Bellaterra (Cerdanyola del Vallès), Barcelona, Spain H I G H L I G H T S G R A P H I C A L A B S T R A C T • Spatial gradients in carbon and wave in- tensity were compared across a Zostera marina (L.) meadow. • Higher near-bed turbulence was corre- lated with lower sediment and carbon ac- cretion rates. • Storm events may dominate resuspension and therefore carbon storage. • Wave dissipation by the seagrass meadow influenced the spatial variation in carbon. A R T I C L E I N F O A B S T R A C T Editor: Jan Vymazal The majority of the carbon stored in seagrass sediments originates outside the meadow, such that the carbon storage capacity within a meadow is strongly dependent on hydrodynamic conditions that favor deposition and retention of Keywords: fine organic matter within the meadow. By extension, if hydrodynamic conditions vary across a meadow, they may Zostera marina give rise to spatial gradients in carbon. This study considered whether the spatial gradients in sediment and carbon ac- Carbon accretion cretion rates correlated with the spatial variation in hydrodynamic intensity within a single meadow. Field measure- Turbulent kinetic energy Resuspension ments were conducted in three depth zones across a Zostera marina L. (eelgrass) meadow in Nahant Harbor, Massachusetts. Four sediment cores were collected in each zone, including one outside the meadow (control) and three within the meadow at increasing distances from the nearest meadow edge. Sedimentation and carbon accretion rates were estimated by combining the measurements of dry bulk density, organic carbon fraction (%OC), 210Pb, and 226 Ra. Tilt current meters measured wave velocities within each zone, which were used to estimate turbulent kinetic energy (TKE). Both sediment and carbon accretion rates exhibited spatial heterogeneity across the meadow, which were correlated with the spatial variation in near-bed TKE. Specifically, both accretion rates increased with decreasing ⁎ Corresponding author. E-mail addresses: jlei@nus.edu.sg (J. Lei), rbs@mit.edu (R. Schaefer), colarusso.phil@epa.gov (P. Colarusso), abnovak@bu.edu (A. Novak), simpsonj@mit.edu (J.C. Simpson), p.masque@ecu.edu.au (P. Masqué), hmnepf@mit.edu (H. Nepf). http://dx.doi.org/10.1016/j.scitotenv.2022.158685 Received 20 April 2022; Received in revised form 17 August 2022; Accepted 7 September 2022 Available online 13 September 2022 0048-9697/© 2022 The Authors. Published by Elsevier B.V. This is an open access article under the CC BY-NC license (http://creativecommons.org/licenses/by-nc/4.0/). J. Lei et al. Science of the Total Environment 854 (2023) 158685 TKE, which was consistent with diminished resuspension associated with lower TKE. A method is proposed for using spatial gradients in hydrodynamic intensity to improve the estimation of total meadow accretion rates. 1. Introduction 2. Methods Seagrass meadows are global hotspots for carbon storage (e.g., Kennedy 2.1. Meadow and sediment characteristics et al., 2010; Duarte et al., 2010; Fourqurean et al., 2012). However, carbon stocks vary significantly between different seagrass sites. For example, This study was conducted in a meadow of Zostera marina in Nahant Har- Lavery et al. (2013) reported an 18-fold range in carbon stock measured bor, Massachusetts (Fig. 1), which is a wave-dominated site. Based on isoto- across seventeen different seagrass habitats (260 to 4800 g C m−2). Simi- pic signature analysis, the carbon in this meadow has a significant larly, Serrano et al. (2014) observed a 14-fold difference in carbon stock allochthonous fraction (Fig. 8 in Novak et al., 2020), so that the accretion over a 10-fold difference in water depth, which they attributed to gradients rate of organic carbon was expected to be influenced, at least in part, by hy- in light availability. The variability in carbon stock is a major source of un- drodynamic conditions favorable for its retention. Meadow boundaries certainty in assessing global seagrass carbon stocks; therefore, there is a were delineated through a boat survey. In Fig. 1, the outer edge of the need to understand the factors driving this variability (Lavery et al., 2013; meadow is shown with pink line segments, which enclose 14.3 acres Serrano et al., 2014). (0.06 km2). A denser inner region (5.2 acres, 0.02 km2) is shown with More than 50 % of the carbon stored in seagrass sediments originates out- green line segments. The meadow was divided into three depth zones, side the meadow (known as allochthonous carbon), so that the carbon storage which were expected to experience different hydrodynamic conditions: capacity within a meadow is dependent on hydrodynamic conditions that Shallow (mean water depth = 2.2 m), Mid (mean depth = 3.2 m), and favor deposition and retention of fine organic matter (e.g., Gacia and Deep (mean depth = 6.2 m). Mean water depths were measured with Duarte, 2001; Gacia et al., 1999; Kennedy et al., 2010). Therefore, hydrody- diver depth gauges and corrected for tidal phase. In July 2020, within namic conditions should help to explain differences in measured seagrass car- each of the three depth zones, three cores were extracted within the bon stocks and burial rates (Granata et al., 2001; Gruber and Kemp, 2010). For meadow and a fourth reference core was collected from the adjacent bare example, Oreska et al. (2017) observed that sediment organic carbon in- bed, resulting in a total of 12 cores. A 50-cm core barrel with 7-cm diameter creased with distance from the meadow edge and attributed this to the atten- was manually driven to the point of refusal, which ranged from 24 cm to uation of current by meadow drag forces. Ricart et al. (2020) attributed an 30 cm. The core locations are shown in Fig. 1. increase in seagrass meadow carbon stocks between the lower and upper re- After the cores were extracted, they were capped underwater and kept gions of an estuary to a decrease in hydrodynamic intensity. Similarly, vertical during transport from the extraction sites to the beach. Depending Novak et al. (2020) observed higher carbon stocks at sites with lower wave on the core length, the cores were extruded into 22 to 25 sections, with 1- and current exposure. cm increments from 0 to 20 cm and 2-cm increments from 20 to 30 cm. The correlation between hydrodynamic intensity and carbon stock is The core subsamples were stored in Ziploc bags on ice in coolers and likely due to sediment resuspension. For example, Dahl et al. (2018) transported to the Novak Lab at Boston University in Boston, Massachu- showed that sediment resuspension within a meadow led to a loss of carbon setts. Samples were put into a drying oven at 60 °C for seven days or until from the surface sediment. Previous studies have linked resuspension they attained a constant weight. The dry bulk density of each slice was cal- within a meadow to near-bed turbulent kinetic energy (TKE). Resuspension culated by dividing the dry mass of the sample by the volume. Next, the sed- occurs when TKE increases above a critical magnitude, which depends on iment samples were divided in half with a sediment splitter. One half of the median sediment grain size (e.g., Tinoco and Coco, 2018; Tang et al., each sample was placed in a muffle furnace at 450 °C for 16 h to determine 2019). For a unidirectional current interacting with a submerged meadow, the percent organic matter content (%OM) using the loss-on-ignition spatial gradients in TKE have been correlated with spatial gradients of net method (Howard et al., 2014). Percent organic carbon (%OC) was calcu- deposition. Specifically, when a unidirectional current encounters a sub- lated assuming a linear dependence on %OM based on data from sites merged meadow, the hydrodynamic drag generated by the meadow across New England reported in Novak et al. (2020). %OM, %OC, and causes the within-meadow current velocity and TKE to decrease with bulk dry density were averaged over the top 15 cm of each core. These pa- distance from the leading edge (e.g., Chen et al., 2013; Lei and Nepf, rameters showed no trends with depth into the cores. The grain size distri- 2021). Zhang et al. (2020) observed that net deposition increased bution was measured by combining samples from the top 5 cm and sieving with distance from the meadow leading edge, correlated with the pro- from 0.5 mm to 0.0625 mm. The median grain size was calculated from the gressive decrease in TKE. grain size distribution. In addition to reducing current velocity, seagrass meadows also attenuate The second half of each sediment sample was shipped to Edith Cowan wave energy, which may diminish near-bed wave velocity (e.g., Fonseca and University in Joondalup, Australia for 210Pb analysis to determine sediment Cahalan, 1992; Mendez and Losada, 2004; Infantes et al., 2012; Lei and Nepf, accumulation rates (SAR) and mass accumulation rates (MAR). 210Pb was 2019). The interaction between the wave velocity and the seagrass sheath gen- determined by analyzing 210Po using alpha spectrometry after the addition erates turbulence near the bed (Zhang et al., 2018), which can trigger sedi- of 209Po as an internal tracer and digestion in acid media using an analytical ment resuspension (Tinoco and Coco, 2018; Tang et al., 2019). Combining microwave (Sanchez-Cabeza et al., 1998). Gamma spectrometry measure- these effects, Tang et al. (2019) predicted the meadow size needed to diminish ments determined the concentration of 226Ra in selected samples along wave energy enough to eliminate wave-induced resuspension and enhance each core. The excess 210Pb was defined as the difference between the particle retention. total 210Pb and 226Ra (i.e., supported 210Pb). Using the 210Pb concentration Building on the studies described above, the present study used field ob- profiles, an average MAR was calculated using the Constant Flux and Con- servations at a wave-dominated site to examine the connection between stant Sedimentation (CF: CS, applied piecewise) model (e.g., Krishnaswami sediment and carbon accretion rates, wave velocity, and TKE within a single et al., 1971; Arias-Ortiz et al., 2018). Cores that indicated significant mixing meadow. An understanding of how resuspension may impose spatial gradi- were not analyzed for MAR. The carbon accretion rate (CAR) was estimated ents on carbon accretion could improve the assessment of carbon credit for by combining MAR with %OC (see Table 1). seagrass meadows, as well as provide guidelines to optimize long-term car- On the same day the cores were extracted, the shoot density was mea- bon accretion potential. To this end, a method is proposed for using pre- sured near each core site with three randomly placed quadrats (25 cm x dicted gradients in hydrodynamic exposure to improve estimates of 25 cm). The leaf width, length, sheath length, number of leaves, and meadow-scale carbon accretion rates. shoot biomass were measured for five selected shoots within each quadrat. 2 J. Lei et al. Science of the Total Environment 854 (2023) 158685 Fig. 1. The field site in Nahant Harbor, Massachusetts, on the East Coast of the United States. Four cores were collected in each of three depth zones; three within the meadow at increasing distance from the nearest meadow edge, numbered 1 (diamonds), 2 (circles), 3 (stars), respectively, and one reference core outside of the meadow (“Ref,” blank markers). White, pink, and green denote the Shallow, Mid, and Deep depth zones, respectively. Green and magenta outlines show the regions of dense (shoot density > 200 m−2) and sparse (shoot density < 200 m−2) vegetation mapped using boat surveys in 2014. The bathymetry surveyed in 2014 is indicated by labeled blue contour lines. Leaf area index (LAI) was calculated as the product of the mean leaf area per the cylinder measure the degree and direction of tilt, which is converted shoot (sheath included) and shoot density (see Table 1). through factory calibration coefficients into velocity and heading. Data were recorded at 16 Hz for 60-s bursts every 5 min at all sites except Shal- 2.2. Velocity and turbulence low(Ref). Due to concerns about boat traffic, the instrument at Shallow (Ref) recorded at 16 Hz continuously for 24 h, at which point it was re- In August 2021, six tilt current meters (TCMs) designed and manufac- moved. The remaining TCMs were deployed for a minimum of two tured by Lowell Instruments, LLC (East Falmouth, Massachusetts) were de- weeks. All TCMs were fully submerged during their deployment. The ployed in Nahant Harbor at Shallow(Ref), Shallow1, Mid(Ref), Mid2, Deep total height of each TCM setup was 33 cm, which was less than the observed (Ref), and Deep3. The TCM deployment sites were chosen because the ex- meadow height at all sites, confirmed by measured shoot height and GoPro tracted sediment cores at these sites yielded mass accretion rates from the footage recorded at the Deep3 site. Although TCMs are primarily marketed 210 Pb analysis. Each TCM consisted of a cylindrical housing tethered by a to record unidirectional flow, a comparison with an acoustic Doppler velo- short cord to an anchor, which in this study was a 20-lb, 12-in by 12-in cimeter (Nortek Vectrino) in an experimental channel confirmed their abil- by 2-in square concrete paving stone. In the absence of water motion, the ity to record wave velocity (see Fig. A.1 in the Supporting Information). The long axis of the cylinder is vertical. The cylinder tilts in response to the average root-mean-squared velocity, URMS, was estimated over each 1-h in- wave and current velocities. An accelerometer and magnetometer within terval, from which wave velocity amplitude was defined as Uw = 21/2URMS. Table 1 Summary of measured meadow structure and sediment characteristics averaged over the top 15 cm of the sediment cores, consistent with the depth of cores that was used to assess MAR and SAR (Fig. 2). Standard deviations among quadrat samples (shoot density) and over sediment profiles (dry bulk density, %OC, SAR) are reported. Derived quantities (LAI, organic carbon, CAR) include propagated uncertainty. Distance to Mean water Dry bulk %OC Shoot LAI Organic carbon Median grain MAR SAR CAR edge depth density density size m m g/cm3 % m2 m2/m2 g C/cm3 mm g/cm2 yr mm/yr g C/m2 yr D1 11 6.3 1.54 ± 0.13 0.33 ± 0.08 91 ± 46 0.9 ± 0.5 0.0051 ± 0.0012 0.11 15 ± 5 D2 19 6.2 1.53 ± 0.14 0.35 ± 0.04 59 ± 9 0.69 ± 0.11 0.0054 ± 0.0007 0.11 16 ± 4 D3 29 6.0 1.45 ± 0.20 0.42 ± 0.07 85 ± 18 0.74 ± 0.16 0.0062 ± 0.0013 0.11 0.45 ± 0.10 3.1 ± 0.7 19 ± 5 DR 1 6.5 1.41 ± 0.19 0.29 ± 0.07 0.0040 ± 0.011 0.12 M1 5 3.0 1.21 ± 0.19 0.68 ± 0.25 261 ± 9 3.04 ± 0.11 0.008 ± 0.003 0.10 14 ± 6 M2 54 3.3 1.36 ± 0.17 0.47 ± 0.16 331 ± 76 7.9 ± 1.8 0.007 ± 0.002 0.10 0.21 ± 0.04 1.6 ± 0.3 10 ± 4 M3 72 3.3 1.39 ± 0.16 0.47 ± 0.11 261 ± 65 4.3 ± 1.1 0.0065 ± 0.0017 0.09 10 ± 3 MR 1 3.0 1.36 ± 0.14 0.21 ± 0.03 0.0028 ± 0.0005 0.10 S1 11 2.0 1.41 ± 0.18 0.21 ± 0.08 240 ± 48 3.7 ± 0.7 0.0030 ± 0.0011 0.10 0.12 ± 0.02 0.90 ± 0.18 2.5 ± 1.0 S2 19 2.3 1.41 ± 0.17 0.20 ± 0.04 251 ± 76 4.8 ± 1.4 0.0027 ± 0.0006 0.09 2.3 ± 0.6 S3 29 2.6 1.33 ± 0.39 0.23 ± 0.05 176 ± 55 1.2 ± 0.4 0.0031 ± 0.0011 0.10 2.8 ± 0.8 SR 8 1.7 1.33 ± 0.22 0.29 ± 0.08 0.0039 ± 0.0013 0.07 3 J. Lei et al. Science of the Total Environment 854 (2023) 158685 Velocity records confirmed that this is a wave-dominated site, with wave 3. Results and discussion velocity much larger than time-mean velocity. Tidal elevations at high and low tides were obtained from Nahant US Harbors data. 3.1. Meadow structure and sediment characteristics Because previous studies have correlated resuspension with near-bed turbulence, the measured near-bed wave velocity was used to estimate The measured meadow structure and sediment characteristics are sum- near-bed turbulence. Turbulence is generated by the interaction of the marized in Table 1. The shoot density was highest in the Mid zone and low- wave velocity with the bed and with the relatively stiff sheath of each est in the Deep zone. The measured leaf width averaged over all quadrats seagrass shoot (Zhang et al., 2018; Tang et al., 2019). The vegetation- was b = 0.47 ± 0.09 cm. The median grain size was comparable across generated turbulent kinetic energy per fluid mass (TKE) can be predicted the meadow, increasing by just 20 % between the Shallow and Deep from a model based on the conversion of wave kinetic energy to TKE by zones. The nearly uniform grain size across the meadow can be attributed canopy drag (Tanino and Nepf, 2008; Tang et al., 2019). The bed- to a common source of suspended material, derived from a nearby marsh. generated turbulence is linearly proportional to the bed-generated shear In the Deep zone, the total 210Pb concentrations decreased with depth in stress, τ w ¼ ρf w U w 2 =2, in which fw is a wave friction factor (Tang et al., the meadow cores (green, blue, and red symbols in Fig. 2(a)) from around 2019). The sum of these predicts the total near-bed TKE: 25 Bq/kg at the surface down to a constant of 14.9 ± 1.6 Bq/kg below 10 to 12 cm. Based on the Deep1, Deep2, and Deep3 cores, the mass accre- tion rate was MAR = 0.45 ± 0.01 g cm−2 yr−1, and the sedimentation ac-  23 ns b2 cretion rate was SAR = 3.1 ± 0.7 mm/yr (Fig. 2(a)). In the meadow Mid near  bed TKE ¼ δ CD U 2RMS þ Cb,w U 2RMS , (1) 2ð1  ϕÞ zone, MAR = 1.6 ± 0.3 g cm−2 yr−1 and SAR = 1.6 ± 0.3 mm/yr based vegetation  generated turbulence bed  generated turbulence on the Mid2 core (Fig. 2(b)). In the meadow Shallow zone, MAR = 0.90 ± 0.18 g cm−2 yr−1 and SAR = 0.90 ± 0.18 mm/yr based on the Shallow1 in which the scale factor δ = 1.2 ± 0.2 is an average from two previous core (Fig. 2(c)). studies (Tanino and Nepf, 2008; Tang et al., 2019). ns is the number of shoots per bed area. b is the sheath diameter (see Section 3.2 in Zhang 3.2. Velocity, turbulent kinetic energy, and sediment accumulation rate et al., 2018), which was assumed to be equal to the leaf width. Since the sheath has a nearly cylindrical geometry, the near-bed solid volume frac- Wind speeds observed during the deployment were typical of historical wind conditions (based on records from the weather station at the Boston tion ϕ = π4 ns b2 . CD (=1.2) is the sheath drag coefficient, based on cylinder Logan International Airport), so that the observed spatial trends in hydro- drag measurements in Keulegan and Carpenter (1958). The bare-bed scale fac- dynamic intensity were considered to be representative of the long-term tor Cb,w, defined by measurements over bare bed, is a function of bed rough- conditions that produced the recorded spatial trend in mass accretion. ness. We used a value measured by Tang et al. (2019), Cb,w = 0.02 ± 0.01 Power spectra of velocity records indicated dominant wave periods be- for d50 = 85 μm, which is close to the sediment grain size at Nahant tween 6 and 10 s, which was consistent with historical data from the Na- (d50 = 100 μm). tional Oceanic and Atmospheric Administration National Data Buoy The critical TKE threshold was estimated using the critical velocity for Center (NOAA NDBC) Station 44,013 buoy in Massachusetts Bay. For this bare beds, modified for meadows (Tang et al., 2019). First, Komar and wave period and local depth range, the near-bed wave velocity is sensitive Miller (1973) showed that for a median grain size d50 <0.05 cm, the critical to depth, specifically increasing as depth decreases. As a result, the magni- wave velocity for a bare bed, Uw,c,b, is given by tude of the near-bed URMS was correlated with the tidal phase (Fig. 3), with peak URMS values associated with low tide (vertical gray lines in Fig. 3), rffiffiffiffiffiffiffiffiffiffiffiffiffi which was particularly evident at the Shallow1 site (black curve in ρU w,c,b 2 2Aw,c,b ¼ 0:3 , (2) Fig. 3). Over the two-week deployment, the average RMS velocity was ðρs  ρÞgd 50 d 50 highest at the Shallow1 site (6.6 cm/s), less at the Mid2 site (5.3 cm/s), and smallest at the Deep3 site (4.7 cm/s). Within each depth zone, the in which ρs is the sediment density, g is the gravitational acceleration, RMS velocity measured just outside the meadow differed from the in- Aw,c,b = UwT/(2π) is the critical wave excursion, and T is the wave pe- meadow velocity by <10 % (data not shown). That is, the meadow did pffiffiffi not locally diminish the wave velocity within the meadow compared to riod. The critical RMS wave velocity is then, U RMS,c,b ¼ U w,c,b = 2:The critical velocity threshold within the meadow is reduced, relative to just outside of the meadow. This was consistent with the prediction of the bare bed, due to the contribution from vegetation-generated tur- wave velocity attenuation described in Lowe et al. (2005). Specifically, bulence. Specifically, the critical rms wave velocity in a meadow, based on the meadow density and wave period, the in-canopy wave veloc- URMS,c,v, is (Tang et al., 2019) ity attenuation factor α ≈ 1 for all depth zones. Large wave events on Days 2, 7, and 13 were correlated with stormy weather reported on Logan Inter- national Airport wind speed and precipitation data. U RMS,c,v 1 The near-bed TKE was estimated from measured URMS using Eq. (1). We ¼ rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 23 : (3) U RMS,c,b 2  2  considered the strong wave conditions observed during low tides on Days 2, 1 þ Cδb,w 2CπD 3 1 ϕ ϕ 7, and 13, referred to as Storms 1, 2, and 3. Consistent with the spatial var- iation in URMS in Fig. 3, the Shallow1 site experienced the highest TKE while the Deep3 site experienced the lowest TKE (Fig. 4). The critical TKE The critical RMS wave velocity can be used in Eq. (1) to estimate the threshold for sediment resuspension (based on Eqs. (1), (2), and (3)) are in- critical TKE. We emphasize that this prediction is for near-bed, wave- dicated with horizontal lines in Fig. 4(a). During Storm 3, the TKE at the generated turbulence. This is our focus because near-bed turbulence has Shallow1 and Mid2 sites significantly exceeded the thresholds, while the been shown to correlate with resuspension (Tang et al., 2019). The thresholds were barely crossed at the Shallow1 and Mid2 sites during vertically-averaged, wave-driven turbulence within a meadow is much Storms 1 and 2. Based on this, we infer that events with velocity magnitude smaller, because the turbulence decays sharply with distance from the similar to or greater than that of Storm 3 control the spatial distribution of bed (see Fig. 3(f) in Zhang et al., 2018). This is because bed-generated tur- sediment carbon within the meadow. With this in mind, in the next phases bulence is limited by a thin wave-boundary-layer, and because the seagrass of analysis we focus on Storm 3. leaves can easily move with the waves, reducing relative motion and thus The SAR (Table 1, Fig. 4(b)) varied inversely with near-bed TKE during reducing vegetation-generated turbulence above the sheath region (see dis- Storm 3, which can be explained by resuspension triggered by near-bed tur- cussions in Zhang et al., 2018 and Tang et al., 2019). bulence, i.e., resuspension diminished particle retention and thus sediment 4 J. Lei et al. Science of the Total Environment 854 (2023) 158685 Fig. 2. Total 210Pb and 226Ra versus depth from the surface of the meadow sediment cores for the (a) Deep zone, (b) Mid zone, and (c) Shallow zone. Only profiles without mixing from which SAR and MAR were determined are included. Colors correspond to the site distance from the nearest edge of the meadow. Green, blue, and red denote meadow sites at increasing distances from the nearest edge of the meadow, respectively. Fig. 3. (a) Water depths at the Shallow site at high and low tides, and (b) hourly-averaged URMS over two weeks. The black curve denotes the Shallow zone (Shallow1). The red curve denotes the Mid zone (Mid2). The blue curve denotes the Deep zone (Deep3). Vertical dotted and dashed lines denote times of low and high tide, respectively. accretion rate. This was consistent with previous laboratory studies, which observed between sedimentation rate and near-bed TKE was consistent found sediment resuspension increased as near-bed TKE increased in pure with the impact of resuspension (Fig. 4). wave and pure current conditions (Tinoco and Coco, 2018; Tang et al., 2019; Zhang et al., 2020). The spatial gradient of TKE results in a spatial gradient of exceedance of the critical TKE threshold, and therefore a spatial 3.3. Sediment organic carbon content and carbon accretion rates gradient of resuspension. Consistent with this, previous sediment trap mea- surements within the Nahant meadow indicated significantly higher resus- The measured CAR values in the Nahant meadow were typical of those pension both within and adjacent to the meadow in the Shallow zone observed in other meadows of the same seagrass species. Specifically. pre- compared to the other depth zones (data not shown). Assuming sediment vious studies report a range of CAR in Zostera marina L meadows from 3.1 supply was the same across the meadow, an increase in resuspension to 230 g C m−2 yr−1 (Novak et al., 2020; Poppe and Rybczyk, 2018; would lead to a lower net sedimentation rate, so that the correlation Prentice et al., 2020). The CAR in the Shallow zone is at the low end of 5 J. Lei et al. Science of the Total Environment 854 (2023) 158685 Fig. 4. (a) Turbulent kinetic energy (TKE) at Deep3, Mid2, and Shallow1 for Storms 1, 2, and 3 (Days 2, 7, and 13 events). The horizontal lines indicate the critical TKE thresholds for each meadow site. The thresholds varies between sites due to the different meadow densities, i.e. nb in Eq. (1). (b) Sediment accretion rate (SAR, black bars, left-hand axis) and near-bed turbulent kinetic energy (TKE, red bars, right-hand axis) at Deep3, Mid2, and Shallow1 during Storm 3. this range, while the CAR in the Deep zone is close to the average of re- directly influenced by a spatial variation in shoot density (e.g., see the dis- ported CAR in these studies. cussions in Lowe et al., 2005 and Luhar et al., 2008). Therefore, at current- The mean sediment organic carbon fractions (%OC) were statistically dominated sites with a patchy meadow structure, one would expect spatial significantly different between the Shallow and Mid Zones but not between variation in shoot density to be associated with spatial variation in velocity, the Mid and Deep or Shallow and Deep zones (see Fig. B.1 in the Supporting and thus, potentially, with spatial variation in sediment carbon. Information). However, differences in sedimentation rates generated differ- The distance from the meadow edge is another physical factor that can ences in CAR among all three depth zones. Specifically, CAR were statisti- influence the spatial heterogeneity in carbon within a meadow. For exam- cally significantly different at each depth zone (see Fig. B.2 in the ple, at a current-dominated site, Oreska et al. (2017) observed that sedi- Supporting Information). ment carbon increased with distance from the meadow's leading edge, Organic carbon accretion might differ between depth zones due to dif- and they attributed this to the attenuation of current within the meadow ferences in seagrass productivity. For example, Dennison and Alberte due to meadow drag. Specifically, when a current encounters a submerged (1986) observed that at the deepest sites within a Zostera marina L. meadow meadow with a sufficient shoot density, the velocity within the meadow de- at Woods Hole, Massachusetts, productivity was controlled by light avail- creases with distance from the leading edge over an adjustment length- ability for photosynthesis, whereas light availability was not a major factor scale that is dependent on the meadow height and shoot density (Chen in productivity differences in shallower regions of the meadow. However, et al., 2013). previous measurements in the Nahant meadow (April to October of 1999 At the Nahant site, neither sediment organic carbon nor CAR had a sig- and 2000) did not indicate a dependence of seagrass productivity on nificant correlation with distance to the nearest meadow edge, and the water depth (Colarusso, 2006). highest and lowest CAR occurred at similar distances to the nearest edge Furthermore, the three-fold spatial variation in sediment organic car- of the meadow (Fig. 5(a) and (c)). The meadow sites were not within the bon did not have a significant relationship with spatial variation in shoot expected current adjustment length-scale of the meadow (up to 2.5 m, density, and CAR had a negative correlation with shoot density (see based on Eq. (10) in Chen et al., 2013), such that the distance to the nearest Fig. B.3 in the Supporting Information). Specifically, CAR was highest meadow edge was not expected to be a factor. Furthermore, at this wave- where the shoot density was lowest (Table 1). While surprising for autoch- dominated site, the distance measured along the direction of wave propaga- thonous carbon, this is consistent with allochthonous carbon, which tends tion is more relevant than the distance to the nearest meadow edge. Since to be more significantly influenced by hydrodynamic factors. Recall that the waves transform across the meadow in response to changes in water Nahant has a significant allochthonous fraction (Fig. 8 in Novak et al., depth and to dissipation by meadow drag, one should consider the direction 2020). Specifically, within this wave-dominated meadow, the near-bed ve- of wave propagation and measure the distance from the offshore edge of the locity, which controls resuspension, has little dependence on local shoot meadow. The CAR had a significant relationship with distance to the off- density but instead depends on the integration of shoot density (which de- shore meadow edge (Fig. 5(d)). Specifically, CAR decreased with distance termines meadow drag) along the trajectory of the wave. However, shoot from the offshore edge. This was consistent with the increased near-bed density can be an important metric when comparing carbon accretion wave velocity associated with shoaling, which would tend to produce rates between different meadows. Denser and larger meadows diminish greater resuspension at sites closer to the shore, which in turn influenced wave energy to a greater degree and thus have a greater ability to curtail re- carbon accretion. That is, moving toward shore, the near-bed wave velocity suspension, which would enhance SAR and CAR (Tang et al., 2019). In con- increased, promoting resuspension and diminishing carbon accretion. Con- trast to an oscillating wave velocity, a steady current velocity is more sistent with this, Dahl et al. (2018) observed that resuspension removed 6 J. Lei et al. Science of the Total Environment 854 (2023) 158685 Fig. 5. (a,b) Sediment organic carbon and (c,d) carbon accretion rate (CAR) versus the distance from nearest meadow edge (a and c) and the distance from offshore meadow edge (b and d). Circles, triangles, and diamonds denote the Deep, Mid, and Shallow zones, respectively. The solid curves denote the lines of best fit. organic material from the bed within a meadow. Resuspension can also en- Eq. (17) in Mendez and Losada (2004) was adapted to describe the change hance oxygen within the bed, which would tend to diminish carbon reten- in wave amplitude over the meadow, tion through increased rates of mineralization by aerobic microbes aw ðiÞ (e.g., Ståhlberg et al., 2006). ¼ K v ði  1Þ K s ði  1Þ (5) aw ði  1Þ 3.4. Using hydrodynamic gradients to constrain meadow-scale estimates of in which aw(i), Kv(i) and Ks(i) are the wave amplitude, vegetation dissipation carbon and shoaling coefficients at distance i from the offshore edge of the meadow, respectively, with i corresponding to the progression of 1-m The carbon accretion rates were inversely correlated with hydrody- steps (Δx). The change in wave amplitude due to wave energy dissipation namic exposure expressed in terms of near-bed TKE (Fig. 6), with the fol- by the meadow can be described as (e.g., Dalrymple et al., 1984) lowing linear best-fit based on Storm 3 (90 % CI, Fig. 6): 1    2 K vðiÞ ¼ , (6) gC cm 1 þ K DðiÞ awðiÞ Δx CAR 2 ¼ ð  0:69  0:17ÞTKE þ ð21  3Þ: (4) m yr s2 in which the wave decay coefficient     At this wave-dominated site, the near-bed TKE can be predicted from 2 9 sinh kðiÞ le,m þ sinh 3k ðiÞ le,m K DðiÞ ¼ C D nbðiÞ b k ðiÞ       (7) wave mechanics. Specifically, wave evolution over the meadow can be de- 9π sinh k ðiÞ DðiÞ sinh 2k ðiÞ DðiÞ þ 2kðiÞ DðiÞ scribed by the conservation of wave energy, including the influence of shoaling and vegetation-induced wave dissipation on the wave amplitude, depends on water depth, D, wave number, k, number of leaves per bed area, aw. Bed-induced wave dissipation was assumed to be negligible, compared nb, and average leaf width, b. The drag coefficient CD = 2 for flat leaf mor- to vegetation-induced dissipation. Since both shoaling and wave dissipa- phology (see Luhar and Nepf, 2016). Because seagrass leaves are flexible, tion change with water depth, a marching solution was used to evaluate their tips can move with the waves, reducing the relative velocity between the wave amplitude at Δx = 1-m steps across the meadow. Specifically, the leaf and the water, which reduces the drag. The impact of this 7 J. Lei et al. Science of the Total Environment 854 (2023) 158685 edge (≈250 m in Fig. 7), Eq. (5) was evaluated in 1-m increments over the meadow, using Eqs. (6) through (9) to evaluate Kv (i−1) and Ks (i−1). Linear wave theory (Eq. (10)) was used to convert the spatial evolution of wave pffiffiffi amplitude to spatial variation in near-bed URMS, using U RMS ¼ U w = 2 (Fig. 7(a)). The near-bed URMS was used in Eq. (1) to estimate near-bed vegetation-generated TKE (Fig. 7(b)), assuming a smooth change in shoot density over the meadow. Finally, Eq. (4) was used to predict the spatial variation in carbon accretion rate (Fig. 7(c)). In this way, we used physical models of wave mechanics to interpolate the carbon accretion rate over the meadow. The estimation was made for two wave periods representing the bounds of the measured wave spectral peak (T=6 s and 10 s shown with black and red curves, respectively, in Fig. 7), but the solution was not sen- sitive to wave periods in this range. Differences between predictions and measurements could be attributed to the assumed linearization of the ba- thymetry and smooth variation in the meadow structure. In summary, sim- ple modeling tools performed well in predicting the observed spatial variation in hydrodynamic conditions and CAR. Next, we considered the results from Storm 3 (Fig. 7) in the broader con- text of conditions at the site. The average airport wind speed for August is 4.5 m/s, and the peak hourly airport wind speed on the stormy days (in- cluding Storm 3) was 9 m/s. The dominant wind direction is from the southwest, with a fetch of about 5 km. Assuming that the airport wind speeds were representative of the wind speed 10 m above the water surface, we used the fetch-limited Joint North Sea Wave Project (JONSWAP) param- eter equations (Kamphuis, 2010) to estimate a significant wave height Hm0 of Fig. 6. Carbon accretion rate (CAR) as a function of near-bed turbulent kinetic 16 cm (using a wind speed of 4.5 m/s) and 33 cm (using a wind speed 9 m/s). energy (TKE) during Storm 3. The black curve denotes the linear best fit (Eq. (4)). Wave heights estimated from the TCM velocity data during the stormy periods The shaded region shows the 90 % confidence interval bounds for the fit. were between 20 and 32 cm, confirming that the range of observed conditions were representative of long-term conditions at the site. To examine the impact of the meadow on the hydrodynamic condi- tions and CAR, we repeated the analysis without vegetation, setting reconfiguration is characterized by the effective meadow height, le,m, which K v =1 and n b = 0 (Fig. 7(d–f)). Note that at the offshore edge of the is the height of the meadow that contributes to the drag and wave dissipa- meadow (≈ 250 m in Fig. 7) the wave conditions are the same with and with- tion. The effective meadow height depends on the morphology and mate- out the meadow, because this is the point at which the waves first encounter rial properties of the leaves and on the wave conditions (Luhar and Nepf, the meadow. Moving toward shore (x < 250 m), the wave velocity and near- 2016; Luhar et al., 2017). Over a wide range of pure wave conditions, a sim- bed TKE would be significantly higher without the meadow, with wave veloc- ple but reasonable approximation for le,m is 10 % of the average shoot ity increasing to more than twice the observed value (Fig. 7(d)). The signifi- length (see Fig. 5 in Lei and Nepf, 2019). cant increase in near-bed TKE (Fig. 7(e)), suggests that CAR would decrease The group velocity, Cg, describes the speed of wave energy propagation to zero over most of the meadow (based on Eq. (4)), which emphasizes the im- (Kamphuis, 2010): portant role of wave dissipation by the meadow in creating conditions condu-   12 cive to carbon accretion. We caution that Fig. 7(e) is speculative because 1 2kD g Cg ¼ 1þ tanh ðkDÞ : (8) Eq. (4) is based on only three points and does not resolve the transition to 2 sinh ð2kDÞ k zero CAR at low values of near-bed TKE. The predicted spatial variation in carbon accretion rate CAR(i) (Fig. 7 As the water depth D decreases, Cg decreases which generally results in (c)) was used to calculate the total carbon accretion rate for the meadow, an increase in wave amplitude. This is captured by the shoaling coefficient, TCAR. Again using Δx = 1 m, Ks (Kamphuis, 2010),  1 TCAR ¼ ∑CARðiÞ Δx wðiÞ : (11) C gði  1Þ 2 Ks ¼ : (9) C gðiÞ w(i) denotes the meadow width at distance i m from Shallow1. The meadow width was assumed to increase linearly from 100 m at Shallow1 to 500 m at Based on the bathymetric survey, it was reasonable to assume a linear Deep3, based on meadow contours developed from boat surveys (Fig. 1). bottom slope across the meadow (Fig. 1). The leaf density (nb, number of Using Eq. (11), TCAR = 9.2 x 105 g C yr−1 = 1.0 short tons C·yr−1. For leaves per bed area) was estimated from the shoot density (Table 1) and av- comparison, estimates using individual point measurements at the Deep, erage number of leaves per shoot in each depth zone; both were assumed to Mid, and Shallow zones were, respectively, 1.5 short tons C yr−1 (50 % vary smoothly across the meadow and interpolated between measurement higher), 0.80 short tons C yr−1 (20 % lower), and 0.20 short tons C yr−1 sites. The initial wave amplitude at the outer edge of the meadow, aw,0, was (80 % lower). This comparison demonstrates that individual measurements calculated from linear wave theory and the measured near-bed wave veloc- of CAR do poorly in estimating the meadow average, and that the predict- ity amplitude, Uw, able variation in wave velocity over a meadow can constrain estimates of carbon accretion. 2π coshcosh ðkzÞ Finally, the potential for wave-driven resuspension to generate spatial Uw ¼ aw , (10) T sinhsinh ðkDÞ gradients in CAR at other sites can be assessed using wind data and median grain size (d50). Wind data and fetch provide an estimate of wave climate at in which T is the wave period, and z = 16 cm was the approximate distance the meadow edge, and the model described above can be used to describe above the bed corresponding to half of the vertical height of the TCM. For the gradient in wave velocity across the meadow. For resuspension to this analysis, we focused on Storm 3 conditions. Starting from the offshore drive spatial heterogeneity in sediment carbon, two conditions must be 8 J. Lei et al. Science of the Total Environment 854 (2023) 158685 Fig. 7. Evolution of (a, d) RMS wave velocity, (b, e) near-bed turbulent kinetic energy (TKE), and (c, f) carbon accretion rate (CAR) versus distance from the site Shallow1 (x = 0) during Storm 3. The top row includes the effects of the meadow on wave dissipation and near-bed TKE. The bottom row excludes the impact of the meadow by setting K v =1 and n b = 0. Black circles represent measurements at Shallow1, Mid2, and Deep3. Dashed curves denote predictions based on physical models and the observed correlation between near-bed TKE and carbon accretion (Fig. 6). Two wave periods were considered, T = 6 s (black dashed- dotted curve) and T=10 s (red dashed curve), to span the observed range of dominant wave periods at the site. met. First, there is a gradient in near-bed wave velocity, and thus in near- measurements could advance carbon credit assessment methods and im- bed TKE within the meadow. Second, the range of TKE within the meadow prove the design of targeted seagrass conservation and restoration projects. crosses the critical TKE threshold for resuspension, so that there will be a gradient in resuspension. The grain size, d50, can be used in Eqs. (2) and (3) to infer the critical turbulence level for resuspension. Spatial gradients CRediT authorship contribution statement in resuspension can also be observed with sediment traps. Jiarui Lei: Conceptualization, Methodology, Investigation, Formal analy- sis, Visualization, Writing – original draft, Writing – review & editing. Rachel 4. Conclusions Schaefer: Methodology, Investigation, Formal analysis, Visualization, Writ- ing – original draft, Writing – review & editing. Phil Colarusso: Conceptuali- Measurements in a Zostera marina L. meadow in Nahant Harbor, Massa- zation, Methodology, Resources, Investigation, Writing – review & editing. chusetts, showed that spatial variation in sediment and carbon accretion Alyssa Novak: Conceptualization, Methodology, Resources, Investigation, rates correlated with the spatial variation in hydrodynamic intensity, Writing – review & editing. Juliet C. Simpson: Conceptualization, Methodol- which was characterized by near-bed turbulence. Specifically, sedimenta- ogy, Resources, Investigation, Writing – review & editing. Pere Masqué: tion and carbon accretion rates both decreased as near-bed turbulence in- Methodology, Writing – review & editing. Heidi Nepf: Conceptualization, creased. This suggested that resuspension driven by turbulence provided Methodology, Resources, Funding acquisition, Investigation, Formal analysis, a constraint on the rates of accretion. At wave-dominated sites, near-bed Visualization, Writing – review & editing, Supervision, Project administration. turbulence estimated from wave velocity and meadow morphology may provide a framework for describing spatial heterogeneity in carbon accre- Data availability tion and carbon sequestration within seagrass meadows, which could facil- itate more accurate carbon measurements. More accurate carbon Data will be made available on request. 9 J. Lei et al. Science of the Total Environment 854 (2023) 158685 Declaration of competing interest marshes, and seagrass meadows. Conservation International, Intergovernmental Oceano- graphic Commission of UNSECO, International Union for Conservation of Nature. Infantes, E., Orfila, A., Simarro, G., Terrados, J., Luhar, M., Nepf, H., 2012. Effect of a seagrass The authors declare that they have no known competing financial inter- (Posidonia oceanica) meadow on wave propagation. Mar. Ecol. Prog. Ser. 456, 63–72. ests or personal relationships that could have appeared to influence the https://doi.org/10.3354/meps09754. Kamphuis, J.William, 2010. Introduction to coastal engineering and management. Advanced work reported in this paper. Series on Ocean Engineering. Vol. 30. Word Scientific Co. Pte. Ltd. Kennedy, H., Beggins, J., Duarte, C., Fourqurean, J., Holmer, M., Marbà, N., Middelburg, J., Acknowledgements 2010. Seagrass sediments as a global carbon sink: isotopic constraints. Glob. Biogeochem. Cycles 24, GB4026. https://doi.org/10.1029/2010GB003848. Keulegan, G.H., Carpenter, L.H., 1958. Forces on cylinders and plates in an oscillating fluid. This study was supported by Shell International Exploration and Pro- J. Res. Natl. Bur. Stand. 60 (5), 423–440. https://doi.org/10.6028/jres.060.043. duction through the MIT Energy Initiative. R. Schaefer was supported in Komar, P.D., Miller, C.M., 1973. The threshold of sediment movement under oscillatory water waves. J. Sediment. Petrol. 43 (4), 1101–1110. https://doi.org/10.1306/74D7290A- part by the National Science Foundation Graduate Research Fellowship un- 2B21-11D7-8648000102C1865D. der Grant No. 1745302. We thank the U.S. Environmental Protection Krishnaswami, S., Lal, D., Martin, J.M., Meybeck, M., 1971. Geochronology of lake sediments. Agency dive team for their assistance in the field: Brent England, Jean Earth Planet. Sci. Lett. 11, 407–414. https://doi.org/10.1016/0012-821X(71)90202-0. Brochi, Chuck Protzmann, and Eric Nelson. The bathymetry and meadow Lavery, P., Mateo, M., Serrano, O., Rozaimi, M., 2013. Variability in the carbon storage of seagrass habitats and its implications for global estimates of blue carbon ecosystem ser- mapping were collected by Mike Sacarny and Juliet C. Simpson under vice. PLoS One 8 (9), e73748. https://doi.org/10.1371/journal.pone.0073748. award NA18OAR4170105 from the National Sea Grant College Program Lei, J., Nepf, H., 2019. Wave damping by flexible vegetation: connecting individual blade dy- of the U.S. Department of Commerce's National Oceanic and Atmospheric namics to the meadow scale. Coast. Eng. 147, 138–148. https://doi.org/10.1016/j. coastaleng.2019.01.008. Administration. The International Atomic Energy Agency is grateful for Lei, J., Nepf, H., 2021. Evolution of velocity from leading edge of 2D and 3D submerged can- the support provided by the Government of the Principality of Monaco. opies. J. Fluid Mech. 916, 1–37. https://doi.org/10.1017/jfm.2021.197. Ved Ahuja assisted with processing sediment samples at Boston University. Lowe, R.J., Koseff, J.R., Monismith, S.G., 2005. Oscillatory flow through submerged canopies: 1. Velocity structure. J. Geophys. Res. Oceans 110 (C10). https://doi.org/10.1029/ 2004JC002788. Appendix A. Supplementary data Luhar, M., Nepf, H.M., 2016. Wave-induced dynamics of flexible blades. J. Fluids. Struct. 61, 20–41. https://doi.org/10.1016/j.jfluidstructs.2015.11.007. Luhar, M., Rominger, J., Nepf, H., 2008. Interaction between flow, transport and vegetation Supplementary data to this article can be found online at https://doi. spatial structure. Environ. Fluid Mech. 8, 423–439. https://doi.org/10.1007/s10652- org/10.1016/j.scitotenv.2022.158685. 008-9080-9. Luhar, M., Infantes, E., Nepf, H.M., 2017. Seagrass blade motion under waves and its impact on wave decay. J. Geophys. Res. Ocean 122, 3736–3752. https://doi.org/10.1002/ References 2017JC012731. Mendez, F.J., Losada, I.J., 2004. An empirical model to estimate the propagation of random breaking and nonbreaking waves over vegetation fields. Coast. Eng. 51, 103–118. Arias-Ortiz, A., Masqué, P., Garcia-Orellana, J., Serrano, O., Mazarrasa, I., Marbà, N., https://doi.org/10.1016/j.coastaleng.2003.11.003. Lovelock, C., Lavery, P., Duarte, C.M., 2018. Reviews and syntheses: 210Pb-derived sed- Novak, A.B., Pelletier, M.C., Colarusso, P., Simpson, J., Gutierrez, M.N., Arias-Ortiz, A., iment and carbon accumulation rates in vegetated coastal ecosystems – setting the record Charpentier, M., Masqué, P., Vella, P., 2020. Factors influencing carbon stocks and accu- straight. Biogeosciences 15, 6791–6818. https://doi.org/10.5194/bg-15-6791-2018. mulation rates in eelgrass meadows across New England, USA. Estuar. Coast. 43 (8), Chen, Z., Jiang, C., Nepf, H., 2013. Flow adjustment at the leading edge of a submerged 2076–2091. https://doi.org/10.1007/s12237-020-00754-9. aquatic canopy. Water Resour. Res. 49 (9), 5537–5551. https://doi.org/10.1002/wrcr. Oreska, M.P.J., McGlathery, K.J., Porter, J.H., 2017. Seagrass blue carbon spatial patterns at 20403. the meadow-scale. PLoS One 12. https://doi.org/10.1371/journal.pone.0176630. Colarusso, Phil, 2006. Natural and Stress Induced Changes in Non-structural Carbohydrate Con- Poppe, K.L., Rybczyk, J.M., 2018. Carbon sequestration in a pacific northwest eelgrass centrations in Eelgrass (Zostera marina L.). Northeastern University Doctoral dissertation. (Zostera marina) meadow. Northwest Sci. 92 (2), 80–91. https://doi.org/10.3955/046. Dahl, M., Infantes, E., Clevesjö, R., Linderholm, H.W., Björk, M., Gullström, M., 2018. In- 092.0202. creased current flow enhances the risk of organic carbon loss from Zostera marina sedi- Prentice, C., Poppe, K.L., Lutz, M., Murray, E., Stephens, T.A., Spooner, A., Hessing-Lewis, M., ments: insights from a flume experiment. Limnol. Oceanogr. 63 (6), 2793–2805. Sanders-Smith, R., Rybczyk, J.M., Apple, J., Short, F.T., Gaeckle, J., Helms, A., Mattson, https://doi.org/10.1002/lno.11009. C., Raymond, W.W., Klinger, T., 2020. A synthesis of blue carbon stocks, sources, and ac- Dalrymple, R.A., Kirby, J.T., Hwang, P.A., 1984. Wave diffraction due to areas of energy dis- cumulation rates in eelgrass (Zostera marina) meadows in the Northeast Pacific. Glob. sipation. J. Waterw. Port Coast. Ocean Eng. 110 (1), 67–79. https://doi.org/10.1061/ Biogeochem. Cycles 34 (2), 1–16. https://doi.org/10.1029/2019GB006345. (ASCE)0733-950X(1984)110:1(67). Ricart, A.M., York, P.H., Bryant, C.V., Rasheed, M.A., Ierodiaconou, D., Macreadie, P.I., 2020. Dennison, W.C., Alberte, R.S., 1986. Photoadaptation and growth of Zostera marina L. (eel- High variability of blue carbon storage in seagrass meadows at the estuary scale. Sci. Rep. grass) transplants along a depth gradient. J. Exp. Mar. Biol. Ecol. 98, 265–282. https:// 10, 1–12. https://doi.org/10.1038/s41598-020-62639-y. doi.org/10.1016/0022-0981(86)90217-0. Sanchez-Cabeza, J.A., Masqué, P., Ani-Ragolta, I., 1998. 210Pb and 210Po analysis in sedi- Duarte, C.M., Marbà, N., Gacia, E., Fourqurean, J.W., Beggins, J., Barrón, C., Apostolaki, E.T., ments and soils by microwave acid digestion. J. Radioanal. Nuc. Chem. 227, 19–22. 2010. Seagrass community metabolism: assessing the carbon sink capacity of seagrass https://doi.org/10.1007/BF02386425. meadows. Glob. Biogeochem. Cycles 24 (4), 1–8. https://doi.org/10.1029/ Serrano, O., Lavery, P.S., Rozaimi, M., Mateo, M.Á., 2014. Influence of water depth on the car- 2010GB003793. bon sequestration capacity of seagrasses. Glob. Biogeochem. Cycles 28 (9), 950–961. Fonseca, M.S., Cahalan, J.A., 1992. A preliminary evaluation of wave attenuation by 4 species https://doi.org/10.1002/2014GB004872. of seagrass. Estuar. Coast. Shelf Sci. 35 (6), 565–576. https://doi.org/10.1016/S0272- Ståhlberg, C., Bastviken, D., Svensson, B.H., Rahm, L., 2006. Mineralisation of organic matter 7714(05)80039-3. in coastal sediments at different frequency and duration of resuspension. Estuar. Coast. Fourqurean, J.W., Duarte, C.M., Kennedy, H., Marbà, N., Holmer, M., Mateo, M.A., Shelf Sci. 70 (1–2), 317–325. https://doi.org/10.1016/j.ecss.2006.06.022. Apostolaki, E.T., Kendrick, G.A., Krause-Jensen, D., Serrano, O., 2012. Seagrass ecosys- Tang, C., Lei, J., Nepf, H.M., 2019. Impact of vegetation-generated turbulence on the critical, tems as a globally significant carbon stock. Nat. Geosci. 5, 505–509. https://doi.org/ near-bed, wave-velocity for sediment resuspension. Water Resour. Res. 55 (7), 10.1038/ngeo1477. 5904–5917. https://doi.org/10.1029/2018WR024335. Gacia, E., Duarte, C., 2001. Sediment retention by a Mediterranean Posidonia oceanica Tanino, Y., Nepf, H.M., 2008. Lateral dispersion in random cylinder arrays at high Reynolds meadow: the balance between deposition and resuspension. Estuar. Coast. Shelf Sci. 52, number. J. Fluid Mech. 600, 339–371. https://doi.org/10.1017/S0022112008000505. 505–514. https://doi.org/10.1006/ecss.2000.0753. Tinoco, R.O., Coco, G., 2018. Turbulence as the main driver of resuspension in oscillatory Gacia, E., Granata, T., Duarte, C., 1999. An approach to measurement of particle flux and sed- flow through vegetation. J. Geophys. Res. Earth Surf. 123, 891–904. https://doi.org/ iment retention within seagrass (Posidonia oceanica) meadows. Aquat. Bot. 65, 255–268. 10.1002/2017JF004504. https://doi.org/10.1016/S0304-3770(99)00044-3. Zhang, Y., Tang, C., Nepf, H., 2018. Turbulent kinetic energy in submerged model canopies Granata, T., Serra, T., Colomer, J., Casamitjana, X., Duarte, C., Gacia, E., 2001. Flow and par- under oscillatory flow. Water Resour. Res. 54 (3), 1734–1750. https://doi.org/10. ticle distributions in a nearshore seagrass meadow before and after a storm. Mar. Ecol. 1002/2017WR021732. Prog. Ser. 218, 95–106. https://doi.org/10.3354/meps218095. Zhang, J., Lei, J., Huai, W., Nepf, H., 2020. Turbulence and particle deposition under steady Gruber, R., Kemp, W.M., 2010. Feedback effects in a coastal canopy-forming submersed plant flow along a submerged seagrass meadow. J. Geophys. Res. Oceans 125 (5). https://doi. bed. Limnol. Oceanogr. 55 (6), 2285–2298. https://doi.org/10.4319/lo.2010.55.6.2285. org/10.1029/2019JC015985. Howard, J., Hoyt, S., Isensee, K., Pidgeon, E., Telszewski, M., 2014. Coastal blue carbon: methods for assessing carbon stocks and emissions factors in mangroves, tidal salt 10