Non-breaking Wave Effects on Buoyant Particle Distributions

The dispersal of buoyant particles in the ocean mixed layer is influenced by a variety of physical factors including wind, waves, and turbulence. Microplastics observations are often made at the free surface, which is strongly forced by surface gravity waves. Many studies have used numerical simulations to examine how turbulence and wave effects (e.g., breaking waves, Langmuir circulation) control buoyant particle dispersal at the ocean surface. However these simulations are not wave phase-resolving. Therefore, the effects of an unsteady free surface due to surface gravity waves remain unknown in this context. To address this, we develop an analytical model for the distribution of buoyant particles as a function of wave-phase under wind-wave conditions in deep-water. Using this analytical model and complementary numerical simulations, we quantify the effects of a nonbreaking, monochromatic, progressive wave train on the equilibrium vertical and horizontal distributions of buoyant particles. We find that waves result in non-uniform horizontal distributions of particles with more particles under the wave crests than the troughs. We also find that the waves can stretch or compress the equilibrium vertical distribution. Finally, we consider the effects of waves on the sampling of microplastics with a towed net, and we show that waves have the ability to lower the measured concentrations relative to nets sampling without the influence of waves.


INTRODUCTION
Plastic pollution in the ocean is ubiquitous. It is found across ocean basins, in the deep-sea, and along our coasts. The estimated sources of plastic into the ocean exceed the current amount of plastic measured in the ocean (Jambeck et al., 2015;Geyer et al., 2017); this discrepancy has become a subject of interest in the scientific community. Solving this question of "missing plastic" requires both accurate measurements and a thorough understanding of oceanic microplastics sources and sinks.
Because a large fraction of plastic in the ocean is buoyant, microplastics measurements are frequently conducted at the surface of the ocean. Direct observations are representative of the microplastic concentration at the time of sampling. However, concentrations are a function of the local conditions: on a calm day, the buoyant plastic can aggregate at the surface, but strong winds can disperse plastic throughout the mixed layer. To address this variability, large-eddy simulations of buoyant particles at the ocean surface have been conducted to show how microplastics can redistribute due to turbulent mixing by wind, currents, breaking waves, and Langmuir circulation Kukulka and Brunner, 2015;Liang et al., 2018). While these studies have provided valuable insight into how particles distribute throughout the ocean mixed layer, they are not wave-phase or free surface resolving; rather, they rely on parameterizations to include wave-averaged effects (see Chamecki et al., 2019 for a review). In reality, observations at the surface of the ocean are taken by a sampling apparatus that experiences individual wave events and an unsteady free surface. If waves induce phase variability in microplastics concentrations, and waves are not sampled uniformly, bias in the observations can result. Therefore, by studying buoyant particle distributions without wave-averaging, we are able to see important phenomena that may affect how we interpret data from free-surface net tow sampling.
A large body of work has studied the numerous ways in which waves affect surface transport and mixing processes. Waves induce a net drift in the direction of the waves, referred to as Stokes drift (see van den Bremer and Breivik, 2018 for a review); this drift is known to be important for microplastics transport from coastal shelf to oceanic scales (Isobe et al., 2014;Onink et al., 2019). Furthermore, the properties of the particles, such as inertia and shape, have also been shown to be relevant to their transport in waves (Eames, 2008;DiBenedetto et al., 2018). The ability of a wave field to disperse horizontal tracer particles has been studied as well (Herterich and Hasselmann, 1982). With regard to buoyant particles in waves, simulations similar to the ones reported here were conducted by Boufadel et al. (2006), who demonstrated how waves affect oil droplet plumes. Waves have also been studied in relation to ocean mixing through turbulence induced by wave orbital motion (Babanin, 2006). However, the effects of waves on the equilibrium distributions of buoyant particles has to our knowledge, not yet been reported.
In addition to hydrodynamic conditions, observations can be affected by the microplastic sampling method, e.g., net tows, grab samples, or filtering pumps (Hidalgo-Ruz et al., 2012;Prata et al., 2019). Some work has been done to address variability in measurements across methods (e.g., Setälä et al., 2016;Barrows et al., 2017;Karlsson et al., 2019); these studies have demonstrated that net tows can both under-sample and over-sample relative to other measurement techniques. While net tows are ideally conducted under calm conditions (e.g., sea states with Beaufort scales of 0-2) (Viršek et al., 2016), microplastics observations have been reported from conditions with Beaufort scales up to 5 (Eriksen et al., 2014;Kooi et al., 2016;Poulain et al., 2019). Therefore, moderate waves are often encountered during net tows.
Neuston (surface) net tows are conducted by trawling a net over a transect at the surface of the ocean. Captured particles are identified and counted in order to monitor microplastics surface concentrations. As sampling methods evolve and improve, understanding the variability and uncertainty in historic datasets is still necessary when analyzing temporal trends of microplastic pollution. Many longterm datasets are from observations taken with surface net tows (Law, 2010;Eriksen et al., 2014;van Sebille et al., 2015). Therefore, this study's discussion focuses on the influence of waves on the net tow sampling method.
This work describes and quantifies the effects of nonbreaking surface gravity waves on both the vertical and horizontal distributions of buoyant particles at the free surface, and relates these results to microplastics sampling. In section 2, current models used for interpreting microplastics measurements at the free surface are reported. Next, in section 3.1, an analytical description for the distribution of particles at the free surface including wave effects is developed, and in section 3.2, the numerical simulation methodology is described. The results of the models are shown in section 4. These results are considered in the context of microplastics sampling efforts in section 5. Finally, the article finishes with a discussion and summary of the major findings in sections 6 and 7. The purpose of this work is to inform current and past microplastics observations, as well as to provide context for the design of new microplastic sampling procedures. This study is only a first step in fully understanding how to interpret microplastics data, and wave-resolving direct numerical simulations and controlled laboratory experiments are recommended as future steps.

BACKGROUND
The ocean surface concentration of microplastics is both a function of the baseline average contamination level and the instantaneous hydrodynamics. Transient local forcing, such as wind and waves, introduce temporal and spatial variability to microplastics measurements. This variance disrupts the ability to directly compare observations. To account for this variance, specialized vertical mixing models have been developed; these models allow surface measurements to be extrapolated to total water column concentrations based on knowledge of local conditions. This process most often uses a one-dimensional model which we will refer to as the wind-mixing model as presented by Kukulka et al. (2012).
The wind-mixing model assumes a balance between the upward buoyancy flux and the downward turbulent flux of particles (Kukulka et al., 2012). This is represented as a horizontally-averaged concentration of microplastics c that is a function of vertical position z. At equilibrium, the total flux of particles is equal to zero. Following Kukulka et al. (2012), we assume an eddy dispersivity model for the turbulent flux, where κ t is the eddy dispersivity of the particles and w p is the particle's rise velocity. Setting the buoyancy flux equal to the turbulence flux of particles gives: The ratio of the eddy viscosity to the eddy dispersivity, ν t /κ t , also known as the turbulent Schmidt number Sc t , is assumed to be unity. Further assuming that the particle rise velocity w p and the turbulent eddy dispersivity κ t are constant, Equation 1 can be directly integrated. Thus, the mean concentration of particles as a function of z follows an exponential curve with the surface concentration c 0 at the surface (z = 0): The decay length scale, or mixing scale L m , is defined as L m = κ t /w p . This scale represents the relative depth over which particles are dispersed. Kukulka et al. (2012) assumes that the depth of the mixed layer is much larger than L m for the assumption of a constant κ t to be valid. Figure 1 shows the relationship described by Equation (2). One common use of this model is to extrapolate surface measurements to total water column concentrations (e.g., Cozar et al., 2014;Kooi et al., 2016;Lebreton et al., 2018;Poulain et al., 2019). This extrapolation can be done by reconfiguring the model as a correction factor f defined as where N is the predicted total count of microplastics in the water column, and N tow is the number of microplastics captured in the towed net (Kukulka et al., 2012). The parameter z 0 represents the depth over which plastic was measured, or the average submergence of the towed net.

METHODS
The wind-mixing model is the current standard used to correct surface microplastic measurements, so we start with the assumptions within this model to examine the effects of surface gravity waves. This model includes only uniform turbulence and particle buoyancy. The assumption of constant eddy viscosity is only valid very close to the surface.
In addition, this model assumes that particle inertia effects are negligible. The relative importance of particle inertia can be defined by the particle's Stokes number St, which is a ratio between the particle relaxation timescale and the relevant flow timescale. A small Stokes number implies that the particle behaves as a flow tracer with little to no particle inertia effects. In this case, we assume the particle relaxation time τ p is defined as FIGURE 2 | Example of the wave geometry and the wave induced velocity field shown as scaled vectors. The wave is traveling in the postive x direction, with wavelength λ and wave amplitude A. The free surface is denoted by η, and ζ is the coordinate representing distance from the instantaneous free surface. τ p = w p /g, and the relevant flow time scales are the wave period T and the Kolmogorov timescale of the turbulence τ η . Thus, we can construct a wave Stokes number St w and a turbulence Stokes number St η . A common measured rise velocity for microplastics is w p = 1 cm/s (Kukulka et al., 2012), which corresponds to buoyant microplastics of about 1 cm according to Figure 2 in Poulain et al. (2019). Using τ p = w p /g, we find τ p ≈ 0.001 s. Wind wave periods are on the order of seconds, so St w = O(10 −3 ) ≪ 1 and the assumption of no particle inertia is valid. In terms of the turbulence, a conservative estimate using a high dissipation rate at the ocean surface boundary layer of 10 −3 m 2 /s 3 corresponds to τ η ≈ 0.03 s (Zippel et al., 2018). In this case, St η = O(10 −2 ) ≪ 1 as well. While the assumption of negligible particle inertia theoretically holds in this case for microplastic particles up to approximately 1 cm, finite-size effects may be important for particles at this size and should be studied further.
In order to isolate the wave effects, we consider a progressive train of surface gravity waves in deep-water. We start with linear wave theory which assumes that the wave-induced flow is irrotational, incompressible, and 2-dimensional. It also assumes small amplitude waves with a free surface deflection η, where η = A cos θ and θ = kx − ωt is the wave phase. In this system t represents time, x is the horizontal coordinate, and z is the vertical coordinate where z = 0 is the free surface. We further assume that the waves are traveling over an infinite depth H such that kH ≫ 1, referred to as the deep-water wave limit. The waves are controlled by the wave number k, related to the wave length λ where k = 2π/λ; the wave frequency ω, related to the wave period T where ω = 2π/T; and the wave amplitude A. These parameters are not independent as the deep-water dispersion relationship states ω = gk where g is the acceleration due to gravity. We will refer to the vertical coordinate relative to the instantaneous free surface as ζ = z − η. A schematic of the wave's geometry, its induced velocity field, and the coordinate system is shown in Figure 2.
In this model, without waves, the particles are mixed solely due to turbulence and vertically transported due to buoyancy, and the only length scale in the problem is the mixing scale L m . With the addition of waves, new length scales are introduced. The waves are controlled by ǫ = kA which is the relative wave steepness. Linear wave theory assumes ǫ ≪ 1, and physically ǫ < 0.44 or the waves will break (Stokes, 1880;Perlin et al., 2013). To relate the mixing to the waves, we use the ratios A/L m and kL m .
The exact values of A/L m and ǫ in the ocean will vary considerably. Kukulka et al. (2012) considered a range of L m from about 0.4-4 m. The most extreme sea state that samples are conducted under is Beaufort 5 which can have wave amplitudes on the order of 1.5 m. Therefore, we report findings using A/L m values of 1/3, 1, and 3 across a range of ǫ from 0 to 0.4 to demonstrate the variability of scenarios potentially encountered during sampling.

Wave-Averaged Distribution
Surface gravity waves induce wave-orbital motions and deflect the free surface, both of which can affect buoyant particle distributions. We first address the distortion of the free surface; we model this effect by considering a superposition of the windmixing model shifted to the instantaneous free surface and averaged over a wave period. We assume a linear wave η = A cos θ , and we denote the phase-dependent concentration as c, where in this case We average this concentration over a wave-period, where c is a wave-period averaged quantity, This expression is evaluated exactly for all points under the wave trough, and above the wave trough the integration range is a function of z/A, which can be integrated numerically. The full expression is: where I 0 is an elliptic integral of the first kind. From inspection, we see that the wave-averaged vertical distribution c (z) is a function of the original solution without waves, c(z), and a factor dependent on the ratio A/L m . As A/L m → 0, the wave solution approaches the no-wave solution, or c (z) → c(z).

Wave-Phase Variability
Waves also alter the vertical and horizontal particle distributions as a function of wave phase. These effects are largely due to the wave-orbital motion of the particles. We begin our analysis with a study of how the wave kinematics affect the horizontal distribution of particles. At one point in time, horizontal position maps onto wave phase. Thus, we use wave phase distribution as a proxy for the horizontal distribution. We approximate wave phase as a function of time by linearizing the particle equations of motion. Under linear wave theory in deep-water, the horizontal and vertical wave-orbital velocities, respectively, are: Next, we approximate the horizontal particle motion as ξ wherė ξ = u w (x 0 , z, t), and x 0 is a constant value that refers to the center of the particle's horizontal motion. By integration, ξ (t) = −A exp kz sin(kx 0 − ωt). Assuming x ≈ ξ (t), we plug this into our expression for θ (t) leaving: To obtain the probability distribution function (p.d.f.) of θ over t, we use a Jacobian transformation where p(θ ) = |dt/dθ |p(t) for a uniform probability of t over one wave period. We then further approximate θ ≈ kx 0 − ωt to obtain a closed form solution. At this level of approximation, we find the p.d.f. of θ given a value of z to be We account for the phase variations over the vertical distribution of particles by weighting Equation 9 with the no-wave vertical particle distribution from Equation (2), vertically integrating from z = −∞ to z = 0, and normalizing the resultant equation: This integration is solved numerically to find the verticallyaveraged wave-phase p.d.f. of buoyant particles. Note that the p.d.f. of particles over wave phase is controlled by the parameters k, L m and ǫ. For large kL m or small ǫ, the p.d.f. becomes uniform over a wave. While Equation 10 describes how waves affect the total amount of particles under each wave phase, we also find that waves can change the shape of the distributions under a wave. They do this through vertical stretching and compression which has been studied previously in relation to vortical structures and turbulence (Phillips, 1961;Guo and Shen, 2013). The stretching of the distribution is controlled by the vertical divergence of the flow: dw/dz = ωǫ exp kz sin θ . To quantify these effects, we calculate L m as a function of wave phase and account for the vertical stretching with the following relationship, approximating the stretching effect by integrating in time: After integrating and evaluating at z = η − L m , we are left with In this case, L m is measured from the instantaneous free surface. The waves' effects on L m are strongest for large ǫ and small kL m .
Frontiers in Marine Science | www.frontiersin.org We synthesize our results from Equations (10) and (12) and find an expression for the vertical distribution of particles under waves and turbulence: This equation is similar to Equation (4), but it now accounts both for the wave-orbital motion effects as well as the effects of the unsteady free surface. It is also a function of ζ , which is the vertical coordinate relative to the instantaneous free surface. In order to both validate this model and expand on its theory, we also conducted complementary numerical simulations.

Numerical Simulations
To complement the analytical model, we use simple numerical simulations with Lagrangian particle tracking. The turbulent dispersion is modeled using a random walk, which allows us to track individual particles. These Lagrangian simulations contrast the Eulerian analytical model of the concentration field. They also allow for a more accurate representation of the wave field without the approximations made in the derivations.
In the simulations, the wave field is imposed using the analytical form of Stokes third order deep-water waves. Because we are interested in particle behavior near the free surface, we use Stokes third order waves to minimize the errors between the trajectories of particles at the free surface, and the free surface itself. This error is O(ε 2 ) for linear wave theory, but O(ε 4 ) for Stokes waves. At this level of approximation, the free surface η is defined as with wave phase θ = kx − ω ′ t where ω ′ = (1 + 1 2 ǫ 2 )ω. The waveinduced velocity field described in Equation (7) is extrapolated to the wave crests due to the agreement seen in Baldock et al. (1996) and as summarized in Smit et al. (2017).
The particle trajectories are simulated using a fourth order Runge Kutta integration which calculates particle velocity based on the superposition of the analytical wave velocity field (see Equation 7) and a constant rise velocity w p . After the velocity is calculated, a random walk perturbation in the vertical and horizontal directions is applied based on the time-step t and κ t : with normally distributed random variables X, Y ∼ N (0, 1). Particles were constrained to not go above the free surface by setting z = η for all particles with z > η after each integration step using Equation (14). The ratio A/L m and ǫ were varied for a total of 42 unique scenarios. The values of A/L m were varied from 1/10 to 10, and the values of ǫ were varied from 0 to 0.3. Due to similarities among the trials, we only report results from a subset of the runs. The simulations were run in MATLAB and started with 1,000 particles at the surface distributed evenly over one wave period, defined by the range {x | 0 ≤ x ≤ 2π/k}, and run until the vertical distribution and wave phase distribution reached equilibrium. Once at equilibrium, the simulations were run for at least an additional 1,000 wave cycles, at which point the data was sampled periodically resulting in 1,000,000 data points for each run.

MODEL RESULTS
Both the numerical simulations and the analytical model demonstrate how the unsteady free surface associated with surface gravity waves fundamentally alters the equilibrium vertical and horizontal distributions of buoyant particles. This can be seen in Figure 3, where we plot c(z), the wave-period averaged concentration in the fixed coordinate system for different A/L m values from both the numerical simulations and the analytical model. For each case, the peak concentration is at the level of the wave trough (when z/L m = −A/L m ), and below this depth the distribution resembles an exponentially decaying curve. However, above this depth, the distribution deviates from the exponential curve. This is due to the fact that there is an unsteady free surface wetting and drying, and thus there are fewer total particles on average when compared to the case without waves. As A/L m increases, the distribution above the trough level becomes more uniform, and as A/L m decreases, the distribution approaches the case without waves.
The wave-orbital motion of the particles also introduces phase variation in the distributions. In Figure 4 we plot the phaseaveraged p.d.f. of particles from both the numerical simulations and the theory using Equation 10. In this figure, we see that on average, particle concentration is enhanced under the wave crests (θ = 0, 2π) and reduced below the troughs (θ = π). For the steepest waves, the number of particles under the wave crest can be up to 50% larger than under the wave trough. The total concentration of particles is phase-dependent, and so is the shape of their vertical distribution. In Figure 5, we plot the concentration of particles over depth, relative to the instantaneous free surface. The data are from numerical simulations with ǫ = 0.2 and A/L m = 3, and are subdivided into the concentration under the wave trough, wave crest, the average profile over the wave, and the no-waves case for reference. As expected, the wave crest has higher concentrations than the wave trough. The largest difference in concentrations is near the surface, and the concentration profiles approach each other at depth. The surface concentration at the crest can be up to 50% larger than the average surface concentration across the wave. We also see that the average vertical profile, when adjusted to the instantaneous free surface, recovers a similar profile to the no-wave case.
Another way to visualize the phase-dependence of the particle distributions is with a 2-dimensional p.d.f. of particles over wave-phase and relative depth as plotted in Figure 6. In this figure, we again see that the peak concentration at the surface is largest under the wave crests. To characterize the shape of the distribution, we calculate the exponential decay length scale, or mixing length scale L m , which is plotted as a solid line from the numerical simulations and a dashed line using Equation (12). The two solutions agree with each other. The value of L m is smaller under the troughs than it is under the crests, indicating that the vertical distribution is compressed under the troughs and stretched under the crests.
The results from both the analytical model and the numerical simulations demonstrate that waves can affect the distribution of buoyant particles near the free surface. The benefit of the numerical simulations is that they are more accurate, however FIGURE 5 | Vertical distribution of relative concentration of buoyant particles under waves as a function of wave phase from numerical simulations. Data comes from a simulation with A/L m = 3 and ǫ = 0.2. The wave phase was subdivided where the crest is the region defined by θ < 0.04π or θ > 2π − 0.04π and the trough is defined by π − 0.04π < θ < π + 0.04π. The no waves line shows the distribution from Equation (2), and the average line shows the wave-period averaged distribution. they have a higher computational cost associated with them. Nevertheless, the numerical simulation results largely agree with the theory described in section 3.1, but they differ in a few distinct ways. The waves considered in the theory are from linear wave theory, but the numerical simulations use higher order harmonics. These non-linear wave effects result in waves with taller crests and shallower troughs, the effects of which can be seen in Figure 3 where the microplastics reach higher in the water column in the numerical simulations when compared with the analytical model. The surface boundary condition combined with the discrete random walk used in the simulations also results in a concentration discontinuity at the free surface; this is a small effect only seen when analyzing the results very close to the free surface, which is done in the following section. However, the close agreement between the two methods show that the assumptions behind our analytical model capture the major effects of waves in these scenarios.

SAMPLING IMPLICATIONS
In this section we consider how the waves may affect the sampling of microplastics at the ocean surface. Under the wind-mixing model, a constant free surface profile at z = 0 is assumed. If the free surface is deflected due to waves, the total surface area between two fixed points increases. A net trawled over a wavy free surface will not sample around z = 0, but rather around z = η. The waves may also cause the net to sample less effectively; if the net lags the waves, it will not exactly track the free surface with a constant sample volume and the assumption of a single z 0 value may break down. We consider how these effects may alter the particles sampled under various scenarios.
If we let δ be the submerged height of the net opening (where δ corresponds to a varying z 0 ), W be the width of the net opening, and d be the distance the net traveled over the free surface, then the total sampled volume of water is V = dδW. Let the total number of particles counted in the net be N tow . Thus, the measured concentration is c tow = N tow /V. While W should be constant, d is a function of the free surface deflection, and δ can vary due to the net moving relative to the water surface. We apply this framework to our wave-phase distribution results to discern the possible effects of waves on sampling.

Surface Tow Length
A net following a deflected free surface will sample a larger volume than if trawled over a flat free surface. The total distance trawled is a function of the free surface deflections, where d from x = x a to x = x b is defined using the arc-length formula as Over one wave length λ, we integrate Equation (16) for a linear wave which results in an expression using an elliptic integral of the second kind E(m): The effective increase in the length of free surface due to the wave can be defined by the ratio d/λ which is only a function of the wave steepness ǫ: For a very steep wave where ǫ = 0.3, we find d/λ = 1.022 for a linear wave and 1.024 for a 3rd order Stokes wave. These are small differences in total tow length that are not necessarily larger than other sources of error during the sampling processes. However, this analysis assumes that as the sample volume travels one linear distance λ, it samples exactly one wave. In reality, the sample volume and the wave are both traveling. Under a traveling wave, d is a function of the relative speed between the sampler and the wave. For simplicity, we will assume that the wave and the sampler are traveling in the same plane. (This is often the case, as sampling vessels will align into the wind direction in order to minimize net destabilization from cross-winds.) If the sampler travels at speed U, then the free surface at the sampler as a function of time is This relationship yields a new effective wave steepness ǫ ′ : where c p = ω/k is the wave phase speed. Larger values of ǫ result in larger values of d/λ. Thus, the expression for ǫ ′ suggests that the effective wave steepness is very sensitive to the relative speed of the boat and the waves. Even moderate waves traveling against the direction of the boat can dramatically increase the effective surface area trawled. For example, a boat traveling at 2 knots in the opposite direction of a 50 cm, 5 s deep-water wave where ǫ = 0.08 has ǫ ′ = 0.69 and d/λ = 1.11. The wave-enhanced tow length effect is relevant to the way in which concentrations are calculated from net trawls. The volume of water sampled is sometimes measured directly using a flow meter attached to the net (Isobe et al., 2014;Kooi et al., 2016), but it can also be estimated by using the linear distance of the trawl (often calculated using GPS coordinates) (Law, 2010;Lebreton et al., 2018). The former method may better capture the free surface deflections due to waves as compared to the latter. Let the concentration calculated using the total arclength of the free surface be c tow , and c 0 tow be the concentration calculated using only the linear distance traveled. To assess the sensitivity of these two concentration measurement techniques to waves, c tow /c 0 tow is plotted in Figure 7 as a function of c p /U and ǫ. In this case, c tow /c 0 tow = λ/d assuming constant N tow , δ, and W values. As seen in the figure, waves can reduce the measured concentration by up to 50% in the most extreme cases. This difference could potentially account for some of the discrepancies across different measurement techniques in wavy conditions.

Non-uniform Sampling Effects
The ideal trawling scenario has a net that samples at a constant depth below the free surface. If the free surface is unsteady, then the depth of the trawl relative to mean water level is η−δ. In many cases, the net does not necessarily sample at a constant depth. We consider the effects of non-uniform net sampling by assuming that the net lags the free surface resulting in an oscillation with a smaller amplitude: where α is the amplitude reduction factor. Ideal sampling results in a constant δ value with α = 1. A sampler that never moves vertically and constantly samples around z = 0 corresponds to α = 0. The sample volume over a wave as a function of α is shown schematically in Figure 8. To assess the effects for α = 1, we also assume that the net never leaves the water surface such that δ > (1 − α)A. Finally, we only consider α ≤ 1 where the troughs are sampled less than the crests. We argue that this corresponds to the more relevant scenario where the net lags the free surface, rather than leads it.
We apply Equation (21) to the numerical simulation results and the distribution modeled by Equation 13 to evaluate how a non-uniform sampler affects the measured number of particles.
We assume the net is twice the height of δ; this means on-average the net is halfway submerged. Over a wave cycle, the particles above δ and below the top of the net are integrated to determine the effective N tow . The value N tow is then normalized by the reference N 0 tow which is calculated using Equation (3). The value of N 0 tow represents the number of measured particles using the wind-mixing model with no waves. Assuming a constant volume sampled in both scenarios, the resultant c tow /c 0 tow values are plotted against α for various wave scenarios in Figure 9. We also consider two different net sizes by varying the ratio δ/A.
In Figure 9, we see that as α decreases, the amount of particles sampled decreases relative to the no-wave prediction. The effect is the strongest for large A/L m . The concentrations show the most decrease after the net starts to dip below the free surface for part of the wave cycle, denoted by the vertical lines. In Figure 9A, δ/A = 0.2, and in Figure 9B, δ/A = 0.5; by comparing these figures we see that the larger the net is relative to the wave amplitude, the less sensitivity there is to α. However, once the net goes under water in both scenarios, we see a decrease in sampled concentrations. An increase in c tow /c 0 tow is seen in the numerical simulations for α values close to 1 in Figure 9A. This is due to a slight enhancement of particles at the free surface in the simulations as a result of the surface boundary condition.

DISCUSSION
Using a phase-resolving approach, we have analyzed the effects of waves on the equilibrium distribution of buoyant particles in the ocean mixed layer, demonstrating that waves introduce concentration differences across a wave cycle. This work implies that if a towed net does not sample each wave phase equally, the resulting observations will be biased.
The phase variability in the distributions is controlled by the wave kinematics. A tracer particle under a linear wave field spends more time under the crests than the trough. One way to conceptualize this phenomena is by considering a particle traveling under a trough and under a crest. The particle under the trough moves in the opposite direction of the wave, so it is leaving the trough more quickly than if the wave was not traveling. However, under the wave crest, the particle is moving in the direction of the wave, so it has a longer residence time in its relative wave phase. In other words, dθ/dt has a higher magnitude under the troughs than under the crests. This translates to particles on-average, non-uniformly distributed in wave-phase, with more particles under the crests than the troughs.
Another major focus of this study is how waves stretch and compress the particle distributions as a function of wave phase. While these processes have been studied in relation to turbulence and vortices, this work centers on particle distributions. The wave-induced velocities have non-zero vertical divergence. This results in the distribution of particles being maximally stretched under the crests, i.e., having the maximum L m value, and being maximally compressed under the troughs, i.e., having the minimum L m value. Because the time-averaged vertical divergence in the wave field is zero, there is no net effect, and this stretching is only seen as a function of wave phase. These findings demonstrate that particles are distributed the deepest relative to the free surface under the crests, whereas they are closest to the surface under the troughs. This result is illustrated in Figure 6.
Because particles preferentially accumulate under the crests relative to the troughs, if a towed net over-samples the wave crests, it can over-sample the particles. However if the net tends to go under water at the crests, it can under-sample the particles because it will miss the peak concentration at the surface. Another common assumption during net tow sampling is that a net that is on-average halfway submerged will sample the same particles as a net that has constant halfway submergence. A net that comes out of water half the time and submerges underwater half the time, is on-average halfway submerged and corresponds to small α values; we have shown that small α values consistently lead to the under-sampling of particles. Therefore, if waves cause the towed net's submergence to vary over time, under-sampling can result. These effects should be studied more carefully during the time of sampling to further discern any bias in the measurements.
One can imagine there is a range of frequencies to which the sampling net can respond. At the low frequency limit for long waves, i.e., waves which are much longer than the net, the net should be able to effectively track the free surface. In that case, the wind-mixing model is still an appropriate measure to use. However, if the net tends to sample at a constant z-level and move through the waves, which may occur with short-period waves, i.e., waves with higher frequencies than the net's natural response frequency, then the modified wave-period-averaged wind-mixing model presented here may be more appropriate. In that scenario, the net will also under-sample the assumed concentration estimated by the wind-mixing model as seen in Figure 3.
When compared to other sampling methods, there is some evidence that net tows can under-sample. For example, grab samples of microplastics measurements have been reported with high values relative to the net tows at the same location (Barrows et al., 2017). This effect is also due to the ability for grab samples to capture particles smaller than the net's mesh size. However, this does not necessarily explain all the discrepancy in observed concentrations, and non-uniform sampling due to waves could be partly responsible. An advantage of net tows in this scenario is that they encompass average concentrations over large distances which lowers sensitivity and variability to patchiness, whereas grab samples are very sensitive to patchiness. Microplastics patchiness can occur due to Langmuir circulation and other hydrodynamic features (van Sebille et al., 2020).
We find that the accumulation of particles under the wave crests is most extreme under the steepest waves. Therefore, as waves steepen before they break, they will accumulate particles under them. The turbulence due to breaking waves is also strongest under the crests of waves (Rapp and Melville, 1990), and thus enhanced turbulent mixing will occur at the point with the highest concentration of particles at the time of breaking. This coupling could result in greater transport of particles away from the surface and should be studied further.
The ocean surface boundary layer is subjected to a multitude of unsteady forcings. Understanding how they each can bias measurements is important to accurately interpreting data, and this study represents only one mechanism: non-breaking surface gravity waves. In this work, we have only considered a monochromatic wave train aligned with the sampler, and therefore future work should include a wave field with a broadband spectrum as well as scenarios in which waves travel obliquely to the sampling apparatus. We have also neglected the effects of Langmuir turbulence, which has been shown to be extremely important in both controlling patchiness and transporting particles to depth. Finally, understanding the limits of the wind-mixing model is important. When L m is greater than the mixed layer depth, the model's underlying assumption of a constant eddy viscosity no longer holds (Kukulka et al., 2012). Large L m values can result from small particle rise velocities, and thus this model may not be appropriate for the smallest microplastic particles. A more general model that can account for the smallest particles as well as the effects of waves is therefore needed.

CONCLUSIONS
Using both numerical simulations and an analytical model, we have demonstrated how waves affect the vertical and horizontal distributions of buoyant particles at the free surface. In the fixed frame, waves change the wave-period-averaged vertical distribution of particles by lowering the peak concentration of particles to the level of the wave trough. Over a wave cycle, the wave orbital motion redistributes the particles so that on average, more particles are under the wave crest than under the wave trough. Waves also stretch the vertical distribution resulting in a smaller mixing length under the wave trough and a larger mixing length under the wave crest. The free surface is stretched due to waves, and when a sampling volume is traveling relative to the waves, the effective length of the trawl can increase dramatically. Finally, we find that non-uniform sampling of a wave can result in reduced measured particle concentrations relative to predictions based on the wind-mixing model (Kukulka et al., 2012).
While we have shown how a train of monochromatic progressive waves distort buoyant particle distributions, the real ocean surface is subjected to a spectrum of waves. The turbulent mixing at the ocean surface is not necessarily constant nor isotropic. It is also not independent from the waves. Particle inertia can also be important for large particles and low turbulence forcing. However, this work shows the possible range of effects of surface gravity waves on buoyant particle distributions. We recommend a more detailed study of how the trawling process can bias the concentrations observed, and how our findings change under more complicated scenarios.

DATA AVAILABILITY STATEMENT
The datasets generated for this study are available on request to the corresponding author.