Finite-size scaling of Lennard-Jones droplet formation at fixed density

We reaccess the droplet condensation-evaporation transition of a three-dimensional Lennard-Jones system upon a temperature change. With the help of parallel multicanonical simulations we obtain precise estimates of the transition temperature and the width of the transition for systems with up to 2048 particles. This allows us to supplement previous observations of finite-size scaling regimes with a clearer picture also for the case of a continuous particle model.


Introduction
Despite the long-lasting scientific interest in droplet formation, it is still a modern problem with many open questions. This is partially due to the general nature of droplet formation, with relevance ranging from metastable decay over cloud formation to cluster formation in protein solutions. In this work, we consider equilibrium droplet formation which yields a firm basis to study the transition between a homogeneous gas and a mixed phase of a droplet in equilibrium with surrounding vapor [1,2,3,4,5]. The problem is formulated in the canonical ensemble. A common control parameter is the density while the temperature is fixed. This allows one, in principle, to access temperature-dependent quantities such as the isothermal compressibility and the interface tension.
At fixed temperature the theory has been verified by several computational studies, including lattice [6,7,8] and Lennard-Jones [9,10,11] systems. Instead, however, one can fix the density and vary the temperature which yields an equivalent but "orthogonal" finite-size scaling behavior [12,13]. In this case, we recently observed an intermediate scaling regime consistent with finite-size scaling results for polymer aggregation [14]. Here, we present new results for the three-dimensional Lennard-Jones system at fixed density which extend and supplement our previous results [13].

Model and Method
We consider N Lennard-Jones particles in a three-dimensional box of length L with periodic boundary conditions, see Fig. 1. The self-avoidance and short-range attraction is modeled by the pairwise Lennard-Jones potential where r ij is the distance between particle i and j. As in Ref. [13], we set = 1 and σ = 2 −1/6 . The computational demand can be reduced by introducing a cutoff radius r c = 2.5σ above which particles do not interact anymore. The potential is then shifted by V LJ (r c ) in order to be continuous, yielding This is in accordance with the existing literature and enables the application of a domain decomposition, where the periodic box is decomposed into equally large (cubic) domains. These domains have to be at least of the size r c . Then, the interaction of each particle is obtained by evaluating only its domain and the adjacent ones (in three dimensions this adds up to 3 3 = 27 domains). Especially in the gas phase the simulation benefits from this procedure, where the particles are equally distributed in the full box. Following Ref. [13], we fix the density ρ = N/L 3 and vary the temperature T . This allows the application of multicanonical simulations [15,16,17,18], which are well-suited for firstorder phase transitions such as condensation. At a first-order transition two phases are in coexistence with suppressed transition states in between. This is circumvented by replacing the canonical Boltzmann weight exp(−E/k B T ) with an alterable weight function, which is iteratively adapted in order to yield a flat histogram in the energy. Each iteration is in equilibrium, sampling the distribution according to the current weight function. This leads to a straightforward parallelization [19,20] which has been shown to perform very well for the problem at hand [8]. In the end, canonical expectation values are estimated by reweighting the data from a multicanonical production run.
A crucial aspect is the selection of Monte Carlo updates. While in Ref. [13] we restricted ourselves to local particles shifts, we added here particle "jumps" with a larger update range. The update range for a particle shift was set to r = 0.3 and for a particle jump to r = L/2. A new position is then selected with equal probability from a sphere of size r around the old position. This simple non-local update significantly enhanced the sampling of the gaseous phase and increased the performance of the simulation drastically. Instead of 512 we are now able to sample up to 2048 particles.

Theory
Here we only want to briefly recapture the results of many previous works [2,3,4,5,13]. For a supersaturated gas, it was shown that the probability of intermediate-sized droplets vanishes [3,4], which leaves us with the scenario of a homogeneous gas phase on the one side and a droplet in equilibrium with surrounding vapor on the other side, see Fig. 1. At a fixed temperature, the droplet formation may be achieved by adding more particle excess to the already supersaturated gas. For both sides it is possible to formulate a contribution to the free-energy in terms of fluctuations (entropy of the gaseous phase) and surface tension (energy of the droplet). Importantly, the finite-size dependence could be rewritten in terms of the fraction of particle excess in the droplet as a function of a dimensionless density. This in turn allowed us to expand the results around the infinite-size transition temperature T 0 to yield the leading finite-size scaling behavior of the condensation temperature T c at fixed density [13]. Similarly, we showed that an expansion of the free-energy difference at the condensation transition yields an estimate of the leading-order finite-size scaling of the transition rounding ∆T , i.e., the width of the transition. For details we refer to the prior literature and give here only the three-dimensional results to leading order: A crucial observation is now that the size R of the droplet at the condensation transition itself scales non-trivially with system size, namely in three dimensions R ∝ N 1/4 . Thus, the leading-order scaling may be identified in terms of powers of the linear extension of the droplet itself, ∝ R −1 respectively ∝ R −3 . This is consistent with the interpretation that the droplet size at transition is the relevant system size. Then, a virtual subsystem around the droplet would lead to a transition between a homogeneous liquid (droplet) phase to a homogeneous gas phase with open boundary conditions. The competition between finite-size contributions from volume (∝ R 3 ) and surface (∝ R 2 ) would in this case give rise to an intuitive finite-size correction of the order R −1 [21,22,23,13,14].

Results
Due to the non-local Monte Carlo update, we are now able to extend our previous results [13] for the three-dimensional Lennard-Jones system to N = 2048 particles with improved statistics. The temperature scale is fixed by setting k B = 1. The finite-size transition temperature T c is obtained as the location of the largest peak in the specific heat C V = k B β 2 E 2 − E 2 /N and plotted in Fig. 2 (left) versus the expected scaling behavior Eq. (3). Overall, the data qualitatively shows a linear behavior as predicted. A leading-order fit for N ≥ 1024 yields T 0 = 0.7032(3) with goodness-of-fit parameter Q ≈ 0.1, shown in the figure as dashed black line. Including additional higher-order corrections, i.e., T c = T 0 + aN −1/4 + bN −2/4 , yields T 0 = 0.6979 (2) for N ≥ 64 with Q ≈ 0.1. Both estimates differ from our previous results outside the error bars, which is expected because the fit error underestimates systematic uncertainties. This is partially due to the fact that we are not yet fully in the asymptotic scaling regime, which we will discuss below. The general trend of the results is compatible with our previous conclusion, where we already noted that the available three-dimensional Lennard-Jones system sizes were too small for clear results. However, it is worth noting that the current estimates are closer to each other, implying an infinite-size limit of the transition temperature around T 0 ≈ 0.70.
The distance from the asymptotic scaling regime is best visible in the rounding of the transition. This is here obtained from the half width of the specific-heat peak, i.e., the width for which C V (T ) > C max V /2, shown in Fig. 2 (right). Since the width vanishes in the thermodynamic  Figure 2. Finite-size scaling of the transition temperature (left) and rounding (right) as the maximum location and the half width of the specific-heat peak, respectively. Error bars are included but smaller than the symbol size.
limit, a double-logarithmic plot reveals the power-law scaling Eq. (4) as a straight line, here dashed black. As previously noticed [13], an intermediate scaling regime ∝ N −1 is clearly visible and marked with a dashed-dotted gray line. This is consistent with the droplet as relevant system size as in this regime still a large fraction of particles contribute to the transition droplet. However, for N ≥ 500 it appears that the asymptotic scaling behavior slowly manifests itself and leading-order estimates become more reliable. This can be tested by performing a direct fit of the power-law ansatz ∆T ∝ N −α , shown in Fig. 3 for variable lower bound N min . Only for N min ≥ 640 the fit yields plausible goodness-of-fit values Q > 0.1. Still, the distance to the predicted asymptotic scaling regime is noticeable. This is consistent with results obtained for the three-dimensional lattice gas where decent results were only obtained for N min ≈ 2000 [13]. This emphasizes once more that care has to be taken with leading-order finite-size scaling away from the asymptotic scaling regime. While estimates of the thermodynamic limit become better with increasing system size, fit errors are difficult to judge because the underlying ansatz is not accurate enough. The same holds true for additional higher-order corrections. While the results are surely good estimates, the fit errors are a difficult measure, not capturing the uncertainty that comes from an incomplete fit model [24].

Conclusions
We have extended our previous results [13] for the three-dimensional Lennard-Jones system to larger system sizes. This allows for more consistent estimates of the thermodynamic limit when considering leading-order and higher-order fits. Still, the current Lennard-Jones system with up to 2048 particles at density ρ = 10 −2 is quite far away from the asymptotic scaling limit. We notice that the rounding of the transition is a good indicator to visualize this distance from the asymptotic scaling regime. Here, also the emerging intermediate scaling regime is best noticeable.