Scaling relations for finite-time first-order phase transition

Introduction.–Phase transition, traditionally lying in the framework of equilibrium statistical mechanics, is primarily concerned with static quantities. Equilibrium phase transitions are drastic changes of a system’s equilibrium state under quasi-static variation of the external parameters (Callen and Griffiths, 1987; Goldenfeld, 1992). At a phase transition, a system undergoes a nonanalytic change of its properties, for example the density in a temperature-driven liquid-gas transition, or the magnetization in a paramagnet-ferromagnet transition (Chaikin and Lubensky, 1995; Imparato and Peliti, 2007). These are static properties with respect to the external parameters. In the past few decades, owing to the recent development of nonequilibrium thermodynamics (Graham and Tél, 1984; Keizer, 1987; Imparato and Peliti, 2007; Seifert, 2012; Sagawa, 2014; Bertini et al., 2015; Pineda and Stamatakis, 2018; Peliti and Pigolotti, 2021), there is a growing interest in the field of stochastic thermodynamics with phase transitions (Derrida, 1987; Garrahan et al., 2007; Mehl et al., 2008; Lacoste et al., 2008; Jack and Sollich, 2010; Ge and Qian, 2010; Herpich et al., 2018; Nyawo and Touchette, 2018; Byrd et al., 2019; Noa et al., 2019; Fei et al., 2020; Vroylandt et al., 2020; Proesmans et al., 2020a; Nguyen and Seifert, 2020; Meibohm and Esposito, 2022). A central concept in thermodynamics is the work performed on such a system (Speck and Seifert, 2004, 2005; Blickle et al., 2006; Quan et al., 2008; Engel, 2009; Nickelsen and Engel, 2011; Ryabov et al., 2013), when some external parameters of the system are varied with time. The second law of thermodynamics constrains that the work $W$ performed on the system is no less than the free energy difference $\Delta F$ between the initial and the final states. In other words, the excess work done in a finite-rate quench with respect to the quasi-static quench process is non-negative $W-\Delta F\geq 0$. A natural question is how thermodynamic quantities such as the excess work scale with the rate of quench.

In the regular cases without phase transitions, it is well-known that in slow isothermal processes, the excess work is proportional to the quench rate as a result of the linear response theory (Salamon and Berry, 1983; Mazonka and Jarzynski, 1999; Crooks, 2007; Scandi and Perarnau-Llobet, 2019; Ma et al., 2020; Chen et al., 2021). While in slow adiabatic processes, the excess work is proportional to the square of the quench rate according to the adiabatic perturbation theory (Sun, 1988; Polkovnikov, 2005; Rigolin et al., 2008; De Grandi and Polkovnikov, 2010; Polkovnikov et al., 2011; Chen et al., 2019). However, in the presence of phase transitions, the scaling behavior becomes much more involved. For the second-order phase transition, it is found that the excess work done during a quench across the critical point exhibits power-law behavior with the quench rate and the corresponding exponents are fully determined by the dimension of the system and the critical exponents of the transition (Fei et al., 2020; Zhang and Quan, 2022), as in the traditional Kibble-Zurek mechanism (Kibble, 1976, 1980; Zurek, 1985, 1996; del Campo, 2018; Cui et al., 2020). Besides the second-order phase transition, the first-order phase transition is ubiquitous in nature, e.g., in biochemical (Ge and Qian, 2009, 2010; Qian et al., 2016; Pineda and Stamatakis, 2018; Nguyen and Seifert, 2020), ecological (Li et al., 2019) and electronic systems (Freitas et al., 2022). Examples include the coil-globule transition in the process of protein folding (Anfinsen, 1972; Kaiser, 1978; Qian et al., 2016) and the information erasure in electronic devices (Koller and Athas, 1992; Hanninen and Takala, 2010; Bérut et al., 2012; Proesmans et al., 2020b). However, the scaling behavior of the excess work during a finite-time first-order phase transition has not been explored so far.

In this Letter, we study the scaling behavior of the excess work with the quench rate in the first-order phase transition. Using the mean-field Curie-Weiss model with a varying magnetic field as a paradigmatic example (Fernández et al., 2012; Collet and Formentin, 2019; Meibohm and Esposito, 2022), we uncover the scaling relation of the excess work with the quench rate as $W_{\mathrm{ex}}\propto v^{2/3}$ for small quench rates in the finite-time first-order phase transition. This replenishes our understanding about nonequilibrium thermodynamics in driven systems with phase transitions. Besides, we find that both the delay time $\hat{t}_{\mathrm{del}}$ and the transition time $\hat{t}_{\mathrm{trans}}$ scale with the quench rate as $v^{-1/3}$. The scaling laws should be valid for a large class of systems with the first-order phase transition since they do not depend on the details of the systems, but are a manifestation of the nonlinear nature of their dynamics.

The model.–We consider the kinetic 1D mean-field Curie-Weiss model with ferromagnetic interaction introduced in Ref. (Meibohm and Esposito, 2022). The Curie-Weiss model consists of $N$ Ising spins $\sigma_{i}=\pm 1$, labeled $i=1,\dots,N$, with the coupling strength $J/2N$. The system is embedded in a heat reservoir at an inverse temperature $\beta=1/(k_{\mathrm{B}}T)$, and subject to a varying field $H(t)$ controlled by some external agent. The state of the system is characterized by the magnetic moment $M=\sum_{i}\sigma_{i}$. In the thermodynamic limit $N\to\infty$, we use the mean magnetization $m\equiv M/N$ to denote the system state. The internal energy density is given by (Goldenfeld, 1992)

The deterministic dynamics is described by the equation of motion of $m(t)$ (see Supplemental Materials (Sup, ) for details):

with the microscopic relaxation time $\tau_{0}$.

The curve characterizing steady state as a function of $H$ is shown in Fig. 1, where we choose the parameters to be $\tau_{0}=2,J=1,\beta=6$. Depending on the value of $H$, the ordinary differential equation (2) can have one, two or three steady-state solutions. For small and large values of $H$, the system has only one steady state ,i.e., the system is monostable. For other intermediate values of $H$, there are three steady states: $m_{\mathrm{top}}$ (blue line), $m_{\mathrm{bot}}$ (green line), and $m_{\mathrm{mid}}$ (orange line), i.e., the system exhibits bistability. Among the three, one is globally stable; one is metastable and the one in between is unstable. The two black dots are the turning points $(H_{\pm},m_{\pm})$, where the metastability disappears.

During the quench process, the external field $H$ is tuned quasi-statically from $H_{i}=-1.5$ to $H_{f}=1.5$ and the transient state of the system depends on its history. The quasi-static hysteresis loop is shown with the pink dotted curve in Fig. 1. The quasi-static state $m_{\mathrm{qs}}$ gets trapped in the local minimum as $H$ is increased quasi-statically from $H_{i}$ to $H_{f}$, i.e., $m_{\mathrm{qs}}$ follows the lower branch $m_{\mathrm{bot}}$ even when the local minimum is metastable rather than globally stable, until the right turning point $(H_{+},m_{+})$. Then the system state jumps to the upper branch. This is what we refer to as “quasi-static first-order phase transition”. So much the same for the return path. The quasi-static state $m_{\mathrm{qs}}$ remains on the upper branch and jumps down to the lower branch at the left turning point $(H_{-},m_{-})$.

When $H(t)$ is tuned at a finite rate $v$, the system tries to keep pace with —but ultimately lags behind— the continually changing quasi-static state, as is shown by the gray dashed curve in Fig. 1. Although at time $0$, the system is prepared in the equilibrium state of $H_{i}$, i.e., $m(0)=m_{\mathrm{qs}}(H_{i})$, at latter times $m(t)\neq m_{\mathrm{qs}}(H(t))$. This mismatch becomes evident around the turning point.

We illustrate the finite-time first-order phase transition with a particle moving in an evolving potential landscape $U$ (where $U$ satisfies $f(m,H)=-\partial_{m}U/\tau_{0}$) in Fig. 2. Far away from the turning point (see Fig. 2(i-ii)), the particle is able to catch up with the change of the potential landscape and is trapped around the left local minimum, no matter whether it is globally stable or metastable. It becomes different around the turning point (see Fig. 2(iii)) where the system transitions from bistable to monostable and the left local minimum disappears. Thus the right minimum becomes the only attractor in the whole landscape. However, the evolution of the particle lags behind the drastic change of the potential. It takes some time for the particle to “feel the change of its surroundings” and make a transition to the right stable minimum. The non-catching-up obviously depends on the rate of quench. After reaching the new stable state, the particle gets stuck as the deforming of the potential landscape continues (Fig. 2(iv)).

Scaling relations of time.–During the quench process, the trajectory follows a stable local equilibrium, be it metastable or globally stable, until an abrupt transition into a totally different state from the original one. This transition relies on the system’s history and is regarded as a finite-time first-order phase transition. This is caused by the disappearance of metastability. The turning points satisfy

where $f(m,H)$ is defined in Eq. (2). The turning points can be calculated analytically as $(H_{\pm},m_{\pm})=(\pm H^{*},\pm m^{*})$, where

Here we consider the case in which the external magnetic field $H(t)$ is varied from the initial value $H_{i}$ to the final value $H_{f}$ according to the linear protocol

with the quench rate $v=(H_{f}-H_{i})/\tau$ and the time duration $\tau$ . We label the time when $H(t)$ reaches the right turning point $H^{*}$ as the turning time $t^{*}$ (Li et al., 2019).

It is important to note that while the potential landscape has significantly changed from bistable to monostable at the turning time $t^{*}$, in a finite-rate quench the dynamics of the system cannot catch up with the change of the potential landscape, i.e., the system displays a delay in transitioning to the monostable state. The transition occurs later than the turning time $t^{*}$ by a delay time $\hat{t}_{\mathrm{del}}$ defined through $m(t^{*}+\hat{t}_{\mathrm{del}})=m^{*},$ which depends on the quench rate $v$.

We further analyze the scaling relation between the delay time and the quench rate. We shift the variables in the following way: $\hat{H}=H-H^{*},\hat{m}=m-m^{*},\hat{t}=t-t^{*}$. The evolution equation is turned into

For slow quench, the equation of motion around the turning point can be approximated as (see Supplemental Materials (Sup, ) for details)

It can be adimensionalized into

by rescaling the variables

where

Under the asymptotic boundary condition $s\to-\infty,u\to\sqrt{-s},$ Eq. (9) is solved in terms of the Airy function (Jung et al., 1990)

In this way, the delay time can be computed as

where $A_{1}^{\prime}\approx-1.019$ represents the first zero point of the Airy Prime function. As is shown in Fig. 3(a), the analytical expression matches well with the numerical results.

Besides the delay time $\hat{t}_{\mathrm{del}}$, there is another time $\hat{t}_{\mathrm{trans}}$ characterizing the non-catching-up of the dynamics. When the system is quenched at a finite rate, the transition time $\hat{t}_{\mathrm{trans}}$ to reach the monostable state, which is defined through $m(t^{*}+\hat{t}_{\mathrm{trans}})=m_{\mathrm{nd}},$ depends on the quench rate $v$. Here $m_{\mathrm{nd}}$ stands for the non-degenerate root of $f(m,H_{+})=0$, with the expression

For $\beta J\gg 1$, $m_{\mathrm{nd}}$ can be approximated by $1$.

The transition time can be calculated as

where $A_{1}\approx-2.338$ stands for the first zero point of the Airy function (Li et al., 2019).

Equations (15) and (17) show that both the delay time and the transition time scale with the quench rate as $v^{-1/3}$. The matching between analytical expressions and numerical results for $\hat{t}_{\mathrm{del}}$ and $\hat{t}_{\mathrm{trans}}$ is illustrated in Fig. 3(a). The exponent $-$1/3 results from the quadratic leading order $\hat{m}^{2}$ around the turning point in the Curie-Weiss model.

With the information about the delay time, we can estimate the time we have to save the system with phase transitions from catastrophe when the metastability disappears in economic systems and climate science (Lenton et al., 2008; Li et al., 2019). And with the transition time, we are able to estimate the time it takes to reach a new stable state from the turning point.

Scaling relation of work.–In the process of a finite-rate quench, the work performed on the system is greater than that in the quasi-static quench process. In other words, the non-catching-up in the dynamics between the transient state and the quasi-static state results in the excess work, which characterizes the irreversibility of the process. It is desirable to explore the relation between the excess work and the quench rate.

From the microscopic definition of work (Jarzynski, 1997; Sekimoto, 2010), the work performed on the system is

In the quasi-static isothermal process, the work equals the free energy difference between the initial and final equilibrium state

The excess work $W_{\mathrm{ex}}$, defined as

describes the excess amount work that one has to perform when the system is quenched at a finite rate rather than quasi-statically. In a slow quench, the excess work can be approximated by the enclosed rectangular area between the dynamic and quasi-static hysteresis

where we approximate $m_{\mathrm{nd}}$ with $1$ in the second line. From the analytical expression in Eq. (21), we obtain that the excess work scales with the quench rate as $W_{\mathrm{ex}}\propto v^{2/3}$, which is verified numerically in Fig. 3(b).

Last but not least, we would like to compare our results in the first-order phase transition with the second-order phase transition and the case without phase transition in Table. 1. For finite-time isothermal processes without phase transition, the excess work $W_{\mathrm{ex}}$ due to non-quasi-static quench is proportional to the quench rate $v$ as a result of linear response theory (Salamon and Berry, 1983; Mazonka and Jarzynski, 1999; Crooks, 2007; Scandi and Perarnau-Llobet, 2019; Ma et al., 2020; Chen et al., 2021). For slow adiabatic processes, the excess work scales as $v^{2}$ according to the adiabatic perturbation theorem (Sun, 1988; Polkovnikov, 2005; Rigolin et al., 2008; De Grandi and Polkovnikov, 2010; Polkovnikov et al., 2011; Chen et al., 2019). In a second-order phase transition, the excess work performed during a quench across the critical point exhibits the $v^{\delta_{1}}$ scaling behavior (Fei et al., 2020; Zhang and Quan, 2022), is relevant to but different from the Kibble-Zurek scaling for the average density of defects (Kibble, 1976, 1980; Zurek, 1985, 1996). The scaling exponent $\delta_{1}$ is determined by the dimension of the system and the critical exponents of the transition. Here our work uncovers the scaling relation between the excess work and the quench rate in the first-order phase transition, thus replenishing our understanding for the finite-time transition. We emphasize that the scaling relation we obtain here is not restricted to the Curie-Weiss model, but should be universal for a large classs of systems with the first-order phase transition. It can also be derived from the finite-time Landau theory (Meibohm and Esposito, 2023) with a similar equation of motion of switching dynamics in bistable systems (Jung et al., 1990).

Conclusions.–In this Letter, we study the time and the work applied in the quench across a first-order phase transition in the Curie-Weiss model with a varying magnetic field. We uncover a $2/3$ scaling relation between the excess work and the quench rate in the finite-time first-order phase transition. This replenishes our understanding about the nonequilibrium thermodynamics in driven systems with phase transitions. Moreover, we show that both the delay time and the transition time exhibit a $v^{-1/3}$ scaling for small quench rates. The model examined here is generic, and the scaling relations should be universal for a large class of systems with the first-order phase transition. The scaling relations in Eqs. (15,17,21) do not depend on the details the models, but are a manifestation of the quadratic leading order $\hat{m}^{2}$ around the turning point. This may also give us some insights into the electronic computation schemes, where non-adiabatic effects impose a tradeoff between the speed and the energy consumption that can be optimized (Bérut et al., 2012; Freitas et al., 2022).

This work is supported by the National Natural Science Foundations of China (NSFC) under Grants No. 12375028 and No. 11825501, No. 12147157.