Steady state dynamic dependence between local mobility and non-affine fluctuations in two-dimensional aggregates

Motivated by qualitative experimental observations in collective behavior of self-propelled camphor particles at air-water interfaces, we study a generic aggregate forming system in two dimensions using canonical ensemble constant temperature molecular dynamics simulation. The aggregates form due to the competition between short-range attraction and long-range repulsion of pair-wise interactions as a generic proxy for the specific case of short-range capillary attraction competing with long-range Marangoni-assisted repulsion in camphor boat systems. Choosing the appropriate set of interaction parameters, we focus on characterising the local dynamics in two specific limiting morphologies, viz. compact and string-like aggregates. We focus on the temporal evolution of the mobility of an individual particle and the dynamic change in its nearest neighbourhood, measured in terms of the Debye–Waller factor () and the non-affine parameter (), respectively (both defined in the text), and their interrelation over several lengths of observation time . The distribution for both measures are found to follow the relation: for the measured quantity x. The exponent is equal to two and one respectively, for the compact and string-like morphologies following the respective ideal fractal dimension of these aggregates. A functional dependence between these two observables is determined from a detailed statistical analysis of their joint and conditional distributions. The results obtained can readily be used and verified by experiments on aggregate forming systems more generic than the specific camphor boat system that motivated us, such as globular proteins, nanoparticle self-assembly etc. Further, the insights gained from this study might be useful to understand the evolution of collective dynamics in diverse glass-forming systems.

to migration patterns of animals [5] and biological oscillations [6]. The examples are practically endless [7] and the scientific explorations to find, tweak and control such possibilities are ever increasing. Yet we are far from conceptualising any microscopic principle able to predict the steady state outcome of such processes in a generic setting. The diversity of control parameters specific to individual situations and their inherent dynamical nature are the primary stumbling blocks to even find a reasonable starting point for such understanding.
Dynamics of self-assembled structures cannot be completely understood in terms of the dynamics of individual constituents alone. Consider, for example, sprinkling a tiny amount of camphor dust on a bowl of water. As camphor readily undergoes sublimation at room temperature, fluctuating interfacial tension imparts a Marangoni stress [8] which repels the particle from its original position relative to its neighboring particle. This self-propulsion [9] is counter ed by the capillary attraction among the particles. Under the influence of competing forces, the camphor particles show aggregation which is highly dynamic. In a more controlled environment, it is possible to control the local structures of such aggregates. In fact, a continuous structural transition from local hexagonal packing (figure 1(a)) of camphor particles to interconnected string-like aggregates (figure 1(b)) have been observed by changing the interfacial tension and thus changing the repulsive Marangoni stress. As we observe the evolution of such structures (see videos SI1 and SI2 (stacks. iop.org/JPhysCM/32/214004/mmedia)), we notice that the motion of individual particles is expected to be strongly influenced by the motion of their neighbouring particles and vice versa. It is natural to expect that the immediate neighbourhood of a moving particle, which itself is dynamically fluctuating, will get somewhat distorted. This can be thought as an effect of dynamic local strain field caused by the motion of particles. Intuitively, this appears to be one basic dynamic process of all aggregations irrespective of their interactions and preparation protocols. However, no functional quantification of interdepend ence between the individual particle motion and associated local distortion of its neighbourhood is readily available to our knowledge.
Observed morphological variation is not specific to the camphor boat system which is athermal, non-equilibrium and non-Hamiltonian per se. This structural phenomenology is shared by many other thermal, Hamiltonian systems and even in equilibrium as has already been mentioned earlier. One essential ingredient for aggregate formation is the competition between interactions of some sort. To probe the interdependence between particle mobility and its neighbourhood fluctuation, we now choose to analyse the simulated trajectories of a generic particulate system interacting via a competition between short-range attraction and long-range repulsion, in two dimensions (2D). Particles in this system, upon cooling down at a suitable rate to appropriate temperature, naturally form stable clusters of different shapes and sizes, thereby rendering it a reasonably appropriate model to study the correlated particle dynamics in a fluctuating heterogeneous environment. Although it is possible to generate a whole hierarchy of self-assembled structures by tuning the competing interactions, we specifically consider two limiting morphologies, namely, finite-size compact droplet-like structures and extended string-like structures. For characterisation of the local dynamical fluctuations, we focus on studying the time evolution of two specific quantities. The individual particle motion is quantified by the mobility ū 2 i defined as the running average of squared displacement of a particle over certain observation time τ w . The degree of fluctuation of the neighbourhood is measured in terms of non-affine parameter χ which computes the mismatch between the particle neighbourhoods in the beginning and the end of the same observation time. First, we show that these two observables follow the same scaling form with their observation time confirming the collective nature of local dynamics in this aggregate forming system. Next, we present a phenomenological characterisation of collective dynamics by providing a functional relationship between the response χ and its predictor ū 2 i suggested by a master curve obtained from detailed statistical analysis of the distributions of respective quantities.

Model system and simulation details
Formation of aggregates is generally attributed to the competing length scales intrinsic to the particle level interactions [10,11]. We consider a 2D system of identical particles of unit mass interacting pair-wise via an effective potential φ, an algebraic sum of a short-range attraction, φ sA and a long-range repulsion, φ lR . The respective functional forms are as follows: φ sA = 4 {(σ/r) 2n − (σ/r) n } and φ lR = (Aξ/r) exp(−r/ξ). The energy, length and time scales, relevant to the simulation are then expressed in units of , σ and σ 2 / , respectively. Following this definition, the repulsion strength A and screening length ξ are expressed in units of and σ, respectively. This specific choice of potentials is a reliable representation of naturally occurring systems, e.g. globular proteins, [12,13] as well as numerous artificially created, yet highly tuneable systems [14] such as polymer-grafted nanoparticles, [15] colloids in weakly polar solvent [16] etc. We note that pattern formation due to competing interactions has also been observed in more exotic systems such as the pasta phase in neutron stars [2] and vortex cluster in type 1.5 superconductors [17].
The range of attraction is chosen to be 0.2σ by fixing n = 18. For comparison, we state that the attraction ranges up to 2.5σ for the well-known Lennard-Jones system (n = 6). By fixing the repulsion strength equal to that of attraction, A = 4 , the system's morphology at fixed density and temper ature can then be systematically controlled by tuning ξ alone. A comprehensive morphology phase diagram along with extensive structural characterization of the observed hierarchy of aggregates obtained as a function of ξ, A and n has already been shown in a previous study [18]. In that same study, it was shown that the model system undergoes a percolation transition when ξ is lowered below 0.5 and the size of the aggregates becomes considerably smaller if the value of ξ is increased above 0.8 for a fixed density ρ = 0.4 and fixed final temperature T = 0.05k B . We recall that percolation is observed for all values of ξ when the density of the system matches the random percolation critical density, ρ = 0.5 for the specified values of T and A. For the present study, we focus only on the finite-size aggregates.
We consider two specific values, ξ = 0.5 and 0.8, for which the particles arrange themselves into compact finite-size crystalline islands with triangular symmetry and as string-like aggregates, as shown in figures 2(a) and (b), respectively. As described in [18], the typical fractal dimension d f of the compact clusters is d f ∼ 1.6 and for the string-like clusters d f ∼ 1.0, computed employing the scaling relation between the radius of gyration and the typical size of the clusters. The systems are prepared by slowly cooling an equilibrium liquid of density, ρ = 0.4 from initial temperature T i = 1.0k B to final temperature T f = 0.05k B with a linear cooling rate of 10 −4 /τ . For comparison, we mention that our choice of T f and ρ is lower than the critical values T c = 0.18 ± 0.01 (in units) and ρ c = 0.6 ± 0.1 for purely attractive systems [19]. Throughout the canonical constant number-area-temperature (NAT) simulation, temperature T is maintained using a Langevin thermostat as implemented in LAMMPS [20] as already detailed in our previous study [21]. The slow cooling protocol, contrasting the usual rapid quench, is adapted to focus on the role of geometric frustration originating from competing interactions on the anomalous long-time dynamics of the model system. Specifically, the model system never reaches equilibrium even after the total observation time 100τ as the potential energy continues to decay monotonically. Also, we recall that the mean square displacement for both the systems falls out of the ballistic regime within one unit of time, τ , in the scale used for the simulation and continues to evolve sub-diffusively at long time after an intermediate relaxation time ranging roughly up to 20τ . (See [21] for further details.)

Evolution of particle mobility
Given the trajectories, we start our analysis by computing the total displacement fluctuation of a particle with respect to its initial position as is the position of the ith particle at instant t within the measurement period τ w . Figures 2(a) and (b) show the spatial map of this quantity for compact and string-like aggregates, respectively. Considering this quantity as a reliable measure of mobility of an individual particle up to τ w , the distribution of the same P(ū 2 i ) is studied for a set of τ w ranging from τ to 100τ . The effect of the shape of the aggregates is evident for both the studied morphologies, P(ū 2 i ) differs from the normal distribution expected for equilibrium liquids. We find that the holds for both compact and string-like morphologies. For compact clusters, the exponent α is unity and α = 0.5 for string-like aggregates.
mobility of particles forming compact crystalline islands is largely distributed in a scale-free way as plotted in figure 3(a) and the algebraic portion of P(ū 2 i ) gets extended over larger range of ū 2 i with increasing τ w keeping the exponent unchanged. For string-like aggregates, P(ū 2 i ) is distinctly different and at long time, it closely follows log-normal statistics over almost two decades as shown by the solid line in figure 3(b). We have confirmed that the system does not show any overall drift as might be expected from the prominent non-zero peak of the distribution for long τ w . Due to finite observation time, both distributions get truncated at certain values of ū 2 i (τ w ) and decay exponentially as shown by dashed lines in figures 3(a) and (b).
The systematic change of P(ū 2 i ) with respect to τ w can be cast in a simple scaling form: τ α w P(ū 2 i ) ∼ū 2 i /τ α w which holds for both compact and string-like morphologies. With the exponent α = 1.0 and 0.5 for the respective cases, scaled P(ū 2 i ) for all τ w collapses to a single master curve as shown in figure 3(c).
i . γ is a morphology dependent exponent equal to the ideal fractal dimension d f for compact (d f = 2) and string-like (d f = 1) aggregates. This finite-time scaling law ensures that a short time series measurement is sufficient to learn about individual particle motion of our model system at long time. This is important as it might considerably save the huge computational resources and/or experimental efforts often required to study a steady state dynamical property of non-equilibrium systems. Further, this result supports a long term theoretical speculation [22][23][24] of the fractal nature of particle trajectories constituting a fractal object, i.e. aggregates and can be reasoned when mobility is conceived as the volume swept by the particular trajectory up to the measurement time. However, our study is inadequate to comment further on this issue. Reserving this question for further exploration, we move on to investigate the change in the neighbourhood of ith particle due to the collective motion of its neighbouring particles.

Evolution of neighbourhood fluctuation
Unique identification of nearest neighbors in a heterogeneous system is difficult. The traditional way of choosing a cutoff radius equal to the first minimum of the radial distribution function of particles is clearly inadequate to capture the spatial heterogeneity of the system. Purely geometric methods, such as Voronoi analysis [25,26] often overestimates the number of nearest neighbours, especially for an inhomogeneous spatial arrangement of particles such as our model system. Also, it is now well known that Voronoi analysis is inadequate to establish any meaningful relationship between the local structure and microscopic dynamics of particulate systems [27,28]. For the present study, we employ an adaptive neighbor search algorithm based on relative angular distance (RAD) [29]. This RAD algorithm starts with a list of prospective neighbors, within an arbitrarily large cutoff, arranged in an ascending order of radial distance relative to the central particle. Only those particles are then identified as nearest neighbors which are not angularly blocked by any other particle in the list based on the solid angle criterion. We note that this method is both geometric and non-parametric and thus suitable for spatially heterogeneous situations such as our model system. For a simple model liquid system, this method has already been successfully used to reveal interesting interrelation between local structure and dynamics [30,31]. Considering the long-range nature of the used repulsive interaction, first list of prospective neighbors is generated employing a cut-off r c = 10σ where the effect of repulsion asymptotically becomes negligible. This ensures that even for a sparse arrangement of particles, all the energetically relevant particles are considered as nearest neighbors. This is important for the particles at the edge of one cluster and sufficiently far away from the other clusters.
Following the identification of the first coordination shell around each particle, we compute the non-affine parameter [32] defined as χ = j [∆ j (τ w ) − D∆ j (0)] where ∆ j (= r j − r i ) is the relative distance between the ith particle and its j th neighbour. The deformation matrix D attempts to match ∆ j (t) between two particles at any instant t to the ∆ j (0) between the same pair of particles at one arbitrary reference time, t = 0, through minimal affine transformation. Computationally, the nearest neighbors of each particle is determined first, in an instantaneous configuration. Then the position of each particle and its nearest neighbors are traced back in the reference configuration. Summation over all the relevant neighbors, then, essentially returns the error in matching the nearest neighborhood of a central particle at any instant to the space enclosed by the same particles around the same central particle at earlier time. For a fixed measurement time interval τ w , χ is then computed for each of the particles for many different choices of reference configurations.
We note that the restriction of nearest neighbours was not strictly posed in the original definition. However, by doing so, we eliminate the arbitrariness in the choice of neighbours which might have affected the computation. Further, this restriction returns χ with a well-defined geometric quantity, namely, procrustean distance which measures the degree of coincidence between an arbitrary simplex (any nearest neighbour shell, in our case) and a reference simplex (the same shell at initial time). The probability distribution of non-affine parameter P(χ) computed at the end of two different measurement windows, τ w = 5.0, 50.0 is presented in figures 4(a) and (b) for the compact and string-like aggregates, respectively. While stretched exponential feature, shown by solid line, is observed for both cases, compact aggregates tend to show exponentially distributed large χ fluctuations (shown by dashed line) compared to string-like aggregates. It has been attempted recently to express the non-affine fluctuations as the effect of anharmonic fluctuation arising in both finite temperature systems [33,34] and granular media [35] as the systems explore the potential/free energy landscape. Following this view, we expect the stretched exponential behaviour as an outcome of underlying rugged energy landscape available to our model system. Further study in this direction is underway and will be presented elsewhere.
Interestingly, P(χ) computed at different τ w shows temporal dependence very similar to that of P(ū 2 i ). Using the same scaling relation P(χ; τ w ) ∼ τ −γ w χ, excellent data collapse (figure 4(c)) can be achieved for both compact and stringlike aggregates with γ = 2.0 and 1.0 for the respective cases. Mentioning both ū 2 i and χ have the same dimensionality of [L 2 ], we recall that ū 2 i is the total displacement fluctuation of an individual particle over certain observation time; whereas χ is a geometric distance quantifying the overlap between the nearest neighbourhood of the same particle computed at the initial and final instances of the observation period. While the strong interdependence between the mobility of a particle and the evolution of its neighbourhood is evident from the existence of a single scaling form with time, it also confirms the collective nature of the local dynamics of our model systems.
As it is not straightforward to establish an analytic relation between these two observables, we adopt the standard statistical analysis route to do the same.

Relation between mobility and neighbourhood fluctuation
Joint probability distribution is a natural starting point to explore any interrelation between two given random variables. P(ū 2 i , χ), joint distribution of ū 2 i and χ, computed at τ w = 10τ for the two sets of studied morphologies, compact and stringlike aggregates, are presented in figures 5(a) and (b) respectively. From this data, it is visually evident that any small change in particle mobility ū 2 i results in a slight deformation of its nearest neighbour shell, quantified by χ, and vice versa. However, the shape of the distributions suggests that the very nature of interdependence is dictated by aggregates' morphology. To quantify the time evolution of this dependence, we compute the covariance and correlation of the two observables at different observation time τ w . By denoting the expectation value of a random variable (r.v.) x by E(x), the covariance, defined as cov( , provides information of a particular kind of dependence, namely, the degree of linear relationship between two observables. The strength of relationship is provided by the correlation, defined where the variance of a r.v. x is denoted by Var(x). We note that a nonzero definite value of either covariance or correlation is indicative of a linear relationship between two r.v.s but a zero value does not exclude the possibility of nonlinear dependence between them. Following Cauchy-Schwarz inequality, covariance is bounded by the standard deviations of individual r.v.s and correlation, by definition, is limited between values +1 and −1 which are only attained for the cases of strict linear dependence between two r.v.s.
Non-zero values of both the covariance and correlation for all τ w , as plotted in figure 5(c) and its Inset, confirms a strong dependence between ū 2 i and χ over all observation times. For small τ w , cov(ū 2 i , χ) is small and it grows monotonically with increasing τ w . We consider this a signature of emergence of collective dynamics within the system. We recall that after time τ , particles in both of our model aggregate forming systems fall out of the ballistic regime as suggested by the appearance of an intermediate plateau in the mean square displacement data plotted with respect to time. However, over this short timescale, the particles' movement and their fluctuations measured by Var(ū 2 i ) are expected to be very small due to low temperature and inherent frustration originating from competing interactions. The neighbourhood also remains mostly unchanged as the particles only start to interact with their neighbours after this time leading to small values of Var(χ). Consequently, the correlation C(ū 2 i , χ), defined as the covariance over individual variances, shows a large value for τ w ∼ τ. As time progresses, the collective dynamics sets in leading to increase in both covariance and individual variances of observables which results in the decay of correlation within τ < τ w < 10τ . For τ > τ w , C(ū 2 i , χ) reaches a steady state value denoting the stabilisation of collective dynamics, namely, caging for compact aggregates and binding for noncompact aggregates as revealed by our earlier study. We can then safely consider τ w > 10τ as the long time limit for our is plotted from (a) compact and (b) string-like aggregates for two different τ w . The stretched exponential nature of the distribution is shown for both cases with the solid line for τ w = 50.0. The exponential tail, shown by the dashed line, is more prominent for the compact case compared to its string-like counterpart. (c) Shows the collapse of P(χ) upon scaling with τ w for these two cases separately. Interestingly, the scaling form and exponents are exactly the same as used for P(ū 2 i ) data collapse (see text).
purpose. Within this limit, cov(ū 2 i , χ) grows linearly with time for compact clusters and the same shows sub-linear time evolution for non-compact aggregates. The specific reason and implication of such dependence, however, require further detailed investigation which is underway and will be reported elsewhere.
Finally, we attempt to characterise the interdependence between the motion of individual particles and the fluctuation in their neighbourhood caused by the motion of their nearest neighbours. It is easy to imagine instantaneous positional rearrangement of particles that constitute an aggregate due to thermal fluctuations, while the overall size and shape of the aggregate remain preserved over long observation time. It is then natural to consider the neighbourhood fluctuation, quanti fied by χ, as a response to the individual particle motion, characterised by ū 2 i . A functional relation between these two observables is obtained by computing the conditional expectation E(χ|ū 2 i ) as the average value of χ for a given value of ū 2 i . Considering these two observables as discrete r.v.s, formal definition reads as E(χ|ū 2 i ) = χ χf (χ|ū 2 i ) where the conditional probability density f (χ|ū 2 i ) = P(ū 2 i , χ)/P(χ) for the non-zero probabilities of the respective quantities. When these two observables are uncorrelated with each other, E(χ|ū 2 i ) should return only the value of E(χ).  i is shown for a range of τ w for the case of (a) compact and (b) stringlike aggregates. Inset shows the interdependence between the expectation of these two quantities. Now, given the knowledge of ū 2 i , all deterministic functions of ū 2 i are ideally known since any such function should behave as a constant in terms of the conditional expected value of the same variable. Among all these functions, it can be shown that the conditional expectation is the best predictor of the response as it minimises the mean square error of prediction. This is one basic result of statistics which is valid regardless of the type of the distribution. A non-linear dependence between χ and ū 2 i is suggested by their mean values plotted against each other in figure 6(Inset). The conditional expectation, E(χ|ū 2 i ) is found to show a power-law behaviour with respect to the mobility over unit observation time, ū 2 i /τ w for their small values and drops down after a maximum value which is dependent on the waiting time of observation. This observation prompts us to conclude the following simple functional relation: χ ∼ (ū 2 i /τ w ) β with the exponent β following the fractal dimension d f of aggregates very closely. To be specific, β values extracted from fitting the data shown in figures 6(a) and (b), are 1.67 and 0.88, respectively, for compact clusters with ideal d f = 2 and string-like clusters with ideal d f = 1. This functional form might be considered as the best predictor of response χ for any small change in ū 2 i at least in regression sense. While in situ measurements of individual particle mobility is possible, tedious post-processing is often required for quantification of non-affine fluctuations. The functional relationship presented here might speed up the process and help us to better speculate the trend of mechanical behaviour of a system on the fly.

Discussion and conclusion
In summary, we have analysed the simulated particle trajectories of a 2D aggregate forming system for its two different morphologies, namely, compact and string-like aggregates over several different lengths of observation time. While correlated particle movements are expected in such non-equilibrium systems, this is confirmed by our analysis and the very nature of this correlation is characterised in detail. The mobility ū 2 i of particles, quantified by the sum of squared displacement over certain observation time, is found to follow a log-normal distribution for string-like aggregates and is scale-free for compact clusters. Irrespective of this morphology dependent difference in the nature of ū 2 i distribution, time evolution of both distributions follow a common power-law form with observation time with the power-law exponent being close to two and one, the ideal fractal dimension of the compact and stringlike aggregates, respectively. Structural change of the nearest neighbour shell due the individual particle motion, measured by non-affine parameter χ, is observed to evolve over certain observation time in exactly the same way as ū 2 i , although the χ distributions are found to follow a stretched exponential form for all observation times and for both compact and string-like aggregates. This scaling form of the time evolution followed by both ū 2 i of a single particle and χ, resulting from all its neighbouring particles is a definitive proof of the correlated nature of local dynamics of the aggregate forming systems. This is one central result of our study which has not been discussed earlier in literature to our knowledge. The fact that the scaling exponents are the same as the fractal dimension of the respective aggregates is remarkable as it suggests an intimate relation between the local structure, specifically geometry, and the local dynamics. A detailed study of the trajectory geometry might be illuminating in this regard and we plan to carry out the same in future.
In recent years, χ has been proved to be a useful tool to understand the mechanical response of both crystalline and amorphous media with and without external perturbation [36]. The mobility, ū 2 i , on the other hand, is easily measurable [37,38] and widely used to study the heterogeneity in local dynamics of a variety of glass-forming systems [39,40]. This study provides a predictive phenomenological relation between these two observables from detailed statistical analysis of their distributions and thus, taking care of all fluctuations of these variables. Most importantly, the common finite-time scaling law followed by both observables, as revealed by the present study, gives us the opportunity to anticipate the long-time mechanical response from a very short-time measurement of particle mobility alone. As these findings are readily applicable to a variety of soft systems, as already mentioned earlier in the text, they might be crucial to design new engineered soft materials with tuneable mechanical properties.