Surface Waves Enhance Particle Dispersion

Abstract

We study the horizontal dispersion of passive tracer particles on the free surface of gravity waves in deep water. For random linear waves with the JONSWAP spectrum, the Lagrangian particle trajectories are computed using an exact nonlinear model known as the John–Sclavounos equation. We show that the single-particle dispersion exhibits an unusual super-diffusive behavior. In particular, for large times t, the variance of the tracer $⟨|X(t){|}^2⟩$ increases as a quadratic function of time, i.e., ${⟨\left|X(t)\right|}^2⟩\sim t^2$. This dispersion is markedly faster than Taylor’s single-particle dispersion theory which predicts that the variance of passive tracers grows linearly with time for large t. Our results imply that the wave motion significantly enhances the dispersion of fluid particles. We show that this super-diffusive behavior is a result of the long-term correlation of the Lagrangian velocities of fluid parcels on the free surface.

1. Introduction

Water waves cause the material transport of fluid particles on the free surface of the fluid. The waves induce a fluid velocity on the free surface which in turn determines the horizontal motion of fluid particles on the free surface. This phenomena has been known since Stokes [1] who studied the average velocity of fluid parcels transported by a linear monochromatic wave. The resulting displacement is referred to as the Stokes drift.

More specifically, denote the horizontal trajectory of a fluid particle on the free surface by $x(t;t_0,x_0)$ (or $x(t)$, for short). For simplicity, we assume that the waves are unidirectional. The map $x(t;t_0,x_0)$ denotes the horizontal position of a fluid parcel at time t given its initial position $x_0$ at time $t_0$. The surface elevation is assumed to be a graph over the horizontal coordinates and is denoted by $\zeta (x,t)$ (see Figure 1). Since the fluid parcels on the free surface are constrained to it, their vertical position of the parcel is given by $\zeta (x(t),t)$. Therefore, knowledge of the horizontal position $x(t)$ of the fluid particles and the free surface elevation $\zeta$ completely determines the position of the particles.

For a fixed initial condition $x_0$, Stokes’ theory [1] seeks to determine the average displacement $\langle x(t)\rangle$ of the particles at a later time. Here, $\langle ·\rangle$ denotes the ensemble average over many realizations of the random surface waves. In other words, the Stokes’ theory is concerned with the first-order statistics of the fluid displacement on the free surface. The main focus of the present work is the second-order statistics of particle dispersion. More precisely, we examine the temporal evolution of the variance $\langle |X(t){|}^2\rangle$ where $X(t)=x(t)-\langle x(t)\rangle$ denotes the mean-zero displacement of the fluid particles.

1.1. Summary of the Main Results

As we review in Section 1.2, this second-order statistics has also been studied extensively. However, previous studies either assume the velocity field as a linear superposition of the velocity induced by linear waves or use perturbation theory to approximate the velocity field. Here, however, we obtained the fluid particle trajectories $x(t;t_0,x_0)$ from an exact model known as the John–Sclavounos (JS) equation (see Section 3 for a review). This exact nonlinear kinematic model was first derived by John [2] for unidirectional irrotational water waves. Sclavounos [3] generalized the equations to two-dimensional waves and removed the irrotationality assumption. Fedele et al. [4] further analyzed the JS equation discovering its underlying Hamiltonian structure.

Using the JS equation, we compute the temporal evolution of the variance $\langle |X(t){|}^2\rangle$ for fluid particles evolving from a given initial condition $x(t_0)=x_0$. This puts us in the framework of the single-particle dispersion theory. Single-particle dispersion was first studied by Taylor [5] in the context of homogeneous, isotropic turbulence (see Section 2, for a review). Taylor’s dispersion theory predicts that the variance $\langle |X(t){|}^2\rangle$ exhibits a ballistic motion for short times (i.e., ${\langle \left|X(t)\right|}^2\rangle \sim t^2$) and a diffusive motion for large times (i.e., $\langle |X(t){|}^2\rangle \sim t$), so that the variance follows the scaling laws

$${\langle \left|X(t)\right|}^2\rangle \sim \left\{\begin{array}{c}{t}^2,\mathrm{for}\mathrm{small}t,\hfill \\ t,\mathrm{for}\mathrm{large}t.\hfill \end{array}\right.$$

(1)

These scaling laws are often assumed in the study of fluid particles dispersed by surface waves. The main contribution of the present study is to show that Taylor’s dispersion laws (1) may be violated for particles advected by free-surface waves. More precisely, our results can be summarized as follows.

• The variance $\langle |X(t){|}^2\rangle$ of the particles advected by a surface wave follow the scaling laws

  $${\langle \left|X(t)\right|}^2\rangle \sim \left\{\begin{array}{c}{t}^4\mathrm{for}t/T_p<1,\hfill \\ t^2\mathrm{for}t/T_p>1,\hfill \end{array}\right.$$

  (2)

  where $T_p$ denotes the wave period. This is markedly different from the prediction of Taylor’s theory. In particular, the long-term evolution of the variance follows a ballistic motion as opposed to a diffusive motion.
• Central to Taylor’s theory of single-particle dispersion is the autocorrelation function $R(\tau )$ of the Lagrangian fluid velocities (see Section 2, for a precise definition of R). To arrive at the scaling laws (1), Taylor assumes that this autocorrelation function decays to zero rapidly enough so that the integral ${\int }_0^{\infty }R(\tau )\tau$ exists and its value is finite. We show, however, that for particles dispersed by surface waves, the autocorrelation function decays to $1/2$ as $\tau \to \infty$. This observation, in turn, explains the unusual scaling (2).

1.2. Earlier Studies

Previous studies of particle dispersion on the free surface of a fluid can be divided into three general categories: (i) Passive tracers on a flat free surface, (ii) Linear or nonlinear waves with the induced velocity field modeled based on simplifying assumptions, and (iii) Velocity field derived from satellite altimetry data.

The first category is concerned with the Lagrangian motion of fluid particles on the flat free surface of a fluid in a container. Although the three-dimensional flow is incompressible, the velocity field restricted to the surface can be compressible. As a result, the passive tracers can form clusters on a nontrivial subset of the surface. Such no-slip surface flows have been studied experimentally [6,7] and numerically [8,9,10,11] with their main emphasis being on the clustering patterns formed by tracers on the free surface. However, they neglect the effect of waves on the dispersion of the tracers by assuming a flat free surface.

The second category considers the effect of waves on the particle dispersion [12,13,14,15,16,17]. These studies often assume (explicitly or implicitly) that the longterm particle dispersion is diffusive (Taylor’s theory) and aim to approximate the diffusion tensor. However, they use a simplified model for the fluid velocity field on the surface. For instance, Herterich and Hasselmann [12] assume that the Eulerian velocity field is a linear superposition of the velocities obtained from linear wave theory. Weichman and Glazman [13] also assume such a linear wave theory for their study of passive tracer advection (also see Refs. [14,15]).

Bühler and Holmes–Cerfon [16] consider dispersion by random waves in the rotating shallow water framework. They go beyond the linear wave theory by using the wave-mean interaction theory to account for second-order corrections to the linear velocity field. Holmes–Cerfon et al. [17] use a similar approach in the framework of the rotating Boussinesq equation. Although accounting for the second-order nonlinear effects, these studies also assume that the longterm dispersion is diffusive (Taylor’s theory).

The third category derives ocean velocity field from the satellite altimetry data [18,19,20,21,22,23]. These studies are concerned with the large scale mixing in the ocean (on the order of a few kilometers). It is believed that, over such scales, the main contribution to mixing comes from the large ocean eddies with negligible contribution from the wave motion. Nonetheless, the velocity field is derived from the sea surface height measured by altimetry techniques [24]. In order to relate the sea surface height to the fluid velocity field, one makes the so-called quasi-geostrophic assumption, resulting in an approximation of the true velocity field.

In contrast to previous studies, here we used the John–Sclavounos equation to compute the exact Lagrangian trajectory of fluid parcels on the free surface and thereby examine the validity of the underlying assumptions of Taylor’s single-particle dispersion theory.

Before proceeding further, we also refer to the work on pilot-wave hydrodynamics which concerns the motion of droplets bouncing on the surface of a fluid [25,26]. These droplets create a wave when bouncing off the surface which in turn guides the motion of the droplet upon subsequent bounces. Clearly, the pilot-wave phenomena is distinct from the dispersion of fluid parcels that belong to the free surface (as considered here) and the two should not be confused.

1.3. Outline of the Paper

In Section 2, we review Taylor’s single-particle dispersion theory. Section 3 reviews the JS equation for the motion of fluid particles on a surface wave. In Section 4, we describe the set-up of our numerical simulations. Our numerical results are presented in Section 5. Finally, we present our concluding remarks in Section 6.

2. Review of Taylor’s Single-Particle Dispersion Theory

In this section, we briefly review the single-particle dispersion theory of Taylor [5]. This theory is not limited to particle on a free-surface wave and applies more generally to fluid particles advected by turbulent velocity fields which satisfy the simplifying assumptions mentioned below. Fluid particles move according to the ordinary differential equation,

$$\dot{\mathbf{x}}=\mathbf{u}(\mathbf{x},t),\mathbf{x}(t;t_0;{\mathbf{x}}_0)={\mathbf{x}}_0,$$

(3)

where $\mathbf{x}(t;t_0;{\mathbf{x}}_0)\in {\mathbb{R}}^d$ denotes the trajectory of a fluid particle starting from the point ${\mathbf{x}}_0$ at the initial time $t_0$. The time-dependent vector field $\mathbf{u}(\mathbf{x},t)\in {\mathbb{R}}^d$ denotes the Eulerian fluid velocity field.

We also define the Lagrangian velocity $\mathbf{v}(t;t_0,{\mathbf{x}}_0)=\mathbf{u}(\mathbf{x}(t;t_0,{\mathbf{x}}_0),t)$ which measures the fluid velocity along the trajectory $\mathbf{x}(t;t_0,{\mathbf{x}}_0)$. Integrating Equation (3) in time, we obtain

$$\mathbf{x}(t;t_0,{\mathbf{x}}_0)={\mathbf{x}}_0+{\int }_{t_0}^t\mathbf{v}(t^{\prime };t_0,{\mathbf{x}}_0)t^{\prime }.$$

(4)

Taylor [5] views the Lagrangian velocity $\mathbf{v}$ as a stochastic process and seeks to derive the properties of the resulting stochastic process $\mathbf{x}$. In the following, we briefly review Taylor’s argument. For notational convenience, we omit the dependence of the Lagrangian trajectory and velocity on the initial conditions $(t_0,{\mathbf{x}}_0)$ and simply write $\mathbf{x}(t)$ and $\mathbf{v}(t)$.

Taylor [5] also assumes that the flow is homogeneous and isotropic. These assumptions imply that the flow can be considered as a one-dimensional motion ($d=1$). Therefore, we omit the vectorization of the quantities and denote the fluid parcel’s position and velocity by $x(t)$ and $v(t)$, respectively. We also define the mean-zero position and Lagrangian velocity of the particles,

$$X(t):=x(t)-\langle x(t)\rangle ,V(t):=v(t)-\langle v(t)\rangle ,$$

(5)

where $\langle ·\rangle$ denotes the ensemble average. It is straightforward to show that

$$X(t)={\int }_{t_0}^{t}V(t^{\prime })t^{\prime }.$$

(6)

Taylor’s theory predicts the scaling (1) for the variance of the quantity $X(t)$. Note that this prediction implies that, for short times, the particle dispersion is ballistic while, for long times, the particles diffuse as if they are undergoing Brownian motion. In order to derive expression (1), we consider the time derivative of the variance of $X(t)$,

$$\frac{\partial }{\partial t}{\langle \left|X(t)\right|}^2\rangle =2\langle X(t)V(t)\rangle =2{\int }_{t_0}^t\langle V(t)V(t^{\prime })\rangle t^{\prime }=2{\int }_0^t\langle V(t)V(t-\tau )\rangle \tau ,$$

(7)

where for the last identity we used the change of variables $\tau =t-t^{\prime }$ and assumed $t_0=0$ (without loss of generality). Assuming that V is a stationary process, the covariance $\langle V(t)V(t-\tau )\rangle$ only depends on $\tau$. Therefore, the autocorrelation function $R(\tau ):=\langle V(t)V(t-\tau )\rangle /\langle |V(t){|}^2\rangle$ only depends on the delay parameter $\tau$. Integrating Equation (7) in time, we obtain

$$\begin{array}{cc}\hfill \langle |X(t){|}^2\rangle & ={2\langle \left|V(t)\right|}^2\rangle {\int }_0^t{\int }_0^{t^{\prime }}R(\tau )\tau t^{\prime }\hfill \\ & ={2\langle \left|V(t)\right|}^2\rangle {\int }_0^t(t-\tau )R(\tau )\tau .\hfill \end{array}$$

(8)

We note that, for homogeneous and isotropic flows with stationary Lagrangian velocities $v(t)$, Equation (8) is exact. In order to arrive at the scaling laws (1), Taylor [5] makes further simplifying assumptions. In particular, he assumes that, for small t, the autocorrelation function is constant. As a result, for small t, we obtain ${\langle \left|X(t)\right|}^2\rangle \sim 2\langle V^2\rangle {\int }_0^t(t-\tau )\tau =\langle V^2\rangle t^2$ (note that the variance $\langle V^2\rangle$ is time-independent since we assumed the Lagrangian velocity is a stationary stochastic process). Furthermore, assume that ${\int }_0^{\infty }R(\tau )\tau$ and ${\int }_0^{\infty }\tau R(\tau )\tau$ exist and are finite. Then, for large t, we have $\langle |X(t){|}^2\rangle \sim 2Dt$ where $D=\langle V^2\rangle {\int }_0^{\infty }R(\tau )\tau$.

To summarize, for the scaling laws (1) to hold, the stochastic Lagrangian velocities $v(t)$ must be stationary. Furthermore, the autocorrelation function $R(\tau )$ must decay to zero fast enough so that the integrals ${\int }_0^{\infty }R(\tau )\tau$ and ${\int }_0^{\infty }\tau R(\tau )\tau$ are finite. In Section 5, we show that some of these assumptions do not hold generally for particles advected by random surface waves. Our results were obtained by numerically integrating an exact model of the fluid trajectories on the free surface. We introduce this model in the next section.

3. John–Sclavounos Equation

The John–Sclavounos (JS) [2,3] equation describes the horizontal motion of fluid particles on a free surface $\zeta$. Lets denote the horizontal coordinate by x and denote the vertical coordinate (corresponding to the depth) by z. We assume that the free surface is a graph over the horizontal coordinate x so that on the free surface $z=\zeta (x,t)$ (see Figure 1). The position of the particles on the free surface at time t is then given by $(x(t),\zeta (x(t),t))$. The Eulerian velocity field inside the fluid is denoted by $\mathbf{u}(x,z,t)$. The fluid velocity $\mathbf{u}$ and the free surface $\zeta (x,t)$ satisfy the water wave equations.

With this notation, the JS equation reads

$$\ddot{x}=-\frac{{\zeta }_{,xx}{\dot{x}}^2+2{\zeta }_{,xt}\dot{x}+(g+{\zeta }_{,tt})}{1+{\zeta }_{,x}^2}{\zeta }_{,x},$$

(9)

where ${\zeta }_{,x}$ is shorthand notation for the partial derivative $\partial \zeta /\partial x$ and similarly for partial derivatives with respect to time. Supplying Equation (9) with the appropriate initial conditions $(x(0),\dot{x}(0))$ and integrating in time, the horizontal motion of the particles on the free surface can be computed. The JS equation is quite a remarkable result. It implies that if the surface elevation $\zeta$ is known, then the motion of the particles on the free surface can be deduced without knowing the full fluid velocity field $\mathbf{u}(x,z,t)$.

Denoting the horizontal Lagrangian velocity of a fluid parcel by $v=\dot{x}$, we write the JS equation as a set of first-order differential equations

$$\dot{x}=v,\dot{v}=-\frac{Q(x,v,t)}{1+{\zeta }_{,x}^2}{\zeta }_{,x},$$

(10)

where

$$Q(x,v,t):={\zeta }_{,xx}{v}^2+2{\zeta }_{,xt}v+(g+{\zeta }_{,tt})=\frac{{}^2}{t^2}\zeta (x(t),t)+g.$$

(11)

Fedele et al. [4] showed that the JS Equations (10) have a Hamiltonian structure which in the 1D case is given by

$$\dot{x}=\frac{\partial H}{\partial p},\dot{p}=-\frac{\partial H}{\partial x},$$

(12)

where the Hamiltonian H reads

$$H(x,p,t)=\frac{1}{2}\frac{{(p-{\zeta }_{,t}{\zeta }_{,x})}^2}{1+{\zeta }_{,x}^2}+g\zeta -\frac{1}{2}{\zeta }_{,t}^2,$$

(13)

and the generalized momentum p is given by $p=(1+{\zeta }_{,x}^2)v+{\zeta }_{,t}{\zeta }_{,x}$.

In the time-independent case, where ${\zeta }_{,t}\equiv 0$, the Hamiltonian is a conserved quantity. In other words, the quantity,

$$H_0=\frac{1}{2}\frac{p^2}{1+{\zeta }_{,x}^2}+g\zeta ,$$

(14)

is invariant along the trajectories of the JS equation. However, in the realistic situation where the free surface elevation is time-dependent, the energy H is no longer conserved and complex particle motion is possible.

The derivation of Fedele et al. [4] also shows that the JS equation holds more generally for any particle constrained on a surface $\zeta (x,t)$ and moving under the gravitational force. In particular, they show that the JS equation can be derived without making use of the continuity equation.

In the following, we do not leverage the Hamiltonian structure of the JS equation. Instead, for a given surface elevation $\zeta (x,t)$, we numerically integrate the JS Equation (10) and compute the resulting particle dispersion properties. Although here we focused on unidirectional waves, the JS equation is applicable to two-dimensional waves (see Refs. [3,4]) and therefore our results can be generalized in a straightforward fashion.

We point out that a model similar to the JS equation was proposed by Herbers and Janssen [27] (see their Equation (20)). That model however is inaccurate since it neglects the denominator in Equation (9).

4. Set-Up

4.1. Irregular Wave Field

We consider random surface waves in deep water consisting of a finite sum of plane waves,

$$\zeta (x,t)=\sum _{i=1}^n{a}_i\cos (k_{i}x-{\omega }_{i}t+{\varphi }_i),$$

(15)

where $k_i$ is the wave number, ${\omega }_i$ is the wave frequency satisfying the linear dispersion relation ${\omega }_i=\sqrt{g{k}_i}$. The random phases ${\varphi }_i$ are uniformly distributed over the interval $[0,2\pi ]$.

We consider waves that follow the Joint North Sea Wave Project (JONSWAP) spectrum [28]

$$S(\omega )=\frac{\alpha g^2}{{\omega }^5}exp\left[-\frac{5}{4}{\left(\frac{{\omega }_p}{\omega }\right)}^4\right]{\gamma }^{exp\left[-\frac{{(\omega /{\omega }_p-1)}^2}{2{\beta }^2}\right]},$$

(16)

where $\gamma =3.3$ is the enhancement factor and ${\omega }_p$ is the angular frequency corresponding to the peak of the spectrum. The standard deviation $\beta =0.07$ for $\omega \le {\omega }_p$ and $\beta =0.09$ for $\omega >{\omega }_p$. The amplitude $\alpha$ will be specified later according to the desired wave steepness.

The wave amplitudes in Equation (15) are set to $a_i=\sqrt{2S({\omega }_i)\mathsf{\Delta }\omega }$ so that the random wave field $\zeta$ has the JONSWAP spectrum (16). This yields ${\sigma }^2:=\langle {\zeta }^2\rangle =\frac{1}{2}{\sum }_{i=1}^n{a}_i^2={\sum }_{i=1}^{n}S({\omega }_i)\mathsf{\Delta }\omega$, where $\sigma$ is the standard deviation of the wave elevation $\zeta$. The significant wave height $H_s$ is defined as four times this standard deviation, $H_s=4\sigma$. Following [29], we define the average wave steepness as $ϵ=H_s{k}_p/2$.

4.2. Initial Conditions

The JS Equation (10) is a two-dimensional system of first-order differential equations. In order to numerically integrate these equations, we need to supply them with the appropriate initial conditions $(x(0),v(0))=(x_0,v_0)$. Since the surface waves (15) are stochastically homogeneous, the choice of the initial position $x_0$ does not alter the final ensemble averaged quantities. Therefore, we simply set $x_0=0$.

The initial velocity $v_0$ should be consistent with the Eulerian velocity induced by the wave motion [3]. Since this Eulerian velocity is unknown (unless one solves the full water wave equations), we have to make an assumption for $v_0$. Here, we set the initial particle velocity $v_0$ to coincide with the induced velocity of a linear wave. The velocity potential of a monochromatic linear wave at the surface is $\varphi =(ag/\omega )sin(kx-\omega t)$. This implies the horizontal velocity $v=\partial \varphi /\partial x=(agk/\omega )\cos (kx-\omega t)$. Since the surface is a superposition of linear waves, each component of the random wave field contributes differently to the horizontal particle velocity. For simplicity, we assume that the main contribution comes from the wave corresponding to the peak of the spectrum. Therefore, we set $v_0=ag{k}_p/{\omega }_p$. Furthermore, the wave amplitude a satisfies $\langle {\zeta }^2\rangle =\frac{1}{2}{a}^2$. This implies

$$v_0=\frac{\sqrt{2}\sigma g{k}_p}{{\omega }_p},$$

(17)

where $\sigma ={\langle {\zeta }^2\rangle }^{\frac{1}{2}}$ is the standard deviation of the wave elevation. We have checked numerically that the resulting ensemble averaged Lagrangian velocity is time-independent and $\langle v(t)\rangle =v_0$, thereby confirming that the choice of the initial velocity does not affect the longterm ensemble behavior of the particles. Furthermore, we have perturbed this initial velocity $v_0$ and observed no significant change in the results reported in Section 5.

Although the initial velocity $v_0$ is deterministic, $v(t)$ is stochastic for $t>0$. This is because v satisfies the JS Equation (10) whose right-hand-side inherits the stochasticity of the surface elevation (15).

4.3. Non-Dimensional Variables

In order to nondimensionalize the variables, we rescale space and time as $x\mapsto x/{\lambda }_p$ and $t\mapsto t/T_p$, respectively. The length scale ${\lambda }_p$ is the wavelength corresponding to the peak of the JONSWAP spectrum (16) and $T_p=2\pi /{\omega }_p$ is the associated wave period.

The gravitational constant in non-dimensional form becomes $g\mapsto g{T}_p^2/{\lambda }_p$. The linear dispersion relation ${\omega }_p^2=g{k}_p$ implies $g{T}_p^2/{\lambda }_p=2\pi$, i.e., the non-dimensional gravitational constant is $g=2\pi$. The initial velocity (17) in terms of non-dimensional variables becomes $v_0=ϵ\sqrt{2}/2$.

5. Numerical Results

In this section, we present the results obtained from numerically integrating the JS Equation (10). All results are reported in the non-dimensional variables discussed in Section 4.3. First, we generate the random wave fields from (15) with $n=200$ and the JONSWAP spectrum (16). The computational domain is $x\in [0,200\pi ]$ with periodic boundary conditions. The domain is large enough to avoid finite-box effects. The parameter $\alpha$ in the JONSWAP spectrum (16) is set to $\alpha =210{(ϵ/2)}^2$ which guarantees that the resulting waves have average steepness $ϵ$. We report our results for three wave steepness values $ϵ=0.05$, $0.075$ and $0.1$ which are below the threshold for breaking waves (Recall that JS equation is not valid for breaking waves). Figure 2 shows a realization of the random wave field with steepness $ϵ=0.05$.

Given a wave field $\zeta$, we integrated the JS Equation (10) with the embedded Runge–Kutta scheme RK5(4) of MATLAB [30]. The initial conditions were $x(0)=x_0=0$ and $v(0)=v_0$ with the initial velocity $v_0$ discussed in Section 4.2 (see Equation (17)). For each steepness $ϵ$, we generated 50,000 random waves. For each realization of the waves, we computed the particle trajectories and estimated the ensemble averages $\langle ·\rangle$ from these 50,000 trajectories.

First, we examine the average drift $\langle x(t)\rangle$ as shown in Figure 3a. The average drift grows linearly with time t which agrees with Stokes’ prediction [1,31]. However, the variance $\langle |X(t){|}^2\rangle$ of the particle positions exhibited a surprising behavior. Recall that $X(t)=x(t)-\langle x(t)\rangle$ is the mean-zero position of the particles. Figure 3b shows this variance as a function of time. For short times $\langle |X(t){|}^2\rangle$ grows with the power law scaling $t^4$. After roughly one wave period, the growth of the variance changes and scales as $t^2$. These two regimes result in the scaling law (2).

Note that this behavior is very different from Taylor dispersion theory (1) in the absence of waves. Taylor dispersion predicts a scaling $t^2$ for short times while, in the presence of waves, we observe the scaling $t^4$. For large times, Taylor dispersion predicts a Brownian-type diffusion where the variance increases linearly with time while, in the presence of waves, we observe a ballistic dispersion with ${\langle \left|X(t)\right|}^2\rangle \propto t^2$.

This super-diffusive asymptotic behavior has been reported in homogeneous, isotropic turbulence where ${\langle \left|X(t)\right|}^2\rangle \propto t^{\gamma }$ with $\gamma >1$ [32,33,34]. In turbulence, the departure from Taylor dispersion theory is often associated with intermittency which manifests itself as heavy tails in the distributions of the particle positions and velocities. However, for particle dispersion by waves, we do not observe such heavy-tail statistics. In fact, as Figure 4 shows, the PDF of particle positions and velocities are Gaussian. Note that the initial conditions $(x_0,v_0)$ are deterministic and therefore their initial distributions are delta functions. However, as shown in Figure 4, they rapidly converge to a Gaussian for $t>0$ (e.g., see the circles marking the PDFs at time $t/T_p=0.1$). As we show in Appendix A, these Gaussian distributions can be deduced from the JS equations and a central limit theorem.

The unusual scaling (2) observed here can be explained by examining the autocorrelation function $R(\tau )=\langle V(t)V(t-\tau )\rangle /\langle |V(t){|}^2\rangle$. For t large enough, this function is independent of t and only depends on the delay $\tau$. Figure 5 shows $R(\tau )$ computed for $t=100T_p$. One important feature of this autocorrelation function is that it does not decay to zero for large $\tau$. Instead, it decays to $R(\tau )=1/2$, indicating that the Lagrangian particle velocities on the free surface remain correlated indefinitely.

This correlation explains the ballistic motion of particles observed in Figure 3b for large t, i.e., ${\langle \left|X(t)\right|}^2\rangle \propto t^2$. Recall Equation (8) that relates the variance of the particles to the autocorrelation function. For large $\tau$, the autocorrelation tends to $1/2$ and therefore the integral ${\int }_0^t(t-\tau )R(\tau )\tau$ scales as $\int (t-\tau )/2\tau \sim t^2/4$ as $t\to \infty$. Furthermore, the variance $\langle |V(t){|}^2\rangle$ is constant for large t (see Figure 6). These numerical observations, together with Equation (8), imply that ${\langle \left|X(t)\right|}^2\rangle \sim t^2$.

Finally, we turn our attention to the short-term behavior ${\langle \left|X(t)\right|}^2\rangle \sim t^4$ reported in Figure 3b. Figure 6 shows that, for $t<T_p$, the variance $\langle |V(t){|}^2\rangle$ increases quadratically in time. This together with Equation (8) implies that the variance $\langle |X(t){|}^2\rangle$ of the fluid parcels should in fact increase as $t^4$.

6. Conclusions

Here, we investigated the nonlinear dispersion of fluid parcels on the free surface of a random wave in deep ocean. The free surface is assumed to be a superposition of linear waves with random phases and the JONSWAP spectrum. However, the fluid particle trajectories are computed using an exact nonlinear model known as the John–Sclavounos (JS) equation [2,3,4].

Our main finding is the breakdown of Taylor’s dispersion theory [5]. In particular, for large times, the particles disperse ballistically as opposed to diffusively. More precisely, the variance $\langle |X(t){|}^2\rangle$ of the distribution of the particles increases quadratically in time for large t so that ${\langle \left|X(t)\right|}^2\rangle \sim t^2$. For short times, the variance is proportional to $t^4$ giving rise to the scaling laws (2).

We showed that this unusual scaling law is a consequence of longterm correlation of Lagrangian velocities of the fluid particles. This is a clear violation of Taylor’s assumption that these correlations decay to zero asymptotically. Clearly, our results have profound implications for modeling the dispersion of fluid particles (and pollutants) on the ocean surface.

Since the JS equations are valid in two dimensions, extending our results to two-dimensional waves is straightforward and will be presented in future work. Furthermore, future work will investigate the single-particle dispersion on the surface of (weakly) nonlinear waves. This can be accomplished, for instance, by the one-way coupling of the JS equation to an envelop equation for the free surface (e.g., the nonlinear Schrödinger equation).