The KPZ equation for kinetic interface roughening (KIR) [1–4],

where h(x, t) is the interface height and ξ(x, t) a Gaussian noise with

is nowadays a paradigm of systems exhibiting non-equilibrium critical scaling [5]. In fact—besides standing out as a representative of a large and robust class of microscopic KIR models¹, from which the phenomenological parameters ν, λ, D can be computed—it is intimately related to two highly non-trivial problems:

Equation (2) is often attributed to the Burgers equation's Galilean invariance (GI), which in turn translates into KPZ equation's invariance under changes in tilt. This opinion has been repeatedly challenged and in fact, the numerically computed exponents obey (2) in a discrete version of (1) where both Galilean invariance and the (1d–peculiar) fluctuation–dissipation theorem are explicitly broken [11, 12]².

Since the first and third terms in the r.h.s. of (1) were already present in the model by Edwards and Wilkinson (EW), the innovation came from the tilt-dependent local growth velocity in the second term (or rather, by its interplay with the third one). Under additive uncorrelated local Gaussian white noises, plane (eventually moving) interfaces—which are stable for ξ(x, t) = 0 [14, 15]—develop (in both EW and KPZ models and for large enough systems) into statistically self-affine fractals, whose typical roughness width scales as a power β (called the growth exponent) of the elapsed time. β is in turn the ratio between the interface's Hölder (or roughness) exponent α and the dynamic exponent z, governing the growth in time of the correlation length. In the well understood d = 1 case, α turns out to be that of simple random walk (namely α=12) in both models. However, (2) imposes that z (and thus β) be appreciably different (at least after some crossover time, needed for the second term to take over the first). The initial success of (1) was then due to the consistency of β_KPZ with numerical results on microscopic KIR models [2–4].

By about half its lifetime so far, the field was mature enough to approach not only the second moment of the fluctuations in h, but its full statistics [16–20]. In the 1d case, by proposing for h(x, t) the asymptotic behavior

with A :=ν2D, Γ :=A2λ2, and x′ :=Ax2(Γt)-2/3, one can see that the stochastic variable χ obeys the Tracy–Widom (TW) statistics of the largest eigenvalue of a random-matrix ensemble, here determined by the dynamics of the substrate³. In the stationary state however, temporal correlations are governed by the Baik–Rains (BR) F₀ limit distribution [23], not related to random-matrix theory (RMT).

In the last few years, a handful of exact solutions to (1) have arisen: whereas most of them were inspired by RMT [18, 19, 24, 25], one was born right within the field of stochastic differential equations [26]. On the other hand, the field is by now so mature that experiments can decide on the statistics [20, 27–29]. And pretty much the same occurs with the numerics: in a recent review article [21], the statistics of KPZ itself has been found to agree with the results transposed from DPRM [16–20, 30, 31].

It is the purpose of this article to contribute to lighten up the abovementioned developments from the perspective of a recent variational formulation of the KPZ equation [11, 12, 32–34]. In the following section we derive general consequences from some properties of the free-energy-like functional Φ[h]. Next, we undertake a numerical study of the time dependence of Φ[h] (featuring an estimator of the EW-to-KPZ crossover time and its dependence on parameter λ) in d = 1, 2, and 3, and derive some consequences from the forced case. Then we relate the time behavior of Φ˙[h] with the asymptotic form in (3) and finally, we collect our conclusions.

The intuitive power of disposing of a potential landscape in the study of dynamical systems—as provided by energy in conservative systems, or free energy in equilibrium thermodynamics—and even in other fields (given its closeness to concepts such as “information” or “fitness function”), is beyond doubt. In gradient dissipative systems, it is immediate to find a Lyapunov function to play that role. In general relaxational ones, meeting the integrability conditions may be hard and moreover, the obtained potential landscape partially loses its intuitive power. A breakthrough proposal [35] to circumvent the first drawback—relying on Langevin (thus Fokker–Planck) dynamics having stationary states—was made three decades ago. Unfortunately, KIR is but an example of those cases where this proposal cannot be applied.

In the literature, it is commonly assumed that the KPZ equation cannot be directly obtained from a Hamiltonian. However, as shown in Wio [32], there is no trouble in expressing (1) as a stochastic gradient flow,

The functional Φ is defined as

where h₀(x, 0) is an arbitrary initial pattern (usually assumed to be constant, in particular h₀ = 0). The first term is clearly the Landau–Ginzburg free-energy functional

associated to the (equilibrium) EW process. Unfortunately, the second one has not an explicit density. Even though at time t, Φ depends only on the field h(x, t), its evaluation requires knowing the detailed history that led from h₀ to h(x, t). In other words, retrieving information on Φ (such as its landscape at certain time or the time dependence of its mean value) requires averaging not simply over field configurations h(x) at time t, but over (statistically weighted) trajectories of the field configuration.

Being the KPZ equation a stochastic gradient flow, the functional Φ governing its deterministic component provides the landscape where the stochastic dynamics of h(x) takes place at time t, and in the absence of noise fulfills explicitly the Lyapunov property Φ˙[h]=-[δΦ[h]δh(x,t)]2≤0. Unfortunately, this does not make Φ into a Lyapunov functional, since (as shown in the next paragraph) it is unbounded from below. As it is well-known, Equation (1) has the joint symmetry λ ↔ −λ, h ↔ −h. In the case of Equation (5), the change λ ↔ −λ—together with the (∇h)² dependence—implies a specular reflection h ↔ −h with respect to the h axis.

In fact—given a reference configuration h₀(x, t₀)—an expansion of Φ[h] in powers of δh(x, t) = h(x, t) − h₀(x, t₀) exactly terminates at the third order:

wherefrom

Upon discretization (in the 1d case) ΔΦ becomes a polynomial in the variables {δh_i},

with the “local potential” terms

and the definitions (∇h₀)_i : = (h₀)_i+1 − (h₀)_i−1, (∇2h0)i :=(h0)i+1+(h0)i-1-2(h0)i.

Note that δh_i > ν/λ would produce a negative parabolic potential, which implies local instabilities. Such a diffusive instability in [ΔΦ₂]_i tells us that the above procedure makes sense only for short times i.e., it provides an “instantaneous landscape” of the potential. However, this is all we need to extract important information on the qualitative aspects of the dynamics. On one hand—by associating the height δh_i at site i with the coordinate of a particle submitted to the potential [ΔΦ]_i—we retrieve a mechanical picture that helps understanding the local dynamics. We notice that [ΔΦ]_i is at most linear in δh_i. But whereas the coefficient of the first term is a constant, the other two terms are some kind of “mean-field potentials.”

Figures 1, 2 depict different contributions to [ΔΦ]_i, as functions of δh_i−1 and δh_i+1. Figure 1 shows that whereas for EW, [ΔΦ₂]_i has a concave parabolic shape with constant curvature, for KPZ this curvature depends on the sign of δh_i through the contribution from δ³Φ[h(x, t)] in (6). Figure 2 illustrates the breakdown of the “horizontal” symmetry δh_i+1 ↔ δh_i−1 by showing 3d plots of [ΔΦ₃]_i when (∇h₀)_i is either positive or negative.

Rather than the quantitative information these figures may provide, some essential aspects—as the local symmetries involved in each case and the dynamical consequences of their breakdown—are worth discussing. We focus on three facts from (7):

Note that in the EW case these symmetries hold, while in the KPZ all of them are broken. These facts lend themselves to insightful interpretations that we study in what follows.

Before delving however into the thorough analysis of Φ[h] to which this article is devoted, we find it convenient to introduce a genuine Lyapunov functional for KPZ, to see its pros and cons. The functional

of Wio [32], allows the KPZ equation to be written as

with Γ[h] given by

It is straightforward to verify that (i) the Lyapunov property F˙[h]≤0 holds, (ii) the EW free-energy functional is recovered for λ = 0, and (iii) 〈F[h]〉t behaves asymptotically as f(t)exp[λνv∞t], with f(t) far milder than an exponential [2, 3]. Moreover, numerical simulations show that v_∞ ∝ λ, so the asymptotic behavior results ~exp[λ²t].

Since (9) is not of gradient form, F[h] does not provide an intuitive landscape for the stochastic dynamics of h(x, t). So hereafter, we focus on the analysis of Φ[h]. As stated before, its form provokes several thoughts:

Spectral methods pay back the overhead of Fourier transforming h(x) forth and back every time t, with a discrete approximation of ∇h of the order of the system's size. They have proved to be more stable and reliable than pure real-space finite-differences schemes in the integration of some non-linear growth equations [40–42]. Using a spectral method to integrate (1), we have analyzed the time behavior of 〈Φ[h]〉_t (the average of Φ over noise realizations, at time t) as given by (14).

Figure 3 displays the time dependence of the NEP's average over 100 samples, for systems in 1d (size 1,024), 2d (size 128²), and 3d (size 64³), and several values of λ. For any λ > 0 there is a maximum, where the non-linear (KPZ) term overcomes the linear (EW) one. Only for λ = 0 (corresponding to EW) is there true saturation. The λ = 0.01 case shown here (visually resembling a saturation behavior) attains its maximum outside the plotted range.

Figure 4—as well as (16) and (17) below—shows that past this maximum, 〈Φ[h(t)]〉_t ~ A − Bt with B ~ λ². This result shows that—due to the correlations and in an effective way—the NEP behaves as having only a linear dependence on h (just as in the toy model discussed before [37, 38]).

Figure 5 shows the dependence on λ of the time t_xovr at which the maximum occurs, i.e., at which the non-linear (KPZ) contribution overcomes the linear (EW) one, making it a proxy of the EW-to-KPZ crossover time [43–45]. Although roughly compatible with a txovr~λ-α law (α = 3 in Guo et al. [43] and α = 4 in Forrest and Toral, [44] for d = 1) and in agreement with the results of previous studies exploiting the time behavior of the stochastic action [46, 47], this dependence becomes milder as d increases.

As known in the EW case [2–4] and shown in Figure 3, the average 〈 ∫ dx(∇h)²〉 over noise realizations saturates to a constant C₀ at long times. Using in (15) the discrete expression for Φ[h]—namely (12), with the approximation in (13b)—we can write

In the discrete representation, (3) implies δh ~ v_∞τ at long times, where τ corresponds to δt and v_∞ is the asymptotic front velocity. Inserting this expression into the previous one and replacing 〈 ∫ dx(∇h)²〉 by its bound C₀, we get

This result, valid in any dimensionality, completely agrees with the numerical ones presented above. In particular, 〈Φ[h]〉 will decrease at long times only if λC₀ + F ≥ 0. If λC₀ + F < 0, 〈Φ[h]〉 will increase asymptotically, as shown in Figure 6. In addition, it is well-known that v_∞ ∝ λ [3]. Hence (particularly in the F = 0 case) we find the λ² dependence for the slope of 〈Φ[h]〉 vs. t, as obtained in simulations.

We close this section by looking at the asymptotic behavior from another point of view. Here we restrict to d = 1, the continuum representation, and the F = 0 case. As already indicated, the statistics of the stochastic variable χ in (3), is of paramount importance [16, 17]. In addition to the forms of Φ[h] given by (5) and (14), we can also obtain an expression for Φ˙[h] as

We have the following relations

with χ′ :=∂xχ. The time derivative of (3) is

as the 2nd and 3rd terms decay very fast. On the other hand,

Exploiting these relations we obtain

If we look at long times, where the non-linear contribution dominates, we have

Since v_∞ ∝ λ, we get right away the expected result: the constant goes as λ².

We now show why 〈Φ[h]〉 and 〈h〉 are correlated. At long times, when 〈∫dx(∇δh)2〉→C0 and 〈h〉 → v_∞t, (13b)—now in terms of δh—“factorizes” out and thus δ3Φ[h(x,t)]→-λC0〈h〉. Hence the second term in (5)—which already dominates over the first one—behaves in the same way as the term in EW corresponding to an external force:

For any substrate dimensionality d, the KPZ equation can be expressed as a gradient flow by means of the functional Φ[h(x, t)] in (5). This is a state functional, i.e., it depends explicitly only on h(x, t). However, evaluating Φ[h] for a given field configuration h(x, t) requires integrating over its whole story, with initial condition h₀(x). So Φ[h] carries information on the whole growth story. If (5) is discretized according to (∇h)2→12[(hi+1-hi)2+(hi-1-hi)2] [11, 12], a cubic polynomial in the {h_i}—containing hi3 terms—is obtained.

Expansion (6) provides an “instantaneous landscape” of Φ[h(x, t)]. Even using (∇h)2→(hi+1-hi-1)2, its discretization also yields a cubic polynomial (in the {δh_i}) but now the absence of δhi3 terms enables a clearer understanding of “horizontal” and “vertical” symmetry breakdowns when λ ≠ 0—which may in turn be related to the asymmetry of the Tracy–Widom distribution—and brings an explanation to the shift in the dynamical exponent z (a change from diffusive to superdiffusive behavior).

Even though Φ[h] fulfills the Lyapunov property, its functional form is for short times that of a cubic polynomial, and is thus unbounded below. The KPZ equation turns out to be a high-dimensional (non Kramers) escape-like problem, and its asymptotic pdf cannot be normalized on the whole h–space (except for fixed and normalized forms of h₀). We have moreover provided a bridge between the observed asymptotic Tracy–Widom statistics, and the behavior of Φ (which as argued, provides an instantaneous landscape for the stochastic dynamics).

Given an ensemble of stories h(x, t), one can attempt to calculate the time dependence of the lower moments of Φ[h]. Our focus in this work is on 〈Φ[h]〉_t, the average of Φ[h] over noise realizations at time t. We have shown both numerically (Figures 3, 4) and analytically in (17) that in the absence of external forcing (F = 0), 〈Φ[h]〉_t → A − Bt as t → ∞ (in fact, surprisingly soon!), with B ~ λ². We have moreover traced back the origin of this asymptotic linear dependence on t to two known effects: (a) 〈h〉 → v_∞t, (b) 〈∫dx(∇h)2〉→C0, which lead to factorization of the “KPZ term” in 〈Φ[h]〉_t.

The time t_max at which 〈Φ[h]〉_t goes through a maximum, can be regarded as an estimation of the crossover time t_xovr from the EW regime to the KPZ one. We have found that tmax~λ-1, in agreement with previous results obtained by exploiting the time behavior of the stochastic action [46, 47].

Other issues worth remarking are that—because of the abovementioned linear dependence of 〈Φ[h]〉_t on t—the behaviors of 〈Φ[h]〉_t and 〈h〉 are correlated, and that (both with and without external forcing) the NEP's time behavior indicates a (non Kramers) activation-like phenomenon. In a forthcoming work [48] we undertake a thorough analysis of the time behavior of the first few higher moments of Φ[h], namely the NEP's dispersion σ(t), skewness γ₁(t) and kurtosis K(t). An interesting point (left to further work) would be to “interpret” the aging phenomenon [49, 50] from the NEP perspective.

We have discussed some strengths of a novel tool for describing the KPZ dynamics, the NEP approach. This innovative framework is expected to contribute answering questions that remain open by today, and handling the KPZ problem from still another perspective [13, 36].

All authors listed, have made substantial, direct and intellectual contribution to the work, and approved it for publication.