Flexural Rigidity and Shoot Reconfiguration Determine Wake Length Behind Saltmarsh Vegetation Patches

Vegetation patches play an important role in controlling sediment deposition in shallow aquatic environments such as coastal saltmarshes and fluvial systems. However, predicting deposition around vegetation patches is difficult due to the complexity of patch morphology and their dynamic interaction with the flow. Here we incorporate a biomechanical model, parameterized using field data, within a 3‐D computational fluid dynamics model which allows prediction of individual shoot reconfiguration within patches due to flow forcing. The model predicts velocity attenuation and bed shear stresses within the wake of the patch which agree spatially with accretion patterns measured in the field using terrestrial lidar. The model is applied to sparse patches of Suaeda maritima, located in saltmarshes of coastal habitats, to explore the role of (I) shoot distribution, (II) patch geometry, (III) shoot flexural rigidity, and (IV) bulk flow velocity in determining the length of the predicted wake region. We demonstrate that for Suaeda maritima, with intermediate rigidity, the vertical shear layer over the vegetation controls the length of the predicted wake region. Consequently, reconfiguration due to flexural rigidity strongly impacts on wake length, confounding the relationship between patch height and wake length. A simplified model for predicting wake length based on shoot reconfiguration is applied to the simulation data and shows good agreement. The results demonstrate that the observed wake characteristics can be well explained by intraspecific variability in flexural rigidity, thus demonstrating the importance of biomechanical traits in determining flow‐vegetation‐sediment interactions.


Introduction
Vegetation has a major impact on flow, sediment dynamics, and ecological processes within shallow aquatic environments and can control the evolution of fluvial, deltaic, and coastal sedimentary systems including coastal saltmarshes (Neumeier & Amos, 2006;Temmerman et al., 2005). Understanding the role of vegetation in saltmarsh development and evolution is crucial for managing natural sea defenses and mitigating against future sea level rise (Möller, 2006;Shepard et al., 2011;Temmerman et al., 2013). Saltmarsh vegetation communities develop through a process of colonization by pioneer species, which alter flow and sedimentation for succession to occur. Such sedimentation occurs behind isolated vegetation patches and controls the process of succession, enabling the establishment of species more sensitive to higher inundation frequency and shear stresses through facilitation interaction processes (Castellanos et al., 1994;Langlois et al., 2003;Tempest et al., 2015).
Previous research has shown that sediment deposition behind idealized vegetation patches may be explained using the hydraulic principle of flow separation and reattachment ; see Figure 1). Vegetation represents a porous blockage, which extracts energy from the flow as it crosses the patch, causing a wake to form behind, characterized by low bleed velocities and reduced turbulent kinetic energy (Nepf, 1999;Nicolle & Eames, 2011). Here, enhanced deposition may occur (e.g., Figure 2), depending upon local sediment and flow conditions (Chen et al., 2012;Shi et al., 2016;. Where the flow reattaches, turbulent kinetic energy is higher, and vortex shedding may occur, which limits sediment deposition (Chen et al., 2012;. Previous studies have defined the spatial extent of the wake sedimentation region using threshold values of both velocity and turbulent kinetic energy (Chen et al., 2012;de Lima et al., 2015;Yang et al., 2016). To date, most research has focused on the role of the horizontal shear layer and wake structure, formed by the horizontal flow separation around the patch, in determining sediment deposition extent (Figure 1a). In these cases, for given patch porosity and flow conditions, patch width determines wake length . However, for submerged vegetation, the vertical shear layer formed at the top of the patch (Figure 1b) may be of equal or greater importance, depending on patch geometry Liu et al., 2018). If wake length is controlled by the vertical shear layer instead of the horizontal shear layer, then shoot height, rather than patch width, will determine the wake length. Moreover, if the wake is controlled by shoot height, then shoot reconfiguration, the adjustment of the shoot into a more streamlined position to reduce drag (Luhar & Nepf, 2011;Nepf & Vivoni, 2000), will also affect the wake length. Therefore, we hypothesize that flexural rigidity is an important parameter in defining wake length for flexible, submerged vegetation patches.
Testing this hypothesis requires investigating the impact of flow and patch characteristics on wake length in an environment where such parameters may be quantified and controlled. Computational fluid dynamics (CFD) modeling provides an effective tool for this as it allows control of flow, patch geometry, and biomechanical properties which may be difficult to accurately constrain within field or flume environments. This permits testing over a wide parameter space within a controlled setup. However, currently such models either allow for accurate prediction of the motion of single shoots (e.g., Dijkstra & Uittenbogaard, 2010;Marjoribanks et al., 2014) or rely on a broader-scale representation of vegetation using a field-parameterized rigid cylinder drag force approach (Baptist et al., 2007;Temmerman et al., 2005;Tempest et al., 2015;Wu et al., 2017). Hence, there is a need for models which are field parameterized but also allow for prediction of flow-vegetation interactions at the shoot scale. Importantly, existing models of flow and deposition processes in the wake of vegetation patches cannot currently account for the impacts of shoot reconfiguration which, through altering the canopy height and shear layer characteristics, may affect patch wake length.
The aim of this paper is to test the hypothesis that flexural rigidity is important in determining wake length behind flexible, submerged saltmarsh vegetation (Suaeda maritima) patches. This is achieved through (i) developing a field-parameterized combined CFD-biomechanical model that can represent shoot reconfiguration; (ii) using it to investigate the relative importance of flexural rigidity in determining wake length, compared to the impact of patch geometry (patch width/height and shoot distribution) and flow conditions (velocity); and (iii) evaluating the performance of a simplified analytical model for vertical shear layers, adapted from Zong and Nepf's (2012) model for horizontal shear layers ( Figure 1a) and based on flexural rigidity and bulk flow characteristics, in predicting wake lengths observed in the field data and CFD simulations.

Field Data Collection for Model Parameterization and Validation
Our model plant species is Suaeda maritima, an annual halophyte plant common in the pioneer saltmarshes of the Mont Saint-Michel Bay (France, 48.631°N, −1.504°W; Figures 2a and 2c). Suaeda maritima grows from June to October with a peak height reached around September to October (Tessier et al., 2000). Our study area is the seaward part of a tidal meander where sparse patches of Suaeda maritima growing on the inner bar are directly exposed to flow crossing the inner bar during flooding. Wake deposits of several centimeter thickness develop downstream of each patch (Figures 2b and 2c) with respect to the flood direction during the spring tides of September and October when the vegetation is at its peak development. The In each case, the wake length is controlled by patch width (w) and deflected patch height (h d ), respectively. L Z&N is the width-dependent wake length, calculated , and L V is the equivalent height-dependent wake length calculated geometrically as wake deposits are composed of very fine sandy to muddy sediment, with a median grain size ranging from 30 to 100 μm. Owing to the meander configuration, ebb flow during the spring tides is concentrated in the outer bend and tends to circumvent the vegetated inner bar (Leroux, 2013). Consequently, flow velocity and suspended sediment concentration (SSC) remain much lower in the studied area during ebbing compared to flooding and no wake deposits develop in the opposite direction of the patch (see section 2.1.2 and Figures 2 and 3). Hence, apart from a uniform component of accretion related to ebbing that we account  (Lague et al., 2013). Green points are automatically classified vegetation patches using CANUPO (Brodu & Lague, 2012). The coordinate system indicated on the picture is used to compute the vegetation patch characteristics and wake deposit length. The location of the reference ADV is indicated (ADV ref  for in data processing, the development of wake deposits occurs solely during flooding, under unidirectional flow. Here, we used a combination of high-resolution Terrestrial Laser Scanner (TLS or lidar) surveys, in situ measurement of flow velocity, SSC, and plant sampling to both parameterize and validate the numerical model.

3-D TLS Data Analysis
Terrestrial laser scanning provides an ideal method for the large-scale, high-resolution survey of temporally variable sediment and vegetation characteristics that are required to parameterize and validate a patchscale numerical model. Daily TLS surveys took place during spring tides in October 2010 using a Leica Scanstation 2 placed at three different scanning positions. The setup and data processing were as used by Lague et al. (2013) using four targets anchored on consolidated marsh for local georeferencing between daily surveys. Scanning time available between tides and difficulties in scanning in the unconsolidated muddy channel deposits meant that scanning positions did not fully cover the area resulting in partial shadow effects (i.e., no data) in the area of densest vegetation patches. The point density varied spatially between 1.0 × 10 3 and 1.0 × 10 5 points/m 2 ( Figure 2).
We used the first survey on 4 October before the first spring tides to extract vegetation characteristics from the point cloud using various algorithms implemented in the open source 3-D point cloud processing software Cloudcompare V2.5 (General Public Licence, 2014). We used the CANUPO classification algorithm (Brodu & Lague, 2012) to automatically separate vegetation patches from the ground. The CANUPO algorithm is a supervised semantic 3-D classification method operating directly on irregular point clouds. It computes the 3-D geometry of the point cloud at various scales using a principal component analysis of the points coordinates. Once trained with manually selected samples, it results in a classification for each point associated with a metric of classification confidence. On the studied vegetation patches, we used a range of scales from 3 to 6 cm with a 0.5-cm interval. This resulted in classification accuracy better than 99% (Brodu & Lague, 2012). The classified vegetation point cloud was subsequently segmented into individual patches using a connected component algorithm that groups subsets of the 3-D point cloud which have a nearest point distance smaller than a given threshold (5 cm for this data set). This resulted in 216 individual vegetation patches (Figure 2b). We assumed a uniform flow direction for the study site and considered a fixed coordinate system in which the y axis was parallel to flow direction computed as the mean direction of the tail deposits ( Figure 2b). Patch width w, length l, and height h were measured along the x, y, and z directions to ensure a consistent approach for all vegetation patches with respect to the mean flow direction. The size characteristics w, l, and h correspond to the range between the 2nd and 98th percentiles of the coordinates (x, y, z) of all the points of a given patch. Using percentiles reduces the dependency on the occasional single shoot sticking out of the vegetation patch and better describes the geometry of the densest part. The observed patches of Suaeda maritima tended to develop mostly perpendicular to the flow direction ( Figure 2c). Consequently, the width was generally larger than the length. Wake sediment characteristics were extracted from the topographic change between 4 and 11 October of the ground points (supporting information Figure S1). The topographic change was measured vertically using the M3C2 algorithm dedicated to accurate change detection on 3-D point clouds (Lague et al., 2013). M3C2 operates directly on the two 3-D point clouds by averaging the point cloud position within a disk of specified diameter (chosen here as 20 cm) resulting in a map of topographic change (Figure 2b). By looking at topographic change we can define the wake as the zone where flow velocity and turbulent kinetic energy were low enough for sediment to settle. The map of topographic change shows that the ambient accretion rate away from wake deposits slightly varies along the flood direction such that a uniform threshold of accretion cannot simply be used (Figures 2d and 2e; Leroux, 2013). To account for the spatial variation of the ambient accretion rate, we first conducted a manual, coarse removal of the wake deposits highlighted by the M3C2 measurement on the 11 October data (Figure 2b). We used the resulting point cloud to interpolate the 11 October ground data in the location of wake deposits. We used M3C2 with a projection scale of 20 cm to measure the distance between the 11 October ground data and the interpolated point cloud. A distance threshold of 3 mm was chosen as it best distinguished the wake deposits behind vegetation patches. Due to the ground-based viewpoint of the TLS data, the vegetation patches generated significant shadow in the densest patch area and not all wake deposits could be properly detected. We thus manually selected wake deposits and associated vegetation patches in the sparsest area ensuring that the wake deposits could be clearly separated in the case of two nearby vegetation patches. This choice is consistent with the setup of the numerical simulation which explores a single patch. This restricted the analysis to 45 vegetation patches and their associated wake deposits (supporting information Figure S1). Mean accretion on the chosen wake deposits was 40 mm over the survey period with local peak accretion up to~80 mm (Figures 2d and  2e). The length of each patch was measured along the y axis from the center of the corresponding vegetation patch. The mean accretion rate outside the wake deposits was 20 mm.

Flow and Sediment Measurement
The field data were originally collected as part of a separate project, and therefore velocity data were not collected with detailed numerical modeling in mind. Consequently, flow and sediment measurements within the study area were limited. Two 6-MHz Nortek Vector acoustic doppler velocimeters (ADV) measuring the three components of the velocity field were installed between 6 and 10 October and set up to record 1-min burst measurements at 64 Hz every 10 min. One ADV was installed immediately behind a very large patch of vegetation (ADV wake ), and the other directly exposed to the incoming flood and away from the influence of any vegetation wake (ADV ref ; supporting information Figures S1 and S2). The ADVs were fixed on rigid frames and set up to measure the velocity at 0.06 m above ground. After the survey period, this distance was 0.05 m for ADV ref and 0.005 m for ADV wake illustrating the large difference in accretion between the two locations. ADV raw data were filtered to remove low signal to noise ratio data. SSC was estimated using the ADVs backscatter signal averaged over each 1-min burst measurement (Salehi & Strom, 2011) and previously calibrated against water sampling during one tide event (Leroux, 2013). Figure 3 shows a typical record during a tidal cycle. Peak velocity (~0.6 m/s) and SSC (~2.8 g/L) away from the wake of vegetation occurs during flooding. At the same time, the SSC in the wake zone is slightly reduced (~1 g/L), and the flow velocity strongly damped (~0.1-0.15 m/s). The combination of reduced velocity and high SSC is conducive of high accretion rates in the wake zone. During ebbing, the two ADVs exhibit similar characteristics with reduced peak flow velocity (~0.2 m/s) and a much lower SSC (~0.2-0.4 g/L). Leroux (2013) showed that the peak flood velocity away from the vegetation wake increases linearly with the tide amplitude (supporting information Figure S3), while the peak flood velocity in the vegetation wake is nearly constant, around 0.15 m/s.

Stem Geometrical and Mechanical Properties
As the biomechanical analysis of the plants was not initially planned when the TLS surveys started in 2010, it occurred 2 years later in July 2012 on vegetation patches located on the same area. The density of vegetation patches was higher, but the samples were chosen in the more exposed and sparse area to be as close as possible to the configuration of October 2010. The flexural rigidity of 196 stems sampled in two 3.5 × 3.5-m 2 squares was measured with a universal testing machine (Instron 5942 Canton, MA, USA). Basal stem samples (50 mm long) were tested as cantilever beams (one fixed end bending test; Hamann & Puijalon, 2013). Each sample was clamped horizontally at its basal end while a force was applied at its midpoint by lowering a probe at a rate of 10 mm/min. The following biomechanical traits were calculated: 1. The bending Young's modulus, E (Pa), quantifies the material stiffness and is calculated as the slope of the stress-strain curve in the elastic deformation region; 2. The second moment of area, I (m 4 ), accounts for the effect of the cross-sectional geometry of a structure on its bending stress. As stem cross section was circular, I was calculated as I = (πr 4 )/4, where r is the radius of stem cross section (Niklas, 1992); 3. The flexural rigidity, EI (Nm 2 ), quantifies the stiffness of the stem fragment and was calculated by multiplying E and I.
A second sampling took place in September 2013 to evaluate the shoot density and length in a vegetation state more similar to that when measurements were made in October 2010 and to constrain the equivalent porosity of the shoots for the numerical modeling setup.

CFD Modeling Methods 2.2.1. Vegetation Model Conceptualization and Parameterization From Field Data
Each vegetation patch was represented by a given planform area, scaled according to the field data. The patch characteristics measured from the TLS data ( Figure 4) show a linear relationship between patch width and patch length, such that l = 0.6w (r 2 = 0.87), while the relationship between the patch height and width exhibits a weaker power law trend, h = 0.25w 0.6 (r 2 = 0.59). Although there is some scatter, the data capture the mean geometrical characteristics of the patches. From the TLS data, these characteristics appear independent of the density of vegetation patches.
For each simulation, the given patch area (l × w, see Figure 5a) was populated with individual shoots matching the mean field-observed shoot density of 505 ± 202 shoots/m 2 measured over 15 patches of width greater than 0.1 m. The location of each shoot within the patch area was randomly selected from a uniform distribution in both the streamwise and spanwise directions, with the additional criteria that shoots could not overlap and must be separated in each direction by 5 mm (one numerical grid cell) of unoccupied space. Shoot properties were assigned using the characteristics observed in the field (Figure 5b). For shoot height, an average spanwise slope of 0.26 in shoot height (Δ Shoot height) was observed away from the patch center and therefore this was applied within the model such that h = 0.25w 0.6 was the maximum shoot height at the patch center. The stem flexural rigidity, calculated using samples from the field site, was EI=1.56 × 10 −4 ± 1.5 × 10 −4 Nm 2 . However, the distribution was not Gaussian and spanned approximately 2 orders of magnitude from 1.0 × 10 −5 to 1.0 × 10 −3 Nm 2 .
Each vegetation shoot was represented within the model as a dynamic porous blockage which reconfigured to the mean flow at each model iteration. Within each shoot, the individual main stem and leaf features were  ) and an additional drag force term within the Navier-Stokes equations which represented the mass and momentum effects of the vegetation respectively. Each vegetation shoot was conceptualized as a vertical stack of porous arrays, 5 × 5 numerical cells wide (0.025 × 0.025 m, 5Δx × 5Δx; see Figure 5c). This width was chosen based on the average width measured from plants sampled in the field. The center of each array corresponded to a single node on a beam which was used to calculate shoot reconfiguration and was fixed to the bed. Static reconfiguration of each shoot to the time-averaged flow was calculated by solving the steady Euler-Bernoulli beam equation (Marjoribanks et al., 2014), applying the total force acting over each porous array at each single vertical node. A relaxed iterative procedure, whereby shoot posture was updated at every flow calculation iteration, was used to ensure that the shoot positions converged simultaneously to the flow. As shoot position was directly coupled with the flow, the flow convergence criterion was sufficient to imply shoot position convergence. Shoot collisions were not treated within the model.
In order to parameterize the shoot porosity (ϕ) for the mass flux scaling algorithm, images of three characteristic shoots, taken in the field, were analyzed as follows ( Figure 6a). Within the image, the main stem of each shoot was manually identified and interpolated into a continuous profile. The image was then converted into a binary mask using RGB thresholding (Figure 6b; Boothroyd et al., 2017). This binary mask was converted into a shoot solid volume fraction (1 − ϕ) map ( Figure 6c) for each pixel, by sampling over a 0.005-m (Δx) moving window. Using this map, a vertical profile for each shoot was created (Figure 6d), by extracting a 0.025-m horizontal window centered on the interpolated main stem at each pixel up-shoot. The results of this process (Figure 6d) visually demonstrate that using a window width of 0.025 m captures the majority of each shoot with minimal overlap. Nevertheless, there is some inherent error in the process due to the close proximity of the shoots in the original photograph. The vertical profile was then resampled horizontally into five cells, each 0.005 m wide (Δx; Figure 6e). Using these profiles, the mean ( Figure 6f) and standard deviation (Figure 6g) of the shoot solid volume fraction were calculated over all three shoots for the three different relative cell positions shown in Figure 5c, representing grid cells at the shoot center (1), intermediate grid cells (2), and grid cells at the shoot boundary (3). Each porous array was randomly populated according to its position relative to the shoot center assuming a normal distribution of porosity with positiondependent mean and standard deviation values. The porosity (ϕ) was applied at grid cell faces as outlined in Marjoribanks, Hardy, Lane, and Tancock (2017).
The drag force acting within each cell was calculated based upon the simplified assumption that all plant components were cylindrical and therefore for each the drag coefficient (C D ) was equal to 1 (Panton, 1996). Suaeda maritima typically exhibit a range of leaf forms from semiterete to terete (cylindrical; Tison & de Foucault, 2014), and Figure 6a shows that many of the stems and leaves were approximately cylindrical in shape. Drag coefficient will vary between stems/leaves and with orientation, and so choice of drag coefficient may impact upon the accuracy of results. However, through explicitly representing the large-scale impacts of shoot reconfiguration and morphology, we have reduced the impact of drag coefficient to the subgrid leaf and stem scale, compared to other methods where all these effects are represented through a single drag coefficient.
The drag force within every grid cell was calculated by summing up the drag due to each cylindrical stem and leaf within the cell, assuming stems and leaves were perpendicular to the flow: where n is the number of stems and leaves, ρ (kg/m 3 ) is the fluid density, C D is the drag coefficient, h (m) is the stem height, d (m) is the stem diameter, and U (m/s) is the local grid cell velocity.
In order to relate the drag force to porosity, four further simplifying assumptions were made about the stems and leaves. First, it was assumed that stems and leaves had the same diameter (d). Second, it was assumed that they extended across either the full width (h = Δy) or height (h = Δz) of a grid cell, depending on orientation. Third, it was assumed that stems and leaves within a grid cell did not overlap, and fourth, it was assumed that they did not induce any within-cell sheltering effects. This allowed calculation of the number of stems per grid cell volume from the cell face porosity (ϕ) as follows: Here Δy (m) could be used interchangeably with Δz (Δy = Δz) depending on stem/leaf orientation. Therefore, combining equations (1) and (2), the drag force within each grid cell could be expressed as follows: This additional drag term was incorporated within the Navier Stokes equations using a linearized source term approach as outlined by Marjoribanks, Hardy, Lane, and Tancock (2017).

Flow Model
Flow was modeled using the Reynolds-averaged Navier-Stokes (RANS) equations with a Re-Normalisation Group theory derived k − ϵ turbulence closure model to produce a time-averaged, steady state solution. The RANS equations were solved using a staggered grid and a hybrid-upwind differencing scheme. The mass and momentum equations were coupled using the SIMPLEST algorithm (Spalding, 1980), and the solution was 10.1029/2019JF005012 Journal of Geophysical Research: Earth Surface solved iteratively until errors were reduced to 0.1% of the inlet flux. The numerical domain was scaled according to patch width to ensure comparability between scenarios, such that the domain spanned 24w in the downstream direction, 5w in the spanwise direction, and 0.75w 0.6 (3h) in the vertical direction to minimize and standardize blockage ratio and boundary effects. Thus, the blockage ratio of the vegetation patch was 0.2, though this does not consider the porosity of the patches themselves. The vegetation patch was centered at x = 2w, y = 2.5w, where in contrast to the TLS processing, the conventional frame of reference with x as the downstream coordinate and y as the cross-stream coordinate was used. For the section of the domain containing the patch, a regular grid of resolution Δx = Δy = Δz = 0.005 m was used. Downstream of the wake region (x > 4w), the x direction grid spacing was increased gradually for computational efficiency. Thus, the total number of grid cells for each simulation was between 0.46 and 38.6 × 10 6 and was 7.45 × 10 6 for the standard (w = 0.3 m) case. The computational time for each simulation scaled with the number of grid cells and was 30 hr for the w = 0.3-m case. A standard rigid lid treatment was applied at the free surface, while the bed was represented using a no-slip boundary condition and logarithmic wall function (dimensionless wall distance, y þ >34 for all cases). The sides of the domain were treated as frictionless boundaries. A constant inlet velocity was applied as only single-point velocity field data were available, and the outlet was represented as a fixed pressure boundary.

Wake Deposition Length Measurement
The field data presented here consist predominantly of high-resolution sediment deposition records, and therefore in order to compare these with the numerical flow data, we calculated wake lengths which we suggest correspond to potential regions of enhanced deposition. In this study we used three different measures of wake length. Primarily, wake deposition was related to areas of the bed where the bed shear stress did not exceed the critical bed shear stress, τ c (Pa). The bed shear stress, τ 0 (Pa), was calculated using the turbulent kinetic energy per unit mass k (m 2 /s 2 ) and may be written as The turbulent kinetic energy is a function of the fluctuating velocity components (u′, v′, w′) and is time aver- . It is calculated directly by the Re-Normalisation Group k − ϵ turbulence closure model by solving a transport equation, and the near-bed value was evaluated in the first grid cell above the bed (z = 0.0025 m). Recent studies have shown that near-bed turbulent kinetic energy is a good predictor of sediment transport within vegetated channels (Tinoco & Coco, 2018;Yang & Nepf, 2018).
Previous field experiments within oceanic environments have found C 1 = 0.19 − 0.2 (Pope et al., 2006;Soulsby, 1981), whereas a theoretical derivation, assuming approximate local equilibrium between turbulence production and dissipation, gives a value of C 1 = 0.3 (Baranya et al., 2012;Rodi, 1980). Here we used C 1 = 0.3 in line with the CFD model wall parameterization. The critical shear stress was τ c = 0.131 N/m 2 and was chosen to represent the threshold for motion for the largest grain size fraction (0.1 mm) applying the adapted Shields (1936) method of Soulsby (1997). Wake length, L TKE , was determined by the furthest point downstream where the shear stress criterion was met. Sediment cohesion was not considered here as the model was used to predict reentrainment of mobile grains rather than erosion.
In addition to the shear stress method, wake length was also calculated using the theoretical basis illustrated in Figure 1a. Using this method, for horizontal shear layers, wake length may be calculated theoretically as the convergence point of the two shear layers formed at the edge of the patch (Chen et al., 2012;. Here w (m) is the patch width, U (m/s) is the arithmetic mean velocity of the mixing layer, ΔU (m/s) is the velocity difference across the mixing layer, and S δ is the rate of growth of the mixing layer half-width, as measured from the starting location of the shear layer (i.e., any movement of the center of the shear layer is represented within the growth rate). However, considering the vertical shear layer (Figure 1b), equation (5) can be altered to calculate the distance downstream at which the vertical shear layer reaches the bed.
The second term in equation (6) is directly analogous to the Zong and Nepf (2012) model, with patch halfwidth replaced by the deflected shoot height, h d (m), while the first term x d (m) accounts for the fact that the downstream location of the shear layer will alter for flexible vegetation (see Figure 1b). This is similar to the term added by Hu et al. (2018) to account for displacement of the wake limit due to leaf length. For any given patch, both L Z&N and L V can be calculated and comparing the two predicted wake lengths enables classification of each patch wake length as either width dependent (L V > L Z&N ) or height dependent (L V < L Z&N ). In other words, whichever shear layer converges furthest upstream will limit the length of the wake. In order to compute L V and L Z&N , values of S δ were calculated for both the vertical and horizontal shear layers, along with the case-specific U and ΔU values. Shoot height varied within each simulated patch, so the maximum value within the patch was used as this has been shown to represent the effective canopy height (Hamed et al., 2017). To calculate S δ , the shear layer growth rate, the width of the shear layers (δ x , δ z ) were calculated at 0.05-m intervals downstream of the patch. The edge of the shear layer was measured quantitatively as the point z (or y for lateral shear layers), at which U(z) = U 1 +0.1ΔU, where U 1 is the mean velocity of the slower of the two coflowing regions forming the mixing layer (Pope, 2000;. This is a standard method for calculating shear layer width but is not directly related to a specific turbulent kinetic energy threshold. Therefore, it is possible that the spread of the critical threshold in turbulent kinetic energy, which in turn affects the critical shear stress, may occur at a smaller or larger rate. As the simulations were steady state and run at a constant velocity, while the field-site experienced tidal variations, an important consideration was what model velocity should be applied to most effectively reproduce the observed deposition patterns. The study area is a depositional environment, which experiences predominantly unidirectional but varying flow conditions and which exhibits net deposition over individual tide cycles (Leroux, 2013). Therefore, under many velocity conditions, wake sedimentation may not be distinguishable from deposition on the bare unvegetated bed. This phenomenon was studied by Shi et al. (2016), who found both sediment and flow controls on preferential wake deposition. Determining a threshold condition for preferential wake deposition is key to comparing predicted wake deposition extents with field data. SSC field data for vegetated and bare areas provided some detail on the threshold for preferential wake deposition. Figure 3 shows an example of the variation in SSC over an individual tide cycle. The largest difference in SSC between the patch and freestream flow, which we suggest relates to the period of most preferential deposition, occurs at the highest velocities (0.4-0.6 m/s). Below 0.4 m/s, at the end of the tide cycle, the SSC converges. However, the velocity was measured close to the bed (z = 0.06 m) and therefore likely underestimates the mean velocity. Therefore, the primary velocity for investigation was chosen as 0.7 m/s. However, other velocities (0.2-0.6 m/s) were also used for comparison. A single velocity condition, rather than an entire tide cycle, was selected to better highlight, and isolate, the role of velocity and other variables in controlling wake length.

Summary of CFD Simulations
To achieve the aim of this study, simulations were conducted with varying patch (width, shoot configuration, and flexural rigidity) and flow (velocity) conditions. When varying each of the different controlling variables, all other parameters were kept constant. The baseline comparison values used were a patch width of 0.3 m, a shoot flexural rigidity of 1.5 × 10 −4 Nm 2 , and a velocity of 0.7 m/s. The baseline width and flexural rigidity values correspond to the mean values for the field data, while the velocity was chosen as it corresponds approximately to the maximum preferential deposition condition as explained in section 2.2.3. Since the vegetation model randomly assigns shoot positions within each patch, when comparing patches with different widths and therefore different shoot distributions, it is important to understand the influence of subpatch shoot distribution on wake flow patterns which may confound any relationship with patch characteristics. Therefore, 10 different replicate simulations were conducted using identical baseline biomechanical, patch geometry, and velocity conditions with different random shoot distributions to quantify the impact of shoot distribution on wake length. For all other simulations with a patch width of 0.3 m, an identical shoot distribution was used. In total, results from 27 different simulations are presented here, split into four different comparison groups, reflecting investigation of different controlling variables: (I) shoot

Simplified Analytical Model
For patches with horizontal shear layers, Zong and Nepf's (2012) method provides a simple model for predicting wake length. Assuming instead that wake length is controlled by the vertical shear layer and consequently is a function of shoot height, then equation (6) can be used as an equivalent simplified model for predicting wake length. The key difference between this model and that of Zong and Nepf is that, in contrast to patch width which will remain constant for vertically aligned shoots, in the case of the vertical shear layer, both x d and h d will vary with hydrodynamic forcing and flexural rigidity. By simplifying patches to a single shoot with identical length and biomechanical properties to that used in the CFD simulations, a steady state Euler-Bernoulli beam equation was used to calculate h d and x d and consequently predict wake length using equation (6).
For this simplified model, the inlet velocity was constant with depth and set equal to the CFD inlet condition (U = 0.7 m/s). Similarly, the model used a constant expression for drag force along the shoot, based upon the mean inlet velocity.
The drag coefficient was taken as C D = 1, while the area of each section of the shoot, A (m 2 ), was adjusted with shoot inclination to represent the projected shoot frontal area. As with the CFD model, there are limitations with assigning a drag coefficient equal to 1; nevertheless, for the simplified model, while the drag coefficient may affect the exact wake lengths calculated, it should not substantially affect the trends produced by the model for varying rigidity. The spatial resolution of the simplified model was set equal to that used in the CFD model (Δz = 0.005 m). The computational time for the simplified model under baseline conditions was 0.15 s.
As discussed in section 2.2.3, the representative shear layer growth rate S δ for calculating wake length using equation (6) will vary depending on the threshold turbulent kinetic energy and consequently bed shear stress. Therefore, different values of S δ notionally represent different bed shear stress thresholds and consequently different grain size erosion thresholds. Accordingly, sensitivity to this parameter was examined through using four different values for S δ . Values of S δ were chosen based on the relationship between L TKE and L V within the CFD simulations in Group I. Average values of U=ΔU at each rigidity were obtained a priori from the CFD simulations with values varying between 0.77 and 1.03. However, using a single averaged value of 0.98 produced a similar result. The output from the simplified model was compared to wake lengths (L TKE ) calculated from the CFD model. In addition, the model was applied to the field data using the observed flexural rigidity value. However, the field data for flexural rigidity showed considerable uncertainty. Therefore, we also used the model to calculate the expected flexural rigidity value for each patch given its observed wake length. This was done by running the model over a wide parameter range for flexural rigidity (0.01-9.5 × 10 −4 Nm 2 ) for each patch and selecting the flexural rigidity with the lowest wake length error.

General Description of the Flow
Simulation results (Figure 7) demonstrate the typical patch layout and shoot reconfiguration across a range of flexural rigidities for the baseline geometry and velocity conditions (w=0.3 m, U=0.7 m/s). Shoots reconfigure to the time-averaged flow, due to the drag force acting on the patch, demonstrated by the pressure gradient across the patches (Figures 7a-7c). Above the patch, horizontal and vertical shear layers form, which grow downstream, as demonstrated by the turbulent kinetic energy distribution (Figure 7d). Behind the patch, a wake forms, where flow velocity (Figure 7e) is reduced. The extent of this wake is determined by the rate of growth of the shear layers.

Comparison With Field Data
The simulations were designed to replicate the ensemble characteristics of the Mont Saint-Michel vegetation, rather than individual plant patches. Nevertheless, it is important to compare the simulation results,

10.1029/2019JF005012
Journal of Geophysical Research: Earth Surface both velocities and wake lengths, to see whether they fit within the observed parameter space and to validate the assumptions used in modeling. Limited flow data from behind and adjacent to an individual patch were available for one similar but slightly lower flow condition (U z=0.06m =0.61 m/s). These show that observed wake velocities were comparable (blue cross in Figure 8a) to those simulated. The shape of the longitudinal velocity profile in Figure 8a also agrees well with previous studies of flow around porous  obstacles (e.g., Chen et al., 2012;de Lima et al., 2015) showing a decrease in velocity just in front of the patch, a large velocity deficit in the wake and gradual recovery that extends far downstream.
The simulated wake lengths for baseline conditions (Figure 9a, red crosses) lie within the observed data range across a variety of patch widths, suggesting the ensemble parameterization of patch geometry, threshold velocity, and flexural rigidity is broadly consistent. However, as discussed above it is not possible to assess the accuracy of the model for individual patch wakes.

Evaluation of Shoot Distribution Effects (Group I)
The impact of the different shoot distributions on patch centerline velocities is shown by the standard deviation in Figure 8a, which is largest in the wakes, where the shoot-scale heterogeneity introduces variations in the velocity field between replicates. The maximum standard deviation occurs immediately downstream of the patch with a value of 0.12U 0 . However, this drops to 0.045U 0 by two patch widths (2w) downstream of the patch and to 0.03U 0 by three patch widths (3w) downstream. Comparing the trend in standard deviation over different numbers of replicate simulations (Figure 8b) demonstrates that for all three wake measures, the standard deviation levels off by 10 replicates.
The three different wake length metrics (Figure 8c) show several important features. First, the patch wake length is clearly height dependent, with L V < L Z&N for all shoot configurations. Therefore, it is most relevant to compare the theoretical value for the vertical shear layer, L V , with L TKE . There is a difference between the mean values for L V and L TKE with a ratio of approximately 1.5. The two measures also show a difference in range, L V having a range of 0.62w compared to 0.92w for L TKE . The larger range of values for L TKE , and larger standard deviation, may be explained by the lower mean value for wake length, meaning it occurs in a region closer to the patch where flow variability is higher (Figure 8a). Comparison with two other wake length metrics commonly used in the literature can be found in the supporting information (Text S1 and Figure S4).

Geometric (Width) Controls on Wake Length (Group II)
The previous section showed that for the baseline case, the wake length was height dependent. Here we broaden the analysis to the full range of patch widths investigated (w = 0.1 − 0.6 m). Considering the theoretical models for shear layer growth (Figure 1), for a given flow velocity and vegetation porosity the two geometric variables that determine wake length for rigid shoots are patch width and shoot height. From Figure 1 b it can be seen that for the vertical shear layer, the maximum half-width to which the shear layer can grow, before the shear layer interacts with the bed, is equal to the rigid shoot height, δ MAX (h) = h. For the horizontal shear layers, the maximum distance each shear layer can grow laterally before interfering is δ MAX (w) = w/2. The controlling factor for wake length will therefore depend on the ratio of w/h (Hu et al., 2018; Liu et al., 2018). If the vertical and horizontal shear layer characteristics ( U=S δ ΔU ) are identical between the horizontal and vertical shear layers, then it follows that the smallest value of δ MAX will dictate the wake length. Given the geometric relationship between w and h within the CFD model, based on the field data (Figure 4), for all patch widths greater than 0.177 m, the theoretical models suggest that wake length is controlled by shoot height (δ MAX (h) < δ MAX (w); Figure 10a).
For flexible vegetation, the deflected shoot height replaces the rigid height (δ MAX (h d ) = h d ). The height dependence in Figure 10a is accentuated by shoot bending which further reduces the shoot height and therefore maximum shear layer width. When considering the simulated deflected shoot heights, δ MAX (h d ) < δ MAX (w) for all patch widths simulated. The replicates for multiple shoot distributions at w = 0.3 m suggest that this trend is independent of shoot distribution, with a range in maximum deflected shoot height of 0.01 m between different shoot distributions. However, for deflected vegetation, equation (6) shows that the location of the shear layer will be displaced downstream (x d ) and therefore this will counteract the reduction in shear layer width, lengthening the wake. Therefore, the value of x d and the magnitude of the shear layer growth rate may still cause the wake to be width dominated for patches where x d is large.
The preceding analysis assumes that the vertical and horizontal shear layers have a similar growth rate, which may not necessarily be true given the difference in patch boundary roughness, foliage, and shoot reconfiguration (e.g., Hu et al., 2018). However, simulated growth of shear layer width (δ x , δ z ), quantified using the method described in section 2.2.3, shows that for an example case, both the vertical and horizontal shear layer growth rates are similar (Figure 10b). On average, after initial variation, both can be well explained by a growth rate of S δ = 0.1. Zong and Nepf (2012) also found S δ ≈ 0.1, which is within the broad range of theoretical mixing layer spreading rates observed in previous shear layer studies S δ = α= 0.06-0.12 (Dimotakis, 1991;Marjoribanks et al., 2014;Pope, 2000;Sukhodolova & Sukhodolov, 2012). Therefore, it is justifiable to state that for rigid shoots, patches with w> 0.177 m will have height-dependent wake lengths. For flexible shoots, further analysis is required to examine the impact of reconfiguration on h d and x d and therefore wake length.

Effect of Flexural Rigidity on Wake Length (Group III)
Given the dependence of wake length on shoot deflection, it is important to understand how flexural rigidity changes deflected shoot height and consequently wake length. The CFD simulation output (Figure 7) shows the clear impact of three different flexural rigidities on patch posture, with the patch becoming increasingly streamlined with lower flexural rigidity. There is also a change in pressure distribution around the patch with larger zones of high and low pressure at the front and rear of the canopy for the more rigid cases, indicating higher patch-induced drag.

Journal of Geophysical Research: Earth Surface
Simulations across eight different flexural rigidities show a clear relationship between flexural rigidity and wake length ( Figure 11). As the L TKE definition of wake length depends directly on critical shear stress, wake lengths for six different values of critical shear stress, including the field estimated value of τ c = 0.13 N/m 2 , are plotted for comparison (Figure 11b). The trends for all critical shear stress values are similar. At low rigidities (EI<3 × 10 −4 Nm 2 ), wake length increases with rigidity. A maximum wake length then occurs at approximately EI= 3 × 10 −4 Nm 2 , above which the wake length decreases before converging to a constant value. The wake length is constant for high rigidities (EI > 1.5 × 10 −2 Nm 2 ), as most shoots have reached their fully rigid positions (Figure 11d). Critical shear stress has a variable effect on wake length, depending on flexural rigidity, with greater variability evident at higher rigidities. This suggests that the longitudinal gradient of bed shear stress is steeper for highly flexible canopies.

Effect of Velocity on Wake Length (Group IV)
The results with varying velocity (Figure 12) show that velocity in the wake region, downstream of the individual shoot wakes, is relatively independent of inlet velocity, which matches the results from the limited field ADV data (supporting information Figure S3). Varying inlet velocity for a constant flexural rigidity ( Figure 12) has a similar impact on patch reconfiguration and deflected patch height to the case in the previous section with constant velocity and variable rigidity (Figures 11d-11f). This is because both represent an increase in Cauchy number, which is the ratio of drag to rigidity forces and has been shown to determine shoot reconfiguration (see Luhar & Nepf, 2011). However, the impact of velocity on wake length is more complicated as variations in inlet velocity affect the bleed flow and therefore shear layer characteristics ( ΔU; U ). Furthermore, at low velocity cases, turbulent kinetic energy levels are sufficiently reduced (Figure 12b) such that preferential wake deposition may not occur. Wake length, L TKE , decreases with inlet velocity (supporting information Figure S5), but this trend is not discussed further here because of the additional complexities outlined above.

Simplified Analytical Model Results
Based on the observed average relationship between L V and L TKE in the CFD simulation data (Figure 8c), within the simplified analytical model, values of S δ ≈ 0.15 − 0.16 were chosen such that the wake lengths  (Figure 13a) produces a graph similar in shape to Figure 11 and predictions that visually agree well with the simulated wake lengths for the field-estimated critical shear stress, τ = 0.131 N/m 2 (black line in Figure 13a) particularly at lower rigidities (EI<2.0 × 10 −4 Nm 2 ). At higher rigidities, the simulated decrease in wake length is underrepresented by the simplified model. The results with different shear layer spreading rates, S δ , show that the impact of this parameter differs with flexural rigidity, with the difference in wake length between the four cases varying between 0.1w at EI = 1.5 × 10 −5 Nm 2 and 0.44w at EI = 3 × 10 −4 Nm 2 . The simulation data plot within the range S δ = 0.14-0.17 but do not follow the trend of a single value. Instead, the simulation data agree best with S δ = 0.15 at intermediate rigidities and S δ = 0.17 at both high and low rigidities.
Applying the simplified model, with S δ = 0.15, to each of the 45 field-surveyed patches with an identical flexural rigidity of 1.5 × 10 −4 Nm 2 produced an average error in predicted wake length of 0.22 m, which is expected given the local variability in both flexural rigidity and velocity. However, using the analytical model  and observed wake length to predict the flexural rigidity values for each patch (Figure 13c) produces a set of values which have a similar range to those observed in the field. We note that the relationship in Figure 13a can produce multiple possible flexural rigidities for a given wake length. In this analysis, consistent with the majority of the field data (Figure 11c), the lower rigidity solution was used. The simplified model predicts a greater relative occurrence of higher rigidity patches compared to the field measurements ( Figure 13c). This may reflect the dominance of higher rigidity shoots in determining patch wake length, above the influence of lower rigidity shoots also present within the patch. Given that the observed values in flexural rigidity vary over an order of magnitude, the simplified model's reproduction of this range does not provide direct evidence of the model's predictive capability but suggests the model may be capturing the key dynamics determining wake length.

A Simplified Model for Vertical Shear Layers
In contrast to the CFD model, which provides detailed parameterization of shoot, leaf, and patch properties, as well as solution of the three-dimensional flow field, the simplified analytical model for vertical shear layers, based on equation (6), relies upon a single shear layer term, U=S δ ΔU to encapsulate all of the flow and patch complexity. In this study, this parameterization was informed by the CFD simulations, but given approximate flow (velocity) and shoot (rigidity, diameter) data, site-specific models like that shown in Figure 13 could easily be produced following the same methodology and, except in the limiting case of very low velocity, would be expected to show the same overall trend with flexural rigidity.
While the model of Zong and Nepf (2012) produces a linear trend with increasing patch width, when considering the vertical shear layer, the relationship between wake length and shoot height is not linear, as the values of x d and h d are a function of flexural rigidity and determined by solving a beam equation.
The model results (Figure 13a) can be understood by examining the theoretical model ( Figure 1b) and the relative impacts of the terms x d and h d on wake length. It is clear from equation (6) that wake length increases with both x d and h d ; however, these two values are inversely related through the beam deflection model. The model results show a sharp increase in wake length at low rigidities, where the plant is in a more prone position, which is due to a large increase in h d as rigidity increases, compared to a small decrease in x d ( Figure 13b). As flexural rigidity increases, the relative change in h d decreases and x d starts to decrease more significantly, reflecting the fact that in more upright positions, shoot reconfiguration leads to greater changes in x d than h d (Figure 13b). This causes a balancing out between the impacts of h d and x d and therefore the wake length reaches its maximum and then decreases slightly as the decrease in x d begins to dominate the increase in h d for rigidities above EI= 3 × 10 −4 Nm 2 . Due to the shear layer growth term, S δ , in equation (6), the influence of h d on L V is approximately 7 times that of x d ; hence, the reduction in wake length at higher rigidities is relatively small. In contrast, the CFD simulation data (black line in Figure 13a) demonstrates a larger decrease in wake length at higher rigidities, suggesting the influence of other factors which are not accounted for in the simplified analytical model. Such factors may include the additional impacts of streamlining on flow structure or changes to bleed flow and porosity within the patch as the plant reconfigures, but these factors seem to have most significant impact at higher rigidities. While the simplified model cannot capture such complexities, it is able to capture the impact of shoot reconfiguration and requires less than 0.0002% of the computational time compared to the CFD model enabling rapid investigation over a wide parameter space.

Dependence of Wake Length on the Vertical Shear Layer and Flexural Rigidity
The field data and CFD results on shear layer growth demonstrated that for pioneer saltmarsh vegetation (Suaeda maritima) patches, considering the ensemble-averaged patch geometry, patch height controls wake length rather than patch width. This geometric control is similar to what Liu et al. (2018) discovered for submerged rigid vegetation patches. There are, however, key differences with their study. First, their experiments were conducted with a high patch flow blockage, decreasing bleed flow and wake length (Chen et al., 2012;Liu et al., 2018). Second, here we consider flexible vegetation, and therefore as the plant reconfigures, h d and x d become significant. Finally, the variability in flexural rigidity across our data set confounds any direct scaling between patch height, h, and wake deposition length, L TKE .

Journal of Geophysical Research: Earth Surface
This is evident by looking at the individual field data for the observed vegetation patches. Assuming similarity of shear layer growth as observed in the CFD simulations (Figure 10), it is possible to identify the theoretical controlling variable for each patch individually under rigid conditions. From the 45 patches, only five exhibit theoretical dependence on patch width (w < 2h). We note that this classification does not account for reconfiguration as this is not possible to quantify for the individual field patches. For the wakes determined by patch width (Figure 9b) there is a visual positive trend; however, there are too few points to produce a robust correlation. For the wakes determined by patch height (Figure 9c), there is less of an obvious trend with wide scatter. This could be explained by variability in shoot height or density within the patch or local velocity conditions. However, the simplified model demonstrates that it can also be explained by the dependence of deflected patch height on flexural rigidity as shown in Figure 13c.
The field data for flexural rigidity show an order of magnitude of variation in values. Previous studies of both freshwater and saltwater vegetation have also demonstrated large intraspecific variation in values for flexural rigidity between samples on different stems from the same site (Feagin et al., 2011;Miler et al., 2012;Paul et al., 2014;Rupprecht et al., 2015), and it is also likely that flexural rigidity varies depending on patch exposure to the flow and the local age of vegetation patches (Anderson & Smith, 2014;Silinski et al., 2018). The simplified model results (Figure 13c) demonstrate that the range of predicted flexural rigidities agrees well with the observed range. However, they also demonstrate the large variation in wake length caused by flexural rigidity which, over the range of flexural rigidity values observed, represents an uncertainty in wake length of 1.55w. Considering the range of velocities that exhibited preferential deposition (0.4-0.7 m/s), the simulation results suggest uncertainty in the wake length, due to local variations in velocity, of 1.1w (Figure 9a), which is also due in part to the impact of velocity on reconfiguration. By contrast, the impact of shoot distribution on wake length was shown to produce a smaller uncertainty of 0.92w (Figure 8). This highlights the relatively large impact of shoot reconfiguration on wake length, compared to the within-patch structure. Furthermore, it demonstrates the importance of intraspecific variability in flexural rigidity in determining wake length.

CFD-Biomechanical Model: Strengths, Limitations, and Future Directions
The coupled CFD-biomechanical model provides a new methodology for investigating patch reconfiguration and flow-vegetation-sediment interaction, using high-resolution TLS and select field measurements to parameterize a site-specific model. The model is not currently designed to represent the following: (i) the complexity of stem and leaf scale flow-vegetation interactions, (ii) interpatch variation in geometric and biomechanical properties, or (iii) specific patches within a site. Instead, it represents a site-averaged vegetation patch model that enables consideration of complex shoot distribution and flexural rigidity impacts which are either lacking in current models or restricted to single experimental cases (Tempest et al., 2015).
Due to the number of simulations and associated computational time constraints, the simulations were run as steady state RANS solutions. As such, while the model was able to capture static reconfiguration of the vegetation to the time-averaged flow conditions, it could not capture the impacts of time-varying plant motion. In reality, shoots undergo dynamic reconfiguration (Siniscalchi & Nikora, 2013) in response to the turbulent flow which may act to either further dampen turbulence within the wake (Ghisalberti & Nepf, 2006) or introduce higher levels of turbulent kinetic energy and turbulent events which may decrease deposition (Ghisalberti & Nepf, 2009;Marjoribanks, Hardy, Lane, & Parsons et al., 2017).
Furthermore, Mont Saint-Michel is a dynamic environment which undergoes periodic tidal-driven sedimentation events. Therefore, the observed sediment deposition is the product of multiple events and wake deposition will in turn affect the flow field and therefore there may exist positive or negative feedbacks on sediment deposition. Similarly, shoot-scale turbulence-induced resuspension and scour around the edges of the patches may influence the wake flow and sediment processes (Tinoco & Coco, 2018;Yang & Nepf, 2018). Flow-sediment interactions will also drive, and be driven by, longer-term plant dynamics. Individual plant life cycles, changes to patch structure, and colonization processes mean that patch density, porosity, and flexural rigidity, three important variables within the model, will all vary over time. Moreover, patches do not exist in isolation and therefore impacts of wave sheltering and larger-scale configuration will affect flow velocities and lead to trait response.

10.1029/2019JF005012
Journal of Geophysical Research: Earth Surface The current CFD framework provides a mechanism for extending the capability of the model to capture these complexities through time-dependent solution using large eddy simulation, investigation of cumulative spatially distributed erosion and deposition over multiple tide cycles, and through multipatch simulation. Such developments would allow characterization of three-dimensional wake structure, permit investigation of the role of large-scale coherent turbulent structures on plant reconfiguration and sediment deposition, and enable investigation of wave damping and dissipation processes, which have been shown to be strongly dependent on flexural rigidity (Möller et al., 2014;Rupprecht et al., 2015). Nevertheless, the present model, with its representation of shoot reconfiguration and shoot-scale complexity, provides new insight into the impact of patch geometry, shoot distribution, and shoot biomechanics on wake length.

Conclusions
This study introduces a new field-parameterized vegetation patch model which permits the investigation of shoot-scale flow-vegetation interactions, bridging the gap between existing single-stem models and largerscale rigid-cylinder or roughness approaches. The model was applied to investigate flow and potential sediment dynamics around pioneer saltmarsh vegetation patches. The main conclusions of this study are the following: 1. For the Suaeda maritima saltmarsh vegetation patches considered here, the vertical shear layer is dominant in controlling wake length, in contrast to much of the previous research on vegetation patches where width dominance is assumed. 2. Accordingly, flexural rigidity and patch reconfiguration become key controls on wake deposition length, in addition to patch height and width. Therefore, quantifying the intraspecific variability of biomechanical characteristics is critical to understanding the impact of vegetation on saltmarsh accretion. 3. Adapting the model of Zong and Nepf (2012) for height-controlled, flexible vegetation patches produces a simplified model for wake length which matches well with simulation data.