Clustering of floating tracers in weakly divergent velocity fields

TL;DR

This paper investigates how floating tracers cluster in weakly divergent ocean surface flows, extending existing theories to include combined flow components and analyzing the effects on clustering rates and patterns.

Contribution

It extends asymptotic clustering theory to weakly divergent flows by analyzing combined flow components and introduces statistical topography for systematic clustering studies.

Findings

Complete clustering persists even in weakly divergent flows but at a slower rate.

Clustering behavior depends on the combination of solenoidal and potential flow components.

Methodology enables systematic analysis of clustering in complex flow regimes.

Abstract

This work focuses on buoyant tracers floating on the ocean surface and treats the geostrophic and ageostrophic surface velocities as the 2D solenoidal (non-divergent) and potential (divergent) flow components, respectively. We consider a random kinematic flow model and study the process of clustering, that is, aggregation of tracers in localized spatial patches. An asymptotic theory exists only for strongly divergent velocity fields and predicts complete clustering, that is, emergence in the large-time limit of spatial singularities containing all available tracers. To extend the theory, we consider combinations of the solenoidal and potential velocity components and explore the corresponding regimes of the clustering process. We have found that, even for weakly divergent flows complete clustering still persists but occurs at a significantly slower rate. For a wide range of parameters, we have analyzed this process, as well as the other type of clustering, referred to as fragmentation clustering, and the coarse-graining effects on clustering, and interpreted the results. For the analyses we have considered ensembles of Lagrangian particles representing the tracer, then, introduced and applied the methodology of statistical topography, which paves the way for systematic studies of clustering in progressively more complicated flows, and for both passive and floating tracers.

Equations

$$\biggl{(}{\frac{\partial}{\partial t}+{\bf{U}}({\bf{R}},t)\frac{\partial}{{\partial{\bf{R}}}}}\biggr{)}C({\bf{R}},t)=\bm{\kappa}{\Delta_{\bf{R}}}C({\bf{R}},t)\,,\,\,\,C({\bf{R}},0)={C_{0}}({\bf{R}})\,;$$

$$\biggl{(}{\frac{\partial}{{\partial t}}+\frac{\partial}{{\partial{\bf{R}}}}{\bf{U}}({\bf{R}},t)}\biggr{)}\rho({\bf{R}},t)=\bm{\kappa}{\Delta_{\bf{R}}}\rho({\bf{R}},t)\,,\,\,\,\rho({\bf{R}},0)={\rho_{0}}({\bf{R}})\,;$$

$$\mathbf{U}(\mathbf{R},t)=\gamma {\mathbf{U}}^p(\mathbf{R},t)+(1-\gamma ){\mathbf{U}}^s(\mathbf{R},t),$$

$$B_{\alpha \beta }^j({\mathbf{R}}^{\prime },\eta )=⟨U_{\alpha }^j(\mathbf{R},t)U_{\beta }^j(\mathbf{R}+{\mathbf{R}}^{\prime },t+\eta )⟩=\int d\mathbf{k}{E}_{\alpha \beta }^j(\mathbf{k},\eta )e^{i{\mathbf{kR}}^{\prime }},$$

$$E_{\alpha \beta }^p(\mathbf{k},\eta )=E^p(k,\eta )\frac{k_{\alpha }{k}_{\beta }}{k^2},E_{\alpha \beta }^s(\mathbf{k},\eta )=E^s(k,\eta )\left({\delta }_{\alpha \beta }-\frac{k_{\alpha }{k}_{\beta }}{k^2}\right),$$

$$B_{\alpha \beta }^j(0,0)=⟨U_{\alpha }^j(\mathbf{R},t)U_{\beta }^j(\mathbf{R},t)⟩=\frac{1}{2}{\sigma }_U^2{\delta }_{\alpha \beta },$$

$$D^p=\underset{0}{\overset{\infty }{\int }}d\eta \int d\mathbf{k}{k}^2{E}^p(k,\eta )=\underset{0}{\overset{\infty }{\int }}d\eta ⟨{\nabla }_{\mathbf{R}}\mathbf{U}(\mathbf{R},t+\eta ){\nabla }_{\mathbf{R}}\mathbf{U}(\mathbf{R},t)⟩,$$

$$D^s=\underset{0}{\overset{\infty }{\int }}d\eta \int d\mathbf{k}{k}^2{E}^s(k,\eta )=\frac{1}{2}\underset{0}{\overset{\infty }{\int }}d\eta ⟨\nabla \times \mathbf{U}(\mathbf{R},t+\eta )\nabla \times \mathbf{U}(\mathbf{R},t)⟩.$$

$$B_{kl}^{eff}\left(\mathbf{r}\right)=\underset{0}{\overset{\infty }{\int }}d\tau B_{kl}\left(\mathbf{r},\tau \right),B_{kl}^{eff}\left(0\right)=D^0{\delta }_{kl},\frac{\partial }{\partial r_i}{B}_{kl}^{eff}\left(0\right)=0,$$

$$-8\frac{{\partial }^2}{\partial r_i\partial r_j}{B}_{kl}^{eff}\left(0\right)=D^s\left(2{\delta }_{kl}{\delta }_{ij}-{\delta }_{ki}{\delta }_{lj}-{\delta }_{kj}{\delta }_{li}\right)+D^p\left(2{\delta }_{kl}{\delta }_{ij}+{\delta }_{ki}{\delta }_{lj}+{\delta }_{kj}{\delta }_{li}\right).$$

$$P(t;\rho )=\frac{1}{2\rho \sqrt{\pi t/\tau }}\exp \left\{-\frac{{\ln }^2(\rho e^{t/\tau }/{\rho }_0)}{4t/\tau }\right\},$$

$$\Phi (t;\rho )\equiv {\int }_0^{\rho }d{\rho }^{\prime }P(t;{\rho }^{\prime })=\mathrm{Pr}\left(\frac{{\ln }^2(\rho e^{t/\tau }/{\rho }_0)}{2\sqrt{t/\tau }}\right),$$

$$\mathrm{Pr}(z)=\frac{1}{\sqrt{2\pi }}{\int }_{-\infty }^{z}dx\exp \left\{-\frac{x^2}{2}\right\}.$$

$$\phi (\mathbf{R},t;\rho )=\delta (\rho (\mathbf{R},t)-\rho ),$$

$$S(t;\bar{\rho })=\int d\mathbf{R}\theta (\rho (\mathbf{R},t)-\bar{\rho })=\int d\mathbf{R}{\int }_{\bar{\rho }}^{\infty }d{\rho }^{\prime }\phi (\mathbf{R},t;{\rho }^{\prime }),$$

$$M(t;\bar{\rho })=\int d\mathbf{R}\rho (\mathbf{R},t)\theta (\rho (\mathbf{R},t)-\bar{\rho })=\int d\mathbf{R}{\int }_{\bar{\rho }}^{\infty }d{\rho }^{\prime }{\rho }^{\prime }\phi (\mathbf{R},t;{\rho }^{\prime }).$$

$$⟨S(t;\bar{\rho })⟩=\int d\mathbf{R}{\int }_{\bar{\rho }}^{\infty }d{\rho }^{\prime }P(\mathbf{R},t;{\rho }^{\prime }),$$

$$⟨M(t;\bar{\rho })⟩=\int d\mathbf{R}{\int }_{\bar{\rho }}^{\infty }d{\rho }^{\prime }{\rho }^{\prime }P(\mathbf{R},t;{\rho }^{\prime }),$$

$$P(\mathbf{R},t;\rho )=⟨\delta (\rho (\mathbf{R},t)-\rho )⟩,$$

$$⟨s_{hom}(t;\bar{\rho })⟩=⟨\theta (\rho (\mathbf{R},t)-\bar{\rho })⟩=P\{\rho (\mathbf{R},t)>\bar{\rho }\}={\int }_{\bar{\rho }}^{\infty }d{\rho }^{\prime }P(t;{\rho }^{\prime }),$$

$$⟨m_{hom}(t;\bar{\rho })⟩=\frac{1}{{\rho }_0}{\int }_{\bar{\rho }}^{\infty }d{\rho }^{\prime }{\rho }^{\prime }P(t;{\rho }^{\prime }),$$

$$\underset{t\to \infty }{\lim }⟨s_{hom}(t;\bar{\rho })⟩\to 0,\underset{t\to \infty }{\lim }⟨m_{hom}(t;\bar{\rho })⟩\to \frac{1}{{\rho }_0}⟨\rho (t)⟩,$$

$$⟨s_{hom}(t,\bar{\rho })⟩=Pr\left(\frac{\ln ({\rho }_0{e}^{-D^{p}t}/\bar{\rho })}{\sqrt{2D^{p}t}}\right),⟨m_{hom}(t,\bar{\rho })⟩=Pr\left(\frac{\ln ({\rho }_0{e}^{D^{p}t}/\bar{\rho })}{\sqrt{2D^{p}t}}\right),$$

$$⟨s_{\mathrm{hom}}(t,\bar{\rho })⟩=\mathrm{P}\{\rho (\mathbf{R},t)>\bar{\rho }\}\approx \sqrt{\frac{{\rho }_0}{\pi \bar{\rho }{D}^{p}t}}{\mathrm{e}}^{-\frac{1}{4}{D}^{p}t},$$

$$⟨m_{\mathrm{hom}}(t,\bar{\rho })⟩/{\rho }_0\approx 1-\sqrt{\frac{\bar{\rho }}{\pi {\rho }_0{D}^{p}t}}{e}^{-\frac{1}{4}{D}^{p}t}.$$

$$\Phi (t;z^{\ast }(t)={\int }_0^{z^{\ast }(t)}d{z}^{\prime }P(t;z^{\prime })=\frac{1}{2},$$

$$l^{\ast }(t)=\exp \left\{\frac{1}{4}(D^s-D^p)t\right\}.$$

$$\frac{1}{4}{D}^{s}t\ll \ln \frac{l_{corr}}{l_0},$$

$$⟨L(t,\rho )⟩=L_0{e}^{D^{s}t},$$

$$U_{\beta }^p\left(\mathbf{R},t\right)={\sigma }_{\mathbf{U}}\int d\mathbf{k}\left(a\left(k,t\right)+ib\left(k,t\right)\right)\frac{k_{\beta }}{k}\exp \left(i\mathbf{kR}\right),$$

$$U_x^s\left(\mathbf{R},t\right)={\sigma }_{\mathbf{U}}\int d\mathbf{k}\left(a\left(k,t\right)+ib\left(k,t\right)\right)\frac{k_y}{k}\exp \left(i\mathbf{kR}\right),$$

$$U_y^s\left(\mathbf{R},t\right)=-{\sigma }_{\mathbf{U}}\int d\mathbf{k}\left(a\left(k,t\right)+ib\left(k,t\right)\right)\frac{k_x}{k}\exp \left(i\mathbf{kR}\right),$$

$$⟨a\left(k,t\right)⟩=⟨b\left(k,t\right)⟩=⟨a\left(\mathbf{k},t\right)b\left({\mathbf{k}}^{\prime },t^{\prime }\right)⟩=0,$$

Full text

floating tracers, material clustering, stochastic velocity fields

††preprint: APS/123-QED

I Introduction

In turbulent oceanic and atmospheric flows, the phenomenon of clustering, that is, spatial aggregation of various tracers or material objects, manifests itself in different physical situations and on different scales Okubo (1980); McComb (1990); Jacobs et al. (2016). More generally, clustering processes are studied in marine ecosystems Okubo (1980), the dynamics of clouds Kostinski and Shaw (2001), physics of porous media Dagan (1987), palaeontology de Pater and Lissauer (2001), and cosmology Shandarin and Zeldovich (1989), to name a few. Clustering of floating debris and plastic litter in garbage islands Law et al. (2010); Cozar et al. (2014), and clustering of oil spills and sargassum in patches Jacobs et al. (2016) have been actively researched problems of fundamental importance and relevant to our study.

Our focus is on clustering of floating tracers at the ocean surface Thomas et al. (2008); S. (2012); Dencausse et al. (2014); Kalda J. (2014); Haza A. C. (2016); Gaube P. (2017); Väli et al. (2018). Here, being ”floating” is crucial, because, as we show, passive tracers at the ocean surface cluster differently. Apparently, the clustering process is controlled by the multiscale and transient surface velocity, and this motivates systematic theoretical study of clustering in random, progressively more complicated and realistic, kinematic velocity fields, and further in dynamically constrained velocity fields.

At the ocean surface, the dominant (i.e., leading-order) geostrophic velocity can be treated as non-divergent and purely solenoidal. The corresponding ageostrophic velocity corrections to the leading-order fields are 3D, weakly 2D-divergent at the surface, and characterized by small vertical velocities Reznik and Grimshaw (2001); Pallàs-Sanz and Viùdez (2007); Zhmur (2011); Wang et al. (2012); McWilliams (2016); Zhang and Qiu (2018). At the surface, this flow correction can be treated as divergent (potential) component of the combined 2D velocity field responsible for advection and clustering of material tracers. Maxey (1987); Boffetta et al. (1994); Klyatskin (1994); Klyatskin and Saichev (1997); Klyatskin (2003); Cressman and Goldburg (2003); Klyatskin (2004); Monchaux et al. (2010); Bourgoin and Xu (2014); Klyatskin (2015, 2016); Gutierrez and Aumaitre (2016).

Asymptotic theories of clustering are developed only for random and purely divergent kinematic flows Klyatskin (2015, 2016); Klyatskin et al. (1996a, b); Saichev and Woyczynski (1996). Extension to weakly divergent, random kinematic flows was attempted in asymptotic sense, but the corresponding predictions have never been verified, and the underlying clustering mechanisms and other properties remain poorly understood Väli et al. (2018). This paper aims to close this gap.

We investigate clustering in a random kinematic velocity model and use Lagrangian-particle Monte-Carlo approach. The paper is organized as the following. Section $2$ formulates the problem and outlines the mathematical foundations in terms of the statistical topography methodology. Section $3$ considers the kinematic flow model and numerical implementation of the clustering analysis. The main results are in Section $4$, followed by conclusions in Section $5$.

II Problem formulation

The tracer dynamics at the upper ocean surface can be formulated for both passive and floating tracers, and these dynamical descriptions are fundamentally different, thus, resulting in profoundly different clustering properties. First, let us introduce passive-tracer concentration $C({\bf r},t)$ and floating-tracer density $\rho({\bf r},t)$, which are the fields of interest, both varying in space and evolving in time (here, ${\bf r}$ is the position vector and $t$ is time), subject to a velocity field ${\bf u}({\bf r},t)$. Second, let us discuss differences between $C$ and $\rho$ and explain the terminology. The passive tracer is just a marker of fluid (and fluid particles), and its dynamics is subject to the standard continuity equation for material tracer; in turn, this equation can be formulated for the ocean surface. Since our study is limited to incompressible fluids (i.e., with zero 3D velocity divergence), the concentration $C$ is also incompressible, in the sense that it does not experience velocity divergence. The floating tracer is not passive, in the sense that it experiences the buoyancy force that keeps it floating on the surface, therefore, it is not a marker of fluid (and fluid particles), its dynamics is subject to the continuity equation for floating tracer, and its density can be changed on fluid particles due to surface velocity divergence. The latter effect can be viewed as compressibility of the floating tracer density $\rho$, and this justifies our choice of terminology: ”density” for floating (compressible) tracers, versus ”concentration” for passive (incompressible) tracers. Note, that both tracers experience 2D divergence of the surface velocity, but this effect is fundamentally different from the compressibility effect.

On the ocean surface, the flow velocity can be written as ${\bf u}(x,y,0,t)=({\bf U}(x,y,t),w(x,y,0,t))$, where ${\bf U}$ and $w$ are the horizontal and vertical velocity components, respectively. The dynamical equations governing concentration $C$ and density $\rho$ on the ocean surface are derived in Appendix A:

[TABLE]

[TABLE]

where $C_{0}({\bf R})$ and $\rho_{0}({\bf R})$ are the initial concentration and density distributions, respectively; ${\bf R}=(x,y)$ is the 2D coordinate on the surface; $\Delta_{\bf R}$ is the 2D Laplace operator; and $\bm{\kappa}$ is the diffusivity coefficient. For the sake of initial simplicity, in what follows we will consider the tracer dynamics in an unbounded domain and for $\bm{\kappa=0}$, while recognizing that both boundary and diffusive effects can be important.

Although clustering of the passive-tracer concentration is important, in this study we focus on clustering of the floating-tracer density. In order to distinguish between the two types of tracers and, therefore, two different types of clustering, we will refer to the passive-concentration clustering as the C-clustering, and to the floating-density clustering as the D-clustering. Since most of the presented results are about D-clustering, we will refer to it as simply clustering, and the C-clustering terminology will be retained.

In what follows the governing equations are nondimensionalized by the length scale $L_{0}$, the time scale $t_{0}$, and the density scale $\rho^{*}$. The initial area occupied by the tracer is four, the initial density of the tracer is unity, and the time scale $t_{0}$ is the discretization time step, since the random velocity is taken to be delta-correlated in time, and there is no correlation time scale. The velocity scale is $U^{*}=L_{0}/t_{0}$, and the turbulent diffusion scale is $\kappa^{*}=L_{0}^{2}/t_{0}$. Let us now describe and discuss the flow velocity model, which is purely kinematic thus not constrained by the fluid dynamics for the sake of initial theoretical simplicity.

II.1 Random velocity fields: statistical properties

First, we decompose the random velocity as:

[TABLE]

where superscripts $p$ and $s$ indicate the potential (i.e., potential, purely divergent) and solenoidal (i.e., non-divergent, purely solenoidal) velocity field components, respectively; and $0\leq\gamma\leq 1$ is the non-dimensional parameter setting their relative contributions. Second, we define statistical properties of the random velocity following Klyatskin (1994, 2015); Klyatskin and Koshel (2017): ${\bf{U}}^{p}({\bf{R}},t)$ and ${\bf{U}}^{s}({\bf{R}},t)$ are time-correlated, random Gaussian fields with spatially homogeneous, isotropic and stationary statistics. For such random fields the spatio-temporal velocity correlation tensors are:

[TABLE]

where $\bm{k}=\left(k_{x},\,k_{y}\right)$ is the wavevector, $k=|\bm{k}|$ ; $\alpha$ and $\beta$ stand for $x$ and $y$, respectively; $j$ stands for $p,\,s$; ${\bf R}^{\prime}$ is a 2D spatial shift and $\eta$ is a time lag; the angular brackets denote ensemble averaging over many flow realizations. The spectral densities are chosen to be

[TABLE]

where $\delta_{\alpha\beta}$ is the Kronecker delta. The local (in space and time) velocity correlation tensor is

[TABLE]

where $\sigma_{U}^{2}=B_{\alpha\alpha}({\bf 0},0)=\int d{\bf k}E(k,0)$. Let us also introduce the effective diffusivities Klyatskin (1994, 2015, 2016):

[TABLE]

[TABLE]

These diffusivities arise from the time-lag-integrated correlation tensor:

[TABLE]

and its second derivatives:

[TABLE]

In combination with the statistical topography methodology, these quantities will be used in the clustering analyses of Section $4$.

II.2 Statistical topography of random fields

A brief description of the statistical topography methodology, applied here for a quantitative description of the clustering process, follows Klyatskin (2015). Let density field $\rho({\bf R},t)$ be a random in time process with a probability density function (PDF) homogeneous in space but evolving in time. Let us assume that, from a smooth initial density distribution, the dynamics generates a population of clusters, which are localized regions with high tracer density. A solution corresponding to the log-normal process, which describes the clustering, has the PDF

[TABLE]

where $\tau=1/D^{p}$ is the effective diffusion time scale. The corresponding integral distribution function is

[TABLE]

where $\Pr(z)$ is the probability integral

[TABLE]

Let us also consider the indicator function

[TABLE]

which sifts $\rho({\bf R},t)$ at given $\rho$ via the Dirac delta-function. The total area of the regions, where density exceeds some threshold $\bar{\rho}$, is referred to as the clusters area and obtained as

[TABLE]

where $\theta(\cdot)$ is the Heaviside (step) function. The total mass of the tracer within the clusters area is referred to as the clusters mass and obtained as

[TABLE]

Ensemble averaging of (15) and (16) yields:

[TABLE]

because the PDF (11) is completely local (i.e., involves no lags in time and shifts in space) and will be referred to as one-point PDF. The ensemble mean of the indicator function (14) ensures that:

[TABLE]

and the expressions (17) and (18) become:

[TABLE]

where $s_{hom}(t;\bar{\rho})$ and $m_{hom}(t;\bar{\rho})$ are, respectively, the specific clusters area and the specific clusters mass defined by the threshold $\bar{\rho}$ Klyatskin (2016); Klyatskin and Koshel (2017). Since $\rho({\bf R},t)$ is a positive-definite field, the clustering happens with probability one, and the corresponding limits,

[TABLE]

assert that the clusters area shrinks to zero, while the clusters mass incorporates all the available tracer. From (11) and (20), one can derive the specific functions for the purely divergent velocity field (i.e., for $\gamma=0$):

[TABLE]

where $\rho_{0}$ is the initial uniform density. For times much larger than the diffusion time scale $1/D^{p}$, approximate estimates can be obtained from (23) and (13) by using the large-argument asymptotics of the probability integral Olver et al. (2010):

[TABLE]

These relations describe the complete clustering that happens exponentially with probability one, and the one-point PDF and such characteristics as clusters area and mass conveniently describe the process. These characteristics, as well as their two-point extensions Klyatskin (2016), are used below for analysis of clustering in previously unexplored, weakly divergent velocity fields.

We also employ the concept of typical realization curve. Given a random process $z(t)$, this curve is defined as the median of the integral distribution function (12) and found as the solution $z^{*}(t)$ of the equation,

[TABLE]

which implies that for any time $t$:    $P[z(t)>z^{*}(t)]=P[z(t)<z^{*}(t)]=1/2$. Typical realization curve of the distance between two particles Klyatskin (2015) is

[TABLE]

When $D^{s}<D^{p}$, the average distance between particles, as given by $l^{*}(t)$, attenuates exponentially in accord with the process of complete clustering. When $D^{s}>D^{p}$, there is no complete clustering, at least for the times such that the following constraint holds Klyatskin (2003):

[TABLE]

where $l_{0}$ is the average initial distance between the particles, and $l_{corr}$ is the spatial correlation radius of the velocity field.

The statistical topography methodology offers a variety of useful characteristics, such as the average length of the isoline contour corresponding to $\rho({\bf R},t)=\bar{\rho}$, which is estimated as

[TABLE]

where $L_{0}$ is its initial length Klyatskin et al. (1996a); Saichev and Woyczynski (1996). Thus, the contour length in the purely non-divergent velocity field grows exponentially, because the tracer patch is continuously stretched and filamented. We refer to this process as tracer fragmentation. In a weakly divergent velocity field, the tracer fragmentation affects and modulates the process of complete clustering, and we address this issue.

III Kinematic velocity model

In this section we introduce the kinematic model of a random velocity field with a prescribed correlation tensor, and then study tracer clustering in the model solutions.

III.1 Model formulation

Following Roberts and Teubner (1995); Zirbel and Cinlar (1997); Koshel and Alexandrova (1999); Klyatskin and Koshel (2017), we use a spectral representation of the velocity field:

[TABLE]

where $\sigma_{U}$ is the standard deviation of the velocity, index $\beta$ stands for $x$ and $y$; $a({\bf k},t)$ and $b({\bf k},t)$ are Gaussian, $\delta$-correlated in time, spectral coefficients, such that:

[TABLE]

The velocity field (III.1) corresponds to the correlation tensor (4); and in the physical space it is obtained by an inverse Fourier transform with a random phase. The spectral density is prescribed in the form

[TABLE]

where $l=l_{corr}$ is the spatial correlation radius. The effective diffusion coefficients (7) – (8) ensue from (31) and (3) by the integration:

[TABLE]

where $t_{0}$ is the time scale.

III.2 Numerical implementation

Given uniformly gridded velocity fields, equations (1) and (2) are solved numerically by a simulation of the ensemble of Lagrangian particles Klyatskin (1994); Gat et al. (1998); Frisch et al. (1998); Koshel and Alexandrova (1999); Koshel et al. (2013, 2015). We do not consider the effect of molecular diffusion here ($\bm{\kappa}=0$), it means that the equations (1) and (2) are first order linear PDE. The simplest means to solve such equations is the method of characteristics Klyatskin (1994); Koshel et al. (2015). Evolution of each Lagrangian particle is governed by

[TABLE]

Here ${\bm{\xi}}$ is the initial position of the trajectory. These equations are integrated in time using the standard Euler-Ito scheme, and the Eulerian density and concentration fields are obtained by the relations:

[TABLE]

We use 900 tracers in each mesh sells from the initial spot of the impurity to obtain good enough Eulerian density and concentration fields.

It is worse to mention that the effect of the molecular diffusion can be considered by adding the additional vector delta-correlated process in right part of the equations (III.2), see references above.

IV Modelling results

In this section we focus on modelling the nondiffusive clustering in the implemented double-periodic domain ($|x|\leq 10,\hskip 10.0pt|y|\leq 10$) with the uniform numerical grid $2048^{2}$. The initial uniform distribution of Lagrangian particles is in the sub-domain $|x|\leq 1,\hskip 10.0pt|y|\leq 1$, and there are $900$ particles placed initially in each grid cell. Results for the strongly and weakly divergent velocity regimes are discussed and compared.

Let us explain how the particle distributions are obtained further below. The random velocity field within each grid cell is treated as piecewise-constant, and this ensures that the velocity field has discontinuities responsible for the white-noise stochastic process Koshel and Alexandrova (1999); Klyatskin and Koshel (2017). Density and concentration of particles are found by averaging over the grid scale. The averaged concentration is the number of particles in the cell divided by the cell area. The averaged density is the number of particles in the cell divided by the area which these particles occupy. Because each particle conserves its mass, the corresponding area experiences velocity divergence and scales as $\rho_{0}/\rho(t)$. Note, that the averaging can sometimes produce spurious concentration values, therefore, some caution with interpretation is needed (see examples towards the end of the paper).

IV.1 Purely divergent and non-divergent flow regimes

First, let us consider the limiting case of purely divergent velocity ($\gamma=1$). Evolution of particles is computed over the diffusion time scale $1/D^{p}$ and for different values of the spatial correlation radius $l_{corr}$ (Fig. 1,2). Clustering is evident in the concentration field (Fig. 2), but it is more pronounced in the density field (Fig. 1). Fig. 3 shows that the clustering area and mass are very close to the asymptotic exponential prediction, even in individual realizations. We also found that the clustering depends weakly on the spatial correlation radius, provided the time is scaled diffusively. The dependency from spacial correlation radius has two aspects. One is the the speed of clustering by virtue of diffusion time scale 32) which depend from value of the correlation radius. Second one is the space structure of the clusters. We can see qualitatively that particles aggregates in thin lines and the distances between these lines are proportional to spacious correlation radius.

Let us now compare clustering in the purely divergent ($\gamma=1$) and non-divergent ($\gamma=0$) flow regimes (Figs. 4 and 5). In the latter regime there is no difference between the density and concentration clustering, because the density is also materially conserved. Furthermore, in this regime the asymptotic theory (26) predicts no complete clustering at all. Instead, both fields are affected by the stirring; as a result, they develop clusters due to the fragmentation process, which is fundamentally different from the complete clustering process. These clusters are transient and never aggregate large values of density and concentration. This is because the stirring process reshuffles the particles at each time instant because of the kinematic nature of the random velocity field. Presence of dynamical constraints in the velocity field might in fact significantly alter this type of clustering leading to anomalously large values provided the velocity field is induced by coherent structures such as vortices and jets. Investigating this possibility is beyond the scope of the paper and will be addressed elsewhere. As a result, the density excess over threshold value is moderate: the number of grid cells with $\rho({\bf R},\tau)>2$ is small. Moreover, the maximum density value (about 5 after $t=3\tau$) is significantly smaller than in the purely divergent case (where it is more than 100 at $t=\tau$). Strictly speaking, these high density (or concentration) values do not signify the complete clustering in the sense of estimates (26) and (28), but reflect the fragmentation clustering. Fig. 3,5 shows the clustering area and mass based on the density larger than $\bar{\rho}=1$. The clustering mass decreases, but less so for a large correlation radius ($l_{corr}=0.16$), because of the slowed down fragmentation process. The clusters are typically small and uniformly distributed, and there are no clusters on the scale of $l_{corr}$. We hypothesize that fragmentation clustering can dominate over exponential clustering in the presence of diffusion ($\kappa\neq 0$), but this analysis is beyond the scope of the paper.

IV.2 Weakly divergent regime

Let us now consider a more realistic flow regime with equal contributions of the solenoidal (non-divergent) and potential (divergent) components ($\gamma=0.5$; Fig. 6). The asymptotic estimate (26) still predicts no clustering; however, we have found that the clustering not only occurs but is also very efficient. Overall, the clusters develop at a much slower rate (Fig. 7), relative to the purely divergent flow regime, but the complete clustering is eventually reached. The clustering slowdown takes place for all $l_{corr}$ considered.

Let us now consider weakly divergent flow regimes with $\gamma=0.2$ (Fig. 8) and also $\gamma=0.1$ to illustrate the trends. Because of the strong solenoidal component, the particles are faster scattered away from their initial positions, hence we had to increase the model domain. Although the asymptotic theory predicts no complete clustering, we have found that clustering still takes place, but at a significantly slower rate. Figure 9a shows that initial clustering is faster due to the fragmentation process, but then it slows down. The complete clustering limit is reached at about $t=25\tau$. Stirring induced by the solenoidal flow component is significant, and in accord with this: $D^{s}\approx 16D^{p}$. Concurrently, the clustering process evolves over the longer solenoidal diffusion time scale $D^{s}$ (Fig. 9). Clusters on the scale $l_{corr}$ do not emerge due to the shearing effect of the solenoidal component, which facilitates fragmentation of the clusters. Complete clustering still occurs, even for $\gamma=0.1$ (Figs. 10 and 11), but with progressively slower rates. Dependence on $l_{corr}$ is similar to the $\gamma=0.2$ regime, as manifested by the clustering curves.

Let us now discuss effects of spatial averaging (i.e., coarse-graining) on the concentration. If the averaging is too coarse, it may appear that all particles within the averaging region are completely clustered while the corresponding clustering mass curves have not completely saturated yet (see Fig. 9, where the clusters are already red-colored). To explain this, let us consider in detail the regime $\gamma=0.2$ and refine its spatial averaging by $100$ times in each direction (Fig. 12). The finer structure of clusters, otherwise smeared out by the averaging, becomes apparent. Hence, within each coarse-grained cluster there are in fact many finer clusters that can be revealed by refining the averaging scale. The refined distributions show that many particles have not yet clustered, and this is true even for relatively long times. To summarize, some caution is required to interpret coarse-grained concentration and density fields, and in this regard statistical topography diagnostics are useful.

V Conclusions

This study is motivated by clustering of floating tracers and material objects on the ocean surface. We consider both passive and floating (i.e., buoyant) tracers and show that their dynamics at the surface are fundamentally different. This is because passive-tracer concentration is materially conserved (i.e., remains constant on fluid particles), whereas floating-tracer density (we use this term to distinguish it from the passive-tracer concentration) changes due to the compressibility effect induced by the surface flow convergence. Our main focus is on the density but some comparisons with the concentration are made for the sake of better understanding.

We represent tracers by ensembles of Lagrangian particles and consider their evolution in random kinematic velocity fields. The goal is to establish and interpret the clustering properties induced by 2D velocity fields consisting of divergent (potential) and non-divergent (solenoidal) flow components. In the situation when velocity field is purely divergent (i.e., without a solenoidal component) and delta-correlated in time, the existing asymptotic theory (Klyatskin, 2015), which forms the basis of our work, predicts complete clustering, that is, formation of spatial singularities that aggregate all available tracer. The complete clustering develops exponentially in time and depends on the velocity correlation length scale. If the divergent flow component is weaker than the solenoidal one, the asymptotic theory predicts no complete clustering — this regime has never been explored. In this paper we close this gap by direct numerical modelling of the tracer fields and subsequent analysis of the solutions by the means of statistical topography methodology.

The presented results demonstrate the opposite to the asymptotic prediction: the complete clustering is still in place but its development becomes significantly slower in time. The time scale in this case is controlled by the solenoidal flow component and its associated diffusion time scale. We also have found that in weakly divergent flow regimes, there is another type of transient tracer aggregation, referred to as fragmentation clustering, which coexists with the complete clustering and interferes with it. This type of clustering is found to be bounded by finite values of concentration (typically $20$-$30$ times larger than the initial one), and it is controlled by the solenoidal flow component that disperses particles and, thus, inhibits the complete clustering. In this regime clustering depends on the velocity correlation length scale.

The analysis and results of this paper should be of interest to ocean and atmosphere modelers, because realistic ocean surface flows are weakly divergent due to ageostrophic motions. Future research extensions considering effects of large-scale flows, explicit eddy diffusion, finite life time of the tracer and inertial effects due to the buoyancy and finite size of floating objects. Including finite time correlations in the kinematic flow model, as well as considering dynamically constrained and progressively more realistic flows are also obvious avenues for new advances on the clustering problem.

Acknowledgements.

This study was supported by the POI FEB RAS Program ”Mathematical simulation and analysis of dynamical processes in the ocean” ($117030110034-7$) and by the Russian Foundation for Basic Research project $17-05-00035$. EAR and PB were supported by the NERC grant $NE/R011567/1$ and by the Royal Society Exchange Grant $IEC/R2/181033$. KVK work on the statistical topography analysis was supported by the Russian Scientific Foundation project $9-17-0006$. DVS was supported by the Russian Foundation for Basic Research project $19-55-10001-KO-a$. VIK was supported by the Program of Fundamental Studies of the Russian Academy of Sciences ”Advanced methods of mathematical modelling for studying nonlinear dynamical systems”.

Appendix A

Given the acting 3D velocity field ${\bf u}\!=\!(u,v,w)$, the continuity equation for the passive-tracer concentration $C_{3D}(x,y,z,t)$ is

[TABLE]

Let us restrict our dynamical description to the 2D surface of the ocean and assume the rigid-lid approximation: $w(x,y,0,t)\!=\!0,$ with the surface velocity expressed as ${\bf U}\!=\!(U,V).$ In this case, the 2D governing equation for the surface concentration $C(x,y,t)\!=\!C_{3D}(x,y,0,t)$ can be written as

[TABLE]

Now, let us invoke the incompressible fluid density continuity equation $\nabla\!\cdot\!{\bf u}\!=\!0,$ and express

[TABLE]

By taking this into account, the governing equation for evolution of the passive-tracer concentration on the ocean surface becomes

[TABLE]

Now, let us derive a similar equation for the floating-tracer density $\rho$ by writing the standard continuity equation

[TABLE]

and noting that the vertical density flux $w\rho\!=\!0,$ because we assume that the additional buoyancy force completely balances and cancels out any vertical flux of the tracer, and the tracer is always restricted to stay on the ocean surface. By taking this into account, the governing equation for evolution of the floating-tracer density $\rho(x,y,t)$ on the ocean surface becomes

[TABLE]

Note, that the last term, which represents density compression by converging velocity, does not have its counterpart in the governing equation for $C$.