Comments on “The proper view of cluster free energy in nucleation theories” by Runlong Cai and Juha Kangasluoma

Graphical Abstract

Copyright © 2022 American Association for Aerosol Research

Editor:

• Jim Smith

In a recent article by Cai and Kangasluoma (2022), the authors discuss the role of free energy and its interpretation in first-order phase transitions, i.e., nucleation. Their analysis is mainly restricted to basic thermodynamic principles: equilibrium conditions and the second law of thermodynamics in particular. The presented key argument is that the free energy barrier, $\Delta G_{\mathrm{barrier}},$ in the classical nucleation framework (separating the two phases) is fundamentally different from the reaction Gibbs free energy of nucleating clusters, ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}.$ The former being positive and the latter negative. By acknowledging this difference, within the suggested theoretical view the second law of thermodynamics is satisfied and the system is able to evolve spontaneously as the free energy along the transition coordinate is always negative, i.e., the entropy of the system increases. The presented analysis by Cai and Kangasluoma, however, suffers from severe mischaracterizations of both the standard theoretical framework of nucleation and the general nature of the phenomenon.

First, unlike it is stated on multiple occasions in their article, nucleation is not an equilibrium (or even a “quasi-equilibrium”) process. A phase transition designates the jump from a metastable state A to a more stable state B. Due to the free energy difference, $G_{\mathrm{B}}<G_{\mathrm{A}},$ the process is spontaneous. In the proper sense of nucleation, the path between states A and B, however, goes through intermediate states with free energies above the initial phase, hence the barrier. To surmount this barrier, a series of “lucky” out-of-equilibrium fluctuations is required. In the classical view of nucleation, these fluctuations take place in the cluster size space, but the system can also undergo perturbations due to e.g., rapid depletion around the nucleation site (Peters 2011) or deviations in cluster temperatures (Tikkanen et al. 2022). Thus, despite its spontaneous nature, once clusters surmount the free energy barrier, nucleation is a rare event facilitated by local serendipitous breakdowns of the metastable equilibrium, which are in fact local nanoscale violations of the second law. As suggested by the most common nucleation theories, the frequency of observing a nucleation event depends mostly on the free energy barrier height and the characteristic time scales of the system.

The factual breakdown of the metastable equilibrium in nucleation notwithstanding, the concepts adopted from equilibrium thermodynamics have been vital in understanding the fundamentals of the process and developing theoretical models of clustering in a supersaturated medium. In chemistry, free energies are often reported for specific reactants and products (computed at some standard pressure, marked here by the superscript $\theta$), and the favorability of a reaction is interpreted based on the free energy difference ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }.$ This difference is directly related to the relative abundance of reactants and products in an equilibrium state where no net flux exists. Following the notation used by Cai and Kangasluoma (2022), the equilibrium state between two adjacent cluster sizes can be expressed kinetically as (1) $${\beta }_{n-1}{P}_1^{\mathrm{eq}}{P}_{n-1}^{\mathrm{eq}}-{\alpha }_n{P}_n^{\mathrm{eq}}=0,$$(1) where $\beta$ and $\alpha$ are the collision and evaporation rate coefficients, respectively. Note that here superscript eq refers to any hypothetical zero-flux state, and thus pressures $P_n^{\mathrm{eq}}$ are not constrained to correspond to cluster population in a saturated equilibrium vapor.

As customary especially in studies of atmospheric new particle formation, the evaporation rate can be expressed using the standard free energies (given by Equation (3)(3) $${\beta }_{n-1}{P}_1{P}_{n-1}-{\alpha }_n{P}_n=J.$$(3) in Cai and Kangasluoma [2022]): (2) $${\alpha }_n={\beta }_{n-1}{P}^{\theta }\hspace{0.17em}\text{exp}\hspace{0.17em}\left(\frac{{\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n{A}_1\to A_n\right)-{\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left((n-1)A_1\to A_{n-1}\right)}{RT}\right).$$(2) In particle formation processes like nucleation, a non-zero flux $J$ exists and in a steady state this flux is constant throughout cluster sizes, (3) $${\beta }_{n-1}{P}_1{P}_{n-1}-{\alpha }_n{P}_n=J.$$(3) It is important to note that here $P_n$ are the actual “non-equilibrium” partial pressures observable during a process of steady-state nucleation whereas $P_n^{\mathrm{eq}}$ are quantities constrained to equilibrium (direct measuring of $P_n^{\mathrm{eq}}$ is virtually impossible).

Instead of associating $P_n$ and $P_n^{\mathrm{eq}},$ in nucleation theories, the nucleating clusters are just assumed to have their internal equilibrium properties, e.g., shape and molecular structure. According to this assumption, the evaporation rate $\alpha$ is a size-specific constant (depending on temperature, but not pressure), and therefore Equations (2)(2) $${\alpha }_n={\beta }_{n-1}{P}^{\theta }\hspace{0.17em}\text{exp}\hspace{0.17em}\left(\frac{{\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n{A}_1\to A_n\right)-{\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left((n-1)A_1\to A_{n-1}\right)}{RT}\right).$$(2) and (3)(3) $${\beta }_{n-1}{P}_1{P}_{n-1}-{\alpha }_n{P}_n=J.$$(3) can be rearranged as (4) $$\delta {\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n\right)={\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n{A}_1\to A_n\right)-{\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left((n-1)A_1\to A_{n-1}\right)=-RT[\ln \left(P^{\theta }{P}_n\right)-\ln (P_1{P}_{n-1}-\frac{J}{{\beta }_{n-1}})].$$(4) Now, by assuming that the cluster growth is purely monomeric and the system has reached a steady state, one can estimate the reaction free energies directly from the observed cluster distributions using the recursiveness of Equation (4)(4) $$\delta {\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n\right)={\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n{A}_1\to A_n\right)-{\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left((n-1)A_1\to A_{n-1}\right)=-RT[\ln \left(P^{\theta }{P}_n\right)-\ln (P_1{P}_{n-1}-\frac{J}{{\beta }_{n-1}})].$$(4) : ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n{A}_1\to A_n\right)={\sum }_{i=1}^n\delta {\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(i\right).$ By further setting $P^{\theta }=P_1^{\mathrm{eq}}=P_1$ and ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(A_1\to A_1\right)=0,$ one obtains the size-specific reaction free energies, ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}},$ corresponding to the conditions of the initial metastable state.

Equation (4)(4) $$\delta {\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n\right)={\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n{A}_1\to A_n\right)-{\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left((n-1)A_1\to A_{n-1}\right)=-RT[\ln \left(P^{\theta }{P}_n\right)-\ln (P_1{P}_{n-1}-\frac{J}{{\beta }_{n-1}})].$$(4) is often used to analyze obtained data from molecular-level computational simulations (Yasuoka and Matsumoto 1998; Halonen et al. 2021) where tracking even the smallest clusters during nucleation can be carried out rather effortlessly. As a brief demonstration of the actual relation between observable cluster pressures and the reaction free energies, we have performed a classical molecular dynamics (MD) simulation for nucleation in a supersaturated Lennard-Jones (LJ) vapor at a temperature of 0.6 and number density of 0.01 (in the standard LJ units). Such an atomistic simulation acts as a realistic “numerical experiment” as the system only obeys Newtonian dynamics and no external constraints involving either cluster formation or free energies are applied. The simulated system (with a fixed volume) consists of 500,000 nucleating atoms and 50,000 thermalizing carrier gas atoms (the details of the simulations can be found in Tikkanen et al. [2022]). Due to the large size of the system and slow enough nucleation, a sufficiently long steady-state period (at which $J$ is constant) is achieved at nearly constant monomer pressure, $P_1\approx 4.3\times {10}^{-3}$ (see Figure 1a). During this period the average ideal cluster pressures, $P_n,$ are calculated (Figure 1b) and the nucleation rate $J$ is determined using the Yasuoka-Matsumoto threshold method (Yasuoka and Matsumoto 1998) (Figure 1c). Knowing both $P_n$ and $J,$ and estimating $\beta$ using the kinetic gas theory for hard spheres, the reaction free energies can be readily computed using Equation (4)(4) $$\delta {\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n\right)={\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n{A}_1\to A_n\right)-{\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left((n-1)A_1\to A_{n-1}\right)=-RT[\ln \left(P^{\theta }{P}_n\right)-\ln (P_1{P}_{n-1}-\frac{J}{{\beta }_{n-1}})].$$(4) . The obtained ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}$ are shown in Figure 1d for the nucleating Lennard-Jones clusters. The reaction free energies clearly represent the classical parabolic nucleation barrier without any ambiguity. The found values characterizing the critical cluster (i.e., the position of the peak of the reaction free energy profile) are in good agreement with the nucleation barrier calculated at constrained equilibrium using semi-grand canonical Monte Carlo sampling (Halonen, Zapadinsky, and Vehkamäki 2018). The presented inspection of the simulated non-equilibrium steady-state cluster pressures via Equation (4)(4) $$\delta {\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n\right)={\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n{A}_1\to A_n\right)-{\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left((n-1)A_1\to A_{n-1}\right)=-RT[\ln \left(P^{\theta }{P}_n\right)-\ln (P_1{P}_{n-1}-\frac{J}{{\beta }_{n-1}})].$$(4) reveals no difference between ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}$ and $\Delta G_{\mathrm{barrier}}.$

Figure 1. Results from the molecular dynamics simulation of Lennard-Jones (LJ) nucleation. (a) Observed pressures of selected cluster sizes n as a function of time t. (b) Non-equilibrium cluster pressures (red circles) averaged over the steady-state nucleation period: from t = 2500 to 5500 (indicated by gray areas in panels (a) and (c)). The dashed black line shows the corresponding equilibrium cluster pressure curve. (c) Total pressure of clusters above threshold size n = 35 (which should be sufficiently larger than the critical size). The nucleation rate J (here reported in units of pressure/time instead of more conventional clusters/volume/time) corresponds to the slope of the curve at the steady state (dashed red line). (d) Obtained reaction free energies from the steady-state cluster pressure distribution using Equation (4)(4) $$\delta {\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n\right)={\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left(n{A}_1\to A_n\right)-{\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}^{\theta }\left((n-1)A_1\to A_{n-1}\right)=-RT[\ln \left(P^{\theta }{P}_n\right)-\ln (P_1{P}_{n-1}-\frac{J}{{\beta }_{n-1}})].$$(4) (large red circles). The small black circles show the reaction free energies according to the analysis by Cai and Kangasluoma (2022) (Equation (9) in their work). The found critical cluster size, n = 21, and the barrier height, 12.4RT (indicated by a vertical dotted line) are consistent with the data from previous Monte Carlo simulations (Halonen, Zapadinsky, and Vehkamäki 2018). Note that the relevant quantities are presented in the standard LJ units.

The crucial confusion in the free energy analysis by Cai and Kangasluoma (2022) is between the assumed (“quasi-”)equilibrium and the actual nucleation state, and this is overall manifested in their results due to assigning an improper equivalence between $P_n^{\mathrm{eq}}$ and $P_n.$ Specifically, in their article, ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}$ is effectively by mistake set to zero in Equation (7) where the pressures used inside the logarithm should in fact be the equilibrium pressures. While further analyzing the observed/simulated cluster pressures with Equation (9) in Cai and Kangasluoma (2022) and the classical equilibrium distribution function, $\Delta G_{\mathrm{barrier}}=-RT\ln \frac{P_n^{\mathrm{eq}}}{P_1},$ the suggested relation by Cai and Kangasluoma would thus be (5) $${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}\left(n{A}_1\to A_n\right)=RT(\ln \frac{P_n}{P_1}-\ln \frac{P_n^{\mathrm{eq}}}{P_1})=RT\ln \frac{P_n}{P_n^{\mathrm{eq}}}.$$(5)

This erroneous relation quite comprehensively explains the findings presented by Cai and Kangasluoma (see, e.g., Figures 2 and 3 in Cai and Kangasluoma [2022]): At a stable equilibrium state (e.g., unsaturated vapor), Equation (5)(5) $${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}\left(n{A}_1\to A_n\right)=RT(\ln \frac{P_n}{P_1}-\ln \frac{P_n^{\mathrm{eq}}}{P_1})=RT\ln \frac{P_n}{P_n^{\mathrm{eq}}}.$$(5) will result in apparent zero values for ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}$ since by definition $P_n=P_n^{\mathrm{eq}}$ for all values of n. This is also the case for the smallest, pre-critical clusters in a supersaturated vapor, as shown in Figure 1d, but ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}$ starts to diverge from zero as n reaches the critical size. At the critical size, according to theory (Frenkel 1955), the fraction $\frac{P_n^{\mathrm{eq}}}{P_n}$ has a universal value of 2, and indeed the ratio found from the MD simulation is about 1.9. For the over-critical clusters, the analysis by Cai and Kangasluoma will result in quickly decreasing negative values of ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}$ since $P_n^{\mathrm{eq}}\to \infty$ as $n\to \infty$ whereas $P_n$ is expected to decrease monotonically while being finite. This fundamental difference between $P_n^{\mathrm{eq}}$ and $P_n$ is further presented in Figure 1b, where the equilibrium (given by $\Delta G_{\mathrm{barrier}}$) and actual steady-state cluster distributions are depicted by the dashed line and open circles, respectively. Regarding any nucleating system, the considered quantities $P_n$ and $P_n^{\mathrm{eq}}$ are not equivalent, and should not be treated as such, because in the hypothetical equilibrium (described by the reaction free energies) the clusters are considered to be isolated from each other and thus the mutual action between different clusters, essential for the nucleation process, is disregarded (Frenkel 1955). In essence, the presented view of cluster free energies and Equation (9) in particular by Cai and Kangasluoma (2022) are thermodynamically inappropriate and do not challenge the general idea (dating back nearly a century) about the role or the nature of free energy in nucleation.

Yet, it is true, the standard reaction free energies are not directly compatible with the free energy barrier used in the classical nucleation theories (Elm et al. 2020). However, the nature of this incompatibleness is merely technical rather than rooted in some fundamental misconception related to the underlying physics. Recently, Halonen (2022) analyzed the apparent inconsistency between $\Delta G_{\mathrm{barrier}}$ and the routinely reported thermochemical values of ${\Delta }_{\mathrm{r}}{G}_{\mathrm{m}}.$ It was demonstrated that with adequate analysis of the kinetic framework behind classical nucleation theory one can easily calculate the steady-state nucleation rates using the standard free energy data obtained from, e.g., quantum chemical calculations.

To conclude, in equilibrium, all transitions between states are balanced (Equation (1)(1) $${\beta }_{n-1}{P}_1^{\mathrm{eq}}{P}_{n-1}^{\mathrm{eq}}-{\alpha }_n{P}_n^{\mathrm{eq}}=0,$$(1) ) and macroscopically the system stays unchanged; in non-equilibrium, this balance breaks down allowing the system to evolve at some finite rate (Equation (3)(3) $${\beta }_{n-1}{P}_1{P}_{n-1}-{\alpha }_n{P}_n=J.$$(3) ). Despite the debates over possible violations of the second law of thermodynamics in out-of-equilibrium processes at the nanoscale (Evans and Searles 1994; Jarzynski 2011; Spinney and Ford 2012), it is worth emphasizing that in nucleation entropy is not only related to mass transfer but also heat transfer which should be taken into account when assessing the total entropy production of the process. Inevitably, every first-order phase transition involves some degree of latent heat release and this excess energy is eventually dissipated to the surroundings in order to keep the system in thermal equilibrium (in the presented MD simulation, this specific process is emulated by using an explicit carrier gas connected to a heat bath via a thermostat). Thus, to assess the entropy production, the thermodynamic inspection of the nucleating system should not be limited to the formed clusters but to include all the relevant processes initiated by clustering.