9512.net

甜梦文库

甜梦文库

当前位置：首页 >> >> # A Scaling Theory of Bifurcations in the Symmetric Weak-Noise Escape Problem

arXiv:cond-mat/9506097v1 22 Jun 1995

A Scaling Theory of Bifurcations in the Symmetric Weak-Noise Escape Problem

Robert S. Maier?

rsm @ math.arizona.edu Dept. of Mathematics University of Arizona Tucson, AZ 85721, USA

Daniel L. Stein?

dls @ physics.arizona.edu Dept. of Physics University of Arizona Tucson, AZ 85721, USA

Abstract We consider the overdamped limit of two-dimensional double well systems perturbed by weak noise. In the weak noise limit the most probable ?uctuational path leading from either point attractor to the separatrix (the most probable escape path, or MPEP) must terminate on the saddle between the two wells. However, as the parameters of a symmetric double well system are varied, a unique MPEP may bifurcate into two equally likely MPEP’s. At the bifurcation point in parameter space, the activation kinetics of the system become non-Arrhenius. We quantify the non-Arrhenius behavior of a system at the bifurcation point, by using the Maslov-WKB method to construct an approximation to the quasistationary probability distribution of the system that is valid in a boundary layer near the separatrix. The approximation is a formal asymptotic solution of the Smoluchowski equation. Our analysis relies on the construction of a new scaling theory, which yields ‘critical exponents’ describing weak-noise behavior at the bifurcation point, near the saddle.

Keywords : boundary layer, caustics, double well, Fokker-Planck equation, Lagrangian manifold, large deviation theory, large ?uctuations, Maslov-WKB method, non-Arrhenius behavior, nongeneric caustics, Pearcey function, singular perturbation theory, Smoluchowski equation, stochastic escape problem, stochastic exit problem, stochastically perturbed dynamical systems, weak noise, Wentzell-Freidlin theory, WKB approximation. Archive Paper Number : cond-mat/9506097

? ?

Partially supported by the National Science Foundation under grants NCR-90-16211 and DMS-95-00792. Partially supported by the U.S. Department of Energy under contract DE-FG03-93ER25155.

1

1

Introduction

In this paper we build on previous work [28, 29, 30] by analysing an unusual bifurcation phenomenon in the theory of noise-activated transitions. We study its appearance in the overdamped limit of two-dimensional double well systems, with nongradient dynamics. In this context, the new phenomenon is a bifurcation of the most probable transition path (in the limit of weak noise) between the two wells, as a system parameter is varied. In many ways, the behavior of a system whose most probable transition path is just beginning to bifurcate resembles that of a system undergoing a phase transition . In particular, double well systems which are ‘at criticality’ in the bifurcation sense will exhibit non-Arrhenius behavior. This means that the growth of the mean time between inter-well ?uctuations, i.e., the growth of the mean time needed for the system to hop from one well to the other, will not be pure exponential in the weak-noise limit. In double well systems at criticality, relaxation due to activation will proceed (in the limit of weak noise) at an anomalous, in fact anomalously large, rate. To treat the previously unnoticed phenomenon of bifurcation, we need to develop a new approach for treating transitions induced by weak noise, when a ‘soft mode’ appears in the dynamics of transverse ?uctuations around the most probable transition path. Since this is analogous to a phase transition, we introduce a scaling theory . In the context of double well systems, our scaling theory is a theory of behavior near the saddle point between the two wells, since the saddle is where the most probable inter-well transition path begins to bifurcate. We shall demonstrate that the theory explains the weak-noise behavior, at criticality, of a large universality class of double well systems. The scaling theory will reveal a striking feature of the bifurcation phenomenon, which is that in any ‘critical’ double well system there appears (in the weak-noise limit) a nongeneric singularity in the stationary probability distribution, located at the saddle point. As Berry [4] discusses, a singularity is nongeneric if it arises, in an appropriate WKB sense, from a catastrophe of unusual type; i.e., one of in?nite codimension . The stationary distribution near the saddle point is described, in the limit of weak noise, by an unusual (non-canonical) di?raction function. The familiar special functions of WKB theory (Airy functions, Pearcey functions, etc.) do not su?ce. The singularity at the saddle, and the di?raction function with which it is ‘clothed,’ can be viewed as the mathematical source of the non-Arrhenius weak-noise asymptotics. We begin with three largely qualitative sections. In Section 2 we review the physical relevance of overdamped models with non-gradient dynamics, and in Section 3 explain how the weak-noise behavior of any double well model of this type is determined by its ?ow ?eld of instanton trajectories (most probable ?uctuational paths). In Section 4 we sketch the gross features of the bifurcation phenomenon, including features such as further bifurcations and universality. In Section 5, our treatment becomes more quantitative. We ?rst review the matched asymptotic approximations technique we have employed elsewhere [28, 29, 30], and begin extending it to handle models with singularities. In Section 5.3 we explain why the bifurcation transition deserves to be called a phase transition. In particular, we explain 2

how behavior near criticality is described by critical exponents , which characterize the rate of divergence of measurable quantities (e.g., the pre-exponential factor in the weak-noise asymptotics of the mean inter-well ?uctuation time). In Section 5.4 we explain how to determine whether any given double well model is ‘critical.’ The transverse Jacobi operator is the di?erential operator appearing in the second variation of the Onsager-Machlup action functional, when one varies about the most probable inter-well transition path. The onset of bifurcation occurs when this operator acquires a zero eigenvalue. In Section 6 we explain the use we shall make of Maslov’s geometric theory of wave asymptotics [4, 25, 32]. In Section 7 we introduce the concept of a scaling theory, by developing a scaling theory of weak-noise behavior near generic (cusp) singularities. We show how the scaling theory justi?es the Ginzburg-Landau approximation used in this context by Dykman et al. [12]. In Section 8 we develop an analogous scaling theory for the nongeneric singularity associated with the onset of bifurcation. We compare our theory with numerical data, and examine its predictions for non-Arrhenius behavior and the stationary distribution near the saddle point. In Section 9 we discuss our results. The reader may wish to glance ahead at Fig. 14, which is an Arrhenius plot of the inter-well hopping rate of any double well system at criticality. The non-Arrhenius behavior shown there, in particular the ‘logarithmic bend,’ is the key result of this paper.

2

Preliminaries

Statistical physics and chemical physics include many examples of stochastically perturbed dynamical systems. It is often the case that the state of such a system is modelled as a particle moving in an n-dimensional force ?eld F (x), and subject to additional random perturbations (‘noise’). Since our interest is in the modelling of nonequilibrium systems, we shall not assume (as is usually done) that this force ?eld is conservative. If the motion of the particle is isotropically damped, with damping constant γ , in the absence of noise the particle position x would obey the deterministic equation ¨ + γmx ˙ = F (x). mx Adding a random force F random (t) yields the Langevin equation ¨ + γmx ˙ = F (x) + F random (t). mx (2) (1)

In √ physical problems F random (t) is often modelled as Gaussian white noise with amplitude 2γmkB T , where T is the ambient temperature. In this case the associated partial di?erential equation, which describes the time evolution of the probability density of x and its velocity, is known as the (forward) Fokker-Planck equation. A case particularly important in applications is the overdamped , or inertialess case, when 1 ¨ term in (2) can be dropped, γ ? t? 0 , for t0 the physical time scale. In this case the mx 3

and the Langevin equation becomes ?rst order in time. If time is rescaled by a factor γm (i.e., t ← γmt), it may be written in the normalized form ˙ = u(x) + ?1/2 w ˙ (t). x (3)

˙ (t) is a standard n-dimensional Gaussian white noise (the derivative of w(t), a stanHere w dard n-dimensional Wiener process), the ‘drift ?eld’ u equals F , and ? equals 2kB T . The corresponding scalar advection-di?usion equation for the probability density ρ = ρ(x, t) of x, ρ ˙ = (?/2)?2 ρ ? ? · (ρu), (4)

is known as the (forward) Smoluchowski equation. It may be written as ρ ˙ = ?L? ρ, where L? ≡ ?(?/2)?2 + u · ? + ? · u. (5)

It is often necessary to generalize the equation to include the e?ects of anisotropic damping [21], or state-dependent noise [23]. However, in this paper we consider only overdamped systems whose Langevin equation is of the form (3). Since we do not require the deterministic forces to be conservative, we do not require u to be a gradient ?eld . This means that even in stationarity, the system may not display detailed balance. Equivalently, the stationary probability distribution for the system may not be (in the traditional sense) in thermal equilibrium. Attractors of the drift ?eld u, in particular point attractors, correspond to ‘metastable states’: they are stable states of the underlying deterministic dynamics, but the thermal noise may induce transitions between them. Of great physical interest is the time needed for this to occur. For example, how long does it take for the noise in (3) to overcome the drift toward a speci?ed stable point S , and drive the system state x beyond the domain of attraction of S , toward another attractor? The study of such noise-activated transitions is known as the stochastic exit problem , or the escape problem . For general stochastic models only numerical results can be obtained (see, e.g., Ref. [7]). The Smoluchowski equation is particularly di?cult to handle in the ? → 0 limit. This is the weak-noise, or low-temperature limit, in which the mean ?rst passage time (MFPT) τ from S to the boundary of its domain of attraction grows exponentially. In this limit a single escape path (the most probable escape path, or MPEP) usually dominates. Our approach to the weak-noise limit, which does not rely on a numerical simulation of the Smoluchowski equation, will exploit this asymptotic determinism quite heavily.

3

Symmetric Double Well Models

As in two of our earlier papers on the stochastic exit problem [28, 29], we shall focus on two-dimensional ‘double well’ systems, with smooth drift ?eld u = (ux , uy ) of the symmetric form shown in Fig. 1. If x = (x, y ) is the two-dimensional state variable, ux (x, y ) is taken to be odd in x and even in y , while for uy (x, y ) the reverse is true. There is assumed to be 4

y

x

– xs

0

xs

Figure 1: The streamlines of a typical symmetric double well drift ?eld, indicating the path taken by the particle in the absence of noise. There are point attractors at S = (xs , 0) and S ′ = (?xs , 0), and a saddle point at (0, 0). a linearly stable point attractor S = (xs , 0) whose domain of attraction is the entire open right-half plane. By symmetry, its re?ection S ′ = (?xs , 0) attracts the open left-half plane. There is also assumed to be a single saddle , or hyperbolic point, on the y -axis separatrix between the two domains of attraction. It must be at the origin, by symmetry. Nongradient drift ?elds with this topology arise in statistical and chemical physics, and also in theoretical biology, e.g., in stochastic competition models of population dynamics [31]. One expects that as ? → 0, exit from either of the two domains of attraction will occur preferentially over the saddle. The drift ?eld u is assumed to have a nondegenerate linearization at the saddle. So λx = ?ux /?x(0, 0) > 0, and λy = ?uy /?y (0, 0) < 0. We shall see that the character of the abovementioned bifurcation phenomenon depends strongly on the quotient ? ≡ |λy |/λx . A typical (and not necessarily gradient) symmetric double well drift ?eld, which we have used elsewhere for purposes of illustration and shall examine further below, is ux (x, y ) = x ? x3 ? αxy 2 , uy (x, y ) = ??(1 + x2 )y,

(6)

in which ? appears as a parameter. We shall call this drift ?eld the ‘standard’ double well model. For any choice of ? > 0, its structure is that of Fig. 1, with S = (xs , 0) = (1, 0). 5

It is not a gradient ?eld unless the parameter α equals ?. If α > 0 it has a very signi?cant additional property, which we shall require of all our double well models. This is the property that ? 2 ux /?y 2 (x, 0), which by symmetry is an odd function of x, is strictly negative for all x between 0 and xs . If this is the case, the drift from the saddle toward S ‘softens’ as one moves away from the x-axis. The o?-axis softening, for the standard model, increases as α is increased. We remind the reader of our approach to the weak-noise limit of stochastically perturbed dynamical systems. (We review the mathematical aspects in Section 5.1 and 5.2.) Suppose that such a system has a unique stationary probability density ρ0 , which satis?es the timeindependent Smoluchowsk equation L? ρ0 = 0. Typically, as the noise strength ? → 0, ρ0 takes on an asymptotic WKB form. In fact for certain functions W and K whose smoothness properties we shall leave unspeci?ed; K , in particular, may have singularities. In any double well model, by convention W = 0 and K = 1 at x = S and S ′ . Moreover, W > 0 at all points x other than S and S ′ . W is called the nonequilibrium potential of the model [16]. If the drift u equalled the negative gradient of a potential Φ, then W would equal 2Φ, K would reduce to a constant, and the WKB form (7) would reduce to a Maxwell-Boltzmann distribution. For systems with nongradient dynamics, the computation of W and K is more complicated. In general W has an alternative interpretation as a classical action function . As we review in Section 5.1, this is because the WKB approximation (7) is determined by a ?ow ?eld of ‘classical’ trajectories, or WKB characteristics, emanating from the attractors of the deterministic dynamics (e.g., S and S ′ ). These classical trajectories (sometimes called instanton trajectories [6, 28], or optimal trajectories [12]) have a physical interpretation as most probable ?uctuational paths . In the double well case, the exponentially rare ?uctuations from S (resp. S ′ ) to any point x in its domain of attraction become increasingly concentrated around the classical trajectory extending from S (resp. S ′ ) to x. Equivalently, the most probable ‘prehistory’ of any ?uctuation passing through x extends back toward S or S ′ along this trajectory [11]. The trajectories are determined by a classical Lagrangian (the Onsager-Machlup Lagrangian), and W (x) is obtained by integrating this Lagrangian along the classical trajectory terminating at x. W (x) is interpreted as the rate at which ?uctuations to the neighborhood of x are suppressed exponentially, as ? → 0. In symmetric double well models the stationary density ρ0 (and hence W ) must be even in x. In the ? → 0 limit the phenomenon of noise-activated hopping between the two wells is governed by the closely related quasistationary density ρ1 , which is odd rather than even. The quasistationary density is the next lowest lying (i.e., slowest decaying) eigenmode of the Smoluchowski operator L? . For any choice of initial conditions the probability density ρ = ρ(x, t) necessarily satis?es for some constant C , where λ1 is the eigenvalue of ρ1 . The exponential decay of the quasistationary eigenmode is interpreted as describing the equilibration of probability between 6 ρ(x, t) ? ρ0 (x) + Cρ1 (x) exp(?λ1 t), t→∞ (8) ρ0 (x) ? K (x) exp[?W (x)/?], ? → 0, (7)

the two wells due to noise-activated hopping, or the absorption of probability on the separatrix [7]. ρ1 of course satis?es Dirichlet (absorbing) boundary conditions on the separatrix. Its eigenvalue λ1 = λ1 (?) normally falls to zero exponentially as ? → 0. The exponentially small splitting between the ground state eigenvalue λ0 ≡ 0 and the eigenvalue λ1 is analogous to the exponentially small splitting (as h ? → 0) between the ground state and ?rst excited state of a quantum-mechanical Hamiltonian with double well potential. Both are WKB phenomena. In the ? → 0 limit, λ1 is interpreted as the rate at which noise-activated hopping takes place. Equivalently, it is a reciprocal MFPT . The techniques reviewed in Section 5.1 and 5.2 permit a computation of the ? → 0 asymptotics of the eigenvalue λ1 , and hence of the MFPT τ , in most symmetric double well models. Our basic approach is similar to that of Kramers [22]. In the limit of weak noise we approximate ρ1 (x, y ) by ρ0 (x, y ) sgn(x) except in a ‘boundary layer’ of width O (?1/2 ) near the x = 0 separatrix, and compute λ1 as the rate at which probability is absorbed on the separatrix. Performing this computation requires the construction of a boundary layer approximation to ρ1 , valid near the saddle, and matching to the ‘outer’ approximations on either side [8]. Normally, we ?nd τ ? A exp[+W (0, 0)/?], where A ∝ K (0, 0)?1 . So the asymptotic MFPT growth rate in the limit of weak noise is simply ?W = W (0, 0) ? W (S ) = W (0, 0), the height of the ‘action barrier,’ or activation barrier, between the two wells. And the MFPT generally displays a pure exponential (Arrhenius) growth, with an explicitly computable (?-independent) prefactor. We shall see, however, that the bifurcation phenomenon may induce more complicated (non-Arrhenius) weak-noise asymptotics for the MFPT.

4

The Bifurcation Phenomenon: Qualitative Features

We pointed out in Ref. [28] that a bifurcation phenomenon may occur in double well models as their parameters are varied. Figure 2 displays the ?ow of instanton trajectories (i.e., most probable weak-noise ?uctuational paths) emanating from the stable point S = (1, 0) in the standard model (6) with ? = 1, at several values of the parameter α. When 0 < α < 4 the general picture resembles Fig. 2(a): the line segment from (1, 0) to the saddle (0, 0) is the only instanton trajectory from S to the saddle. This line segment is interpreted as the most probable escape path (MPEP). In the weak-noise limit, the (exponentially rare) ?uctuations from the right-half plane to the left-half plane proceed preferentially along it. To leading order, activation kinetics reduce to instanton dynamics. As α is increased, there is a qualitative change, akin to a phase transition, in the behavior of the instanton trajectories. This takes place at the critical value α = 4, as shown in Fig. 2(b). When α > 4 as in Fig. 2(c), they focus at a point (xf , 0) on the x-axis, with xf > 0. xf converges to zero as α → 4+ , so one may speak of the focal point ‘being born’ at criticality, and ‘emerging from the saddle’ as α is increased above its critical value. In geometrical optics the focal point would be called a cusp . From it there extends a fold , or caustic (an envelope 7

0.6 0.4 0.2 0 -0.2 -0.4 -0.6 0

0.6 0.4 0.2 0 -0.2 -0.4 -0.6 0

0.2

0.4

0.6

0.8

1

1.2

1.4

1.6

0.2

0.4

0.6

0.8

1

1.2

1.4

1.6

0.6 0.4 0.2 0 -0.2 -0.4 -0.6 0

0.6 0.4 0.2 0 -0.2 -0.4 -0.6 0

0.2

0.4

0.6

0.8

1

1.2

1.4

1.6

0.2

0.4

0.6

0.8

1

1.2

1.4

1.6

Figure 2: The ?ow ?eld of outgoing instanton trajectories (i.e., most probable ?uctuational paths, in the weak-noise limit) emanating from the stable point S = (1, 0) of the standard double well model (6). Here ? = 1, and parts (a), (b), (c), (d) of the ?gure illustrate the cases α = 1, 4, 5, 10. The α = 4 model is ‘critical’ in the bifurcation sense. Increasing α above 4 causes the instanton trajectories to focus, and the MPEP to bifurcate.

8

of crossing trajectories, with ?y ∝ (?x)3/2 ). Each point in the sharp-tipped region within the fold is reachable from S via three instanton trajectories. Of these three trajectories, only the one(s) with minimum action are ‘physical,’ and can be interpreted as most probable ?uctuational paths. For example, on-axis points (x, 0) with 0 ≤ x < xf are reachable via an on-axis (straight) trajectory, and via two additional symmetrically placed o?-axis (curved) trajectories. Computation shows that the o?-axis trajectories have lesser action, and are dominant. The true (‘least action’) MPEP’s in Fig. 2(c) are accordingly the symmetrically placed pair of o?-axis trajectories, one above and one below the x-axis, that terminate on the saddle. Note that beyond the cusp (i.e., at x < xf ), the physical action W is no longer di?erentiable through the x-axis. This nondi?erentiability arises from di?erent dominant o?-axis trajectories being selected as y → 0+ and y → 0? . The transition at α = 4 can be interpreted as a bifurcation of the MPEP , corresponding to a sort of symmetry breaking. At larger values of α, the drift ?eld u and the Langevin equation (3) remain symmetric about the x-axis, but each of the two MPEP’s is not. The line segment from S to the saddle, formerly the (unique) MPEP, in no way contributes to the leading weak-noise asymptotics for escape. (It remains an extremum of the Onsager-Machlup action functional, but is no longer the minimum.) The occurrence of a bifurcation in the standard model at su?ciently high α (when ? = 1, at α = 4) is due to the fact that by increasing α, one softens the resistance to motion toward the separatrix in the vicinity of the x-axis (though not on the x-axis itself). This enhances the probability of escape trajectories that deviate from the axis. Of course it is only in the limit, as ? → 0, that well-de?ned MPEP’s appear. And the existence of a sharp, well-de?ned transition when α equals some critical value αc is not at all obvious! When α is increased beyond αc , further bifurcations of the on-axis instanton trajectory will occur. In Section 5.4 we explain how the critical values of α are determined by a Jacobi equation , with a classical mechanical interpretation. It turns out that in the standard model (j ) with ? = 1 the j ’th bifurcation occurs at α = αc = (j +1)2 . Figure 2(d) shows the situation (2) at α = 10, when a second focus (xf , 0) has emerged, with its own caustic. Beyond the ?rst focus (xf , 0) each point on the x-axis is reached from S by three instanton trajectories; beyond the second focus, each such point is reached by ?ve. The MPEP’s in Fig. 2(d), however, remain the symmetrically placed pair of o?-axis trajectories that terminate on the saddle. Computation shows that the oscillatory trajectories from S to the saddle arising from the second, third,... bifurcations have higher actions, and are accordingly not physical. That caustics can occur in the ?ow pattern of the most probable ?uctuational paths has been known for some time [9, 19], but our Ref. [28] was the ?rst to consider the e?ects on exit phenomena. We shall see that what occurs at the ?rst critical value of α has much in common with a critical point characterizing a phase transition in a condensed matter system. This is suggested by Fig. 3, which plots the activation barrier ?W = W (0, 0), as determined by the true MPEP or MPEP’s, as a function of α for the standard model with ? = 1. (Recall that ?W = W (0, 0) is the exponential growth rate of the MFPT as the noise strength tends to zero.) W (0, 0) decreases above α = αc = 4 as the bifurcating MPEP’s move away from 9

?W 0.55 0.5 0.45 0.4 0.35 0 2 4 6 8 10 α

Figure 3: The activation barrier ?W = W (0, 0) ? W (S ) = W (0, 0) between the two wells of the standard double well model (6), as a function of the o?-axis softening parameter α. Here ? = 1. The lowering of the activation barrier beyond α = αc = 4 is due to the bifurcation of the MPEP, along which the action di?erence ?W is computed.

10

the x-axis. The WKB prefactor K (0, 0) turns out to be singular at the bifurcation transition; ? in Section 5.3 we note that it diverges as α → αc . As a consequence, in the standard model at least, the weak-noise MFPT asymptotics at criticality cannot be of a pure Arrhenius form . There is in fact a set of critical exponents describing the behavior of the ? → 0 asymptotics (1) of the standard model as α tends to its (?-dependent) ?rst critical value αc = αc , and as x → 0, y → 0. It is a reasonable conjecture that behavior near criticality is universal in the sense that it does not depend on the details of the stochastic model exhibiting the bifurcation phenomenon. To analyse the critical behavior and demonstrate universality, in Sections 7 and 8 we begin the construction of a scaling theory of the bifurcation phenomenon. Our treatment extends from the standard model (6) to any symmetric double well model with a similar ‘o?-axis softening parameter’ α, and a ?rst critical value αc . We ?rst identify the singular behavior, for any double well model at criticality, of the action W and the WKB prefactor K at the saddle point. We then show that at criticality, the stationary density ρ0 and the quasistationary density ρ1 may be approximated on an appropriate (?dependent) length scale near the saddle point by certain ‘di?raction functions,’ which have explicit integral representations. The technique for constructing these representations is due to Maslov [32], and ultimately to Keller [20]. It was Maslov who ?rst worked out, in the context of wave ?elds, the di?raction functions that ‘clothe’ generic singularities other than cusps and folds. A very important discovery, from a mathematical point of view, will be that when α equals the critical value αc where the MPEP begins to bifurcate, the saddle point (0, 0) acquires a certain nonzero singularity index . What this means is best understood by comparing the singularity at the saddle (when α = αc ) with the cusp and fold singularities present when α > αc . The terminology of geometrical optics [4] is appropriate. The cusp at (xf , 0) is a structurally stable singularity (or catastrophe , in the language of Thom), with codimension 2. The fold extending from it, though not ‘physical’ in the above least-action sense, is a catastrophe of codimension 1. For points x in the vicinity of the cusp, the WKB approximation (7) for the value of the stationary density ρ0 (x) breaks down. The proper treatment of points near the cusp and the fold is similar to the short-wavelength treatment of wave ?elds near caustics [4, 12]. The cusp is said to have singularity index 1/4, and points on the fold would (if it were physical) have singularity index 1/6. This means that at these singular points the prefactor in the WKB approximation to ρ0 , which formally diverges, if properly constructed would acquire a factor ??1/4 (resp. ??1/6 ). There is a non-WKB (but uniformly valid) approximation to ρ0 (x) in the vicinity of each such singular point, in terms of canonical di?raction functions. We shall re-derive these facts in Section 7, in terms of scaling functions. We shall show in Section 8 that the singularity index of the point singularity appearing at the saddle, in critical models, depends in a universal way on ?, i.e., on the ratio of the eigenvalues of the linearization of the drift u at the saddle. It turns out to equal (? + 1)/6. Moreover, the approximations to ρ0 and ρ1 near the saddle are given by non-canonical di?raction functions. By using the non-canonical approximation to ρ1 to compute the rate at which probability is absorbed on the separatrix we shall quantify the universal non11

1

0.5

0

-0.5

-1 -1.5

-1

-0.5

0

0.5

1

1.5

Figure 4: The ?ow ?eld of the instanton trajectories emanating from both stable ?xed points S and S ′ , in a critical version of the standard double well model (6). This ?gure reveals that at criticality, a two-sided caustic extends sideways from the saddle point. Although a universal phenomenon, this caustic is nongeneric in the sense of singularity theory. Here ? = 1, and the parameter α is set equal to the corresponding ?rst critical value αc = 4, as in Fig. 2(b). Arrhenius behavior of the weak-noise MFPT asymptotics. We shall show that in symmetric double well models at criticality, τ ? const × ?s exp[+?W/?], ? → 0, (9)

where s = s(?) = (? + 1)/6 is the index of the singularity at the saddle. We shall also derive scaling corrections, in the weak-noise limit, to the normal distribution of exit location points near the saddle. The preceding results will hold for all ? satisfying 3/4 < ? < 3; the weak-noise asymptotics of models with ? ≤ 3/4 and ? ≥ 3 are still under investigation. The point singularity appearing in critical models at the saddle point, which may be termed a nascent cusp , is nongeneric . It is not a member of the well known family of singularities that includes folds, cusps, swallowtails, etc. This becomes clear if one plots the ?ow ?eld of the instanton trajectories emanating from both S and S ′ in the standard model at criticality (? varying, and α set equal to its ?-dependent ?rst critical value). At least 12

when 3/4 < ? < 3/2, one ?nds that at criticality a two-sided caustic extends transversally from the saddle point itself . (Cf. Dykman et al. [12].) Figure 4, which is an extended version of Fig. 2(b), shows the ?ow ?eld when ? = 1 and α = αc = 4. The caustic is clearly visible. It is not ‘physical,’ since it is formed by high-action instanton trajectories that have crossed the separatrix. But the ‘nascent’ cusp is clearly a cusp in its own right, of an unusual sort. Numerically one ?nds that the two-sided caustic extending from it is located at

?1 < |x| ? const × |y |(3/2??) ,

y → 0.

(10)

A conventional (generic) caustic would have an exponent of 3/2. The continuously varying exponent (3/2 ? ?)?1 , which turns out to be universal and which we shall derive in Section 8 from our scaling theory, signals that the two-sided caustic is nongeneric. The nascent cusp from which it extends is itself nongeneric in the sense of singularity theory. As Berry [4] has emphasized, nongeneric singularities arise from catastrophes of in?nite codimension . It is remarkable that a singularity of such complexity is a universal feature of singly parametrized symmetric double well models with non-gradient dynamics.

5

Quantitative Semiclassical Asymptotics

We now begin a quantitative treatment of the weak-noise asymptotics for escape. We ?rst recast our earlier results in a form that facilitates the analysis of singularities. In Section 5.1 we discuss geometric aspects of the WKB approximation, and in Section 5.2 we discuss our matched asymptotic approximations technique for computing MFPT asymptotics. In Section 5.3 we use the standard model (6) to illustrate the nature of the nascent cusp appearing at the saddle at criticality, and the ways in which bifurcation can be viewed as a phase transition. In Section 5.4 we explore the bifurcation phenomenon from a classical mechanical point of view, and relate it to the appearance of a transverse soft mode . We explain how its appearance is governed by a Jacobi equation, and how this equation determines whether or not a given double well model is at criticality.

5.1

The WKB Approximation and Classical Mechanics

The time-independent forward Smoluchowski equation L? ρ0 = 0 may be written as H (x, ???)ρ0 = 0 where H (x, p) = p2 /2 + u(x) · p (12) is the so-called Wentzell-Freidlin Hamiltonian [14], whose dual is the Onsager-Machlup Lagrangian [34] ˙ ) = |x ˙ ? u (x )| 2 / 2 . L(x, x (13) 13 (11)

In equation (11) we have adopted an operator ordering convention according to which the action of ? precedes that of x. In the weak-noise (? → 0) limit the stationary density ρ0 and the quasistationary density ρ1 are given in the interior of each well by a WKB, or semiclassical form. A full WKB expansion for ρ0 would be of the form ρ0 (x) ? [K (0) (x) + ?K (1) (x) + · · ·] exp[?W (x)/?], ? → 0, (14)

as in geometrical optics. By substituting this formal series into (11) and examining the coe?cients of each power of ?, one obtains equations for W and the K (m) . That is what we shall do, though we shall work only to leading order: our WKB Ansatz will be ρ0 (x) ? K (x) exp[?W (x)/?]. Notice that since the eigenvalue λ1 = λ1 (?) of ρ1 is exponentially small as ? → 0, the asymptotic expansions (in powers of ?) for ρ1 and ρ0 will be the same. To see the di?erence between them, which is signi?cant only near the separatrix between the two wells, one would have to go ‘beyond all orders’ in the WKB expansion. The eikonal equation for W is the time-independent Hamilton-Jacobi equation H (x, ?W ) = 0 (15)

so that W is a classical action at zero energy . For any point x in either well, it may be computed by integrating the Lagrangian along the zero-energy classical trajectory extending from S (resp. S ′ ) to x. Each such trajectory, which satis?es the Euler-Lagrange equations, is interpreted as a most probable ?uctuational path in the ? → 0 limit. These trajectories are the ‘instanton trajectories’ of the last section; the term is justi?ed by analogy with the semiclassical limit in quantum mechanics and quantum ?eld theory [3]. In the language of Gutzwiller [18], the points at which the instanton trajectories focus would be called zeroenergy conjugate points. It is convenient to work in the Hamiltonian picture, according to which the classical trajectories of interest lie on a zero-energy surface in a nonphysical phase space, coordinatized by position x and momentum p. The ?ow in this phase space (2n-dimensional, if con?guration space is n-dimensional) is determined by Hamilton’s equations and the WentzellFreidlin Hamiltonian. From this point of view the instanton trajectories of Figs. 2 and 4 are mere images of phase-space trajectories, projected ‘down’ to con?guration space by the map (x, p) → x. The phase-space trajectories emanate from (S, ?) (resp. (S ′ , ?)). In WKB theory the projected trajectories are traditionally called characteristics , and the phase space trajectories bicharacteristics . Characteristics may intersect, as in Figs. 2(c) and 2(d), but bicharacteristics may not. It is easy to verify, using Hamilton’s equations, that (S, ?) and (S ′ , ?) are hyperbolic ?xed points of the Hamiltonian ?ow. And the unstable manifold of (S, ?), for example, comprises all points (x, p) that lie on one of the bicharacteristics emanating from (S, ?). The unstable manifolds of (S, ?) and (S ′ , ?) are Lagrangian [25]: they are invariant under the Hamiltonian ?ow. By the term ‘Lagrangian manifold’ we shall refer to either of these two unstable manifolds, or their union. We denote by M this union, i.e., the set of all points (x, p) that 14

lie on a bicharacteristic emanating from either (S, ?) or (S ′ , ?). If con?guration space is n-dimensional, M will be an n-dimensional manifold. Each point P = (x, p) on M has a value for the zero-energy action W associated with it, computed by W (P ) = p · dx, (16) the line integral being taken along the bicharacteristic terminating at P . If due to intersecting characteristics, or the crossing of characteristics from one well to the other, there are several manifold points Pi = (x, p(i) ) ‘above’ some point x, then W (x) and its gradient p = p(x) will in a mathematical sense be multivalued. As a function of x, W may in fact have branch points, branch lines (cuts), etc. But the physical action W (x) appearing in the WKB approximation will be single-valued: it will equal the minimum of the values W (Pi ) at the manifold points above x. This ‘least action’ computation determines which instanton trajectories are physical. The WKB prefactor K satis?es an easily derived transport equation. (Cf. Talkner [38].) ˙ = p + u(x) (which is one of Hamilton’s equations), the transport If one uses the fact that x equation takes on the comparatively simple form ˙ = ?[? · u + ?2 W/2]K, K (17)

the time derivative referring to instanton transit time, i.e., to motion along a characteristic or bicharacteristic. Similarly to W , K may be regarded as a function on M rather than on con?guration space. Integration of the equation (17) requires knowledge of the second spatial derivatives of W along the characteristic. But (? 2 W/?xi ?xj )(x) equals (?pi /?xj )(x), which is a measure of the ‘slope’ of the manifold above the point x. By di?erentiating Hamilton’s equations it is easy to show that the Hessian matrix Z = (Zij ) whose elements are the partial slopes ?pi /?xj satis?es the matrix Riccati equation ˙ = ?Z 2 ? ZB ? B t Z ? Z pl Y (l)

l

(18)

along any characteristic. (Cf. Ludwig [26].) Here B = (?ui /?xj ) and Y (l) = (?ul /?xi ?xj ) are auxiliary matrices. Since ?2 W = tr Z , the computation of K by numerical integration is straightforward. It is interesting to compare these results with those of Littlejohn [25] on the WKB prefactor for the solutions of the Schr¨ odinger equation in the semiclassical (? h → 0) limit. He introduces a Lagrangian manifold, and a similar integration along characteristics. But because he analyses the time-dependent Schr¨ odinger equation, he ?nds that the transport equation for his analogue of K can be integrated explicitly, yielding a Van Vleck determinant. Matters are not so simple in the time-independent case, for the Schr¨ odinger equation as well as for the Smoluchowski equation. Our WKB analysis of the weak-noise limit of the stationary density actually has more in common with the work of Gutzwiller on the semiclassical approximation of ?xed-energy quantum-mechanical Green’s functions [17, 18, 36] than it does with the semiclassical approximation of time-dependent quantum-mechanical propagators. 15

The prefactor K is analogous to the prefactor of a semiclassical Green’s function (at ?xed energy). It can in fact be related to the density of bicharacteristics on the Lagrangian manifold. This resembles Gutzwiller’s interpretation of the prefactor of a semiclassical Green’s function in terms of the density of classical trajectories on an energy surface [18].

5.2

Matched Asymptotic Approximations and the MFPT

We now specialize to two-dimensional double well models with the structure of Fig. 1. On account of symmetry and smoothness we may expand the drift u = (ux , uy ) thus: ux (x, y ) = v0 (x) + v2 (x)y 2 + · · · uy (x, y ) = u1 (x)y + u3 (x)y 3 + · · · . (19)

By assumption v0 (x) > 0 for all x between 0 and xs , and u1 (x) < 0 for all x between 0 and xs inclusive. If the symmetry through the axis is unbroken, W and K (both of them computed by integration along instanton trajectories emanating from S ) will have similar expansions W (x, y ) = w0 (x) + w2 (x)y 2 /2! + · · · K (x, y ) = k0 (x) + k2 (x)y 2 /2! + · · · . (20) (21)

Here w2m (x) ≡ ? 2m W/?y 2m (x, 0) and k2m (x) ≡ ? 2m K/?y 2m (x, 0). Since W can be viewed as a classical action, the functions w2m can be expressed in terms of the momentum p = p(x) of ′ (x) = px (x, 0) the instanton trajectories passing through near-axis points x. For example, w0 and w2 (x) = ?py /?y (x, 0). Substituting the WKB Ansatz into the Smoluchowski equation L? ρ0 = 0, and examining the coe?cients of each power of ? and y , will yield equations for the ′ = px = ?2v0 , various coe?cient functions in (20) and (21). One ?nds in particular that w0 or xs v0 (x′ ) dx′ . (22) w0 (x) = 2

x

Therefore the Hamilton equation x ˙ = px + v0 (x), which follows from the Wentzell-Freidlin Hamiltonian (12), implies that x ˙ must equal ?v0 (x) at all points between S and the saddle. The instanton trajectory on the x-axis moves with a speed equal to the local value of the drift speed, but in the direction opposite to the drift. Examining coe?cients also yields the two equations ˙ 0 = ?[u1 + w2 /2]k0 k 2 w ˙ 2 = ?w2 ? 2u1w2 + 4v0 v2 (23) (24)

˙ for ?v0 k ′ , and w where we have changed the independent variable from x to t by writing k ˙2 0 ′ for ?v0 w2 . Equations (23) and (24) could equally well be deduced from (17) and (18). For later reference we note that

′ w ˙ 4 = ?v0 w4 ′ 2 ′ = ?4(w2 + u1 )w4 ? 3[(w2 ) + 4v2 w2 + 8u3 w2 ] + 48v0 v4

(25)

16

is the equation satis?ed by the fourth derivative w4 = ? 4 W/?y 4 = ? 3 py /?y 3 on the x-axis. The physical interpretation of the functions k0 and w2 is straightforward. The WKB Ansatz implies that ρ0 (x, y ) ? k0 (x) exp ? w0 (x) + w2 (x) y2 /? , 2! ? → 0. (26)

So when ? is small, w2 governs the small transverse ?uctuations about the x-axis. At any time when the system state x has ?uctuated leftward from xs to x, the distribution of the transverse component y will (provided that w2 (x) > 0) be approximately Gaussian, with variance ? ?/w2 (x). Of course such ?uctuations are exponentially rare, on account of the exp[?w0 (x)/?] factor. The Riccati equation (24) therefore gives the position dependence of the width of the ‘WKB tube’ of probability density surrounding the MPEP, when this MPEP is in fact the line segment between S and the saddle. Moreover this equation captures the essence of the bifurcation phenomenon, as we shall see in Section 5.4. For the moment we note only that it may readily be integrated from t = ?∞ (when the instanton trajectory formally emerges from S ) to t = +∞ (when the trajectory, obeying x ˙ = ?v0 (x), reaches the saddle). Since v0 (x(t)) → 0 as t → ±∞, we see from (24) that w2 must converge as t → +∞ (resp. t → ?∞) to one of the two zeroes of the quadratic polynomial

2 ? w2 ? 2u1 w2 = ?w2 (w2 ? 2|u1 |),

(27)

near S , the same being true of ρ1 . Since νy = u1 (xs ), this will match to the tube approximation if w2 (xs ) = 2|u1(xs )|. Similarly, on the O (?1/2 ) lengthscale near the saddle, ρ0 may be approximated by the inverted Gaussian ρ0 (x, y ) ? const × e+λx x

2 /?

where u1 signi?es u1 (0) (resp. u1 (xs )). On physical grounds one expects that usually (‘generically’) the WKB tube will have a ?nite variance at both endpoints, i.e., as t → ?∞ and x → xs , and as t → +∞ and x → 0. So w2 (xs ) should equal 2|u1 (xs )|, and w2 (0) should equal 2|u1 (0)|. If these endpoint (‘turning point’) conditions hold, it is easy to match the tube approximation (26) to auxiliary, non-WKB approximations valid near the endpoints: the stable points and the saddle. On physical grounds, ρ0 and ρ1 may be approximated on the O (?1/2 ) lengthscale near S by a Gaussian function of the system state x. Let us write νx and νy for ?ux /?x(S ) and ?uy /?y (S ), the two (negative) eigenvalues of the linearization of the drift u at S . Then 2 2 ρ0 (x, y ) ? const × e?|νx |(x?xs ) /? e?|νy |y /? , ? → 0 (28)

e?|λy |y

2 /?

,

? → 0.

(29)

Since λy = u1 (0), the tube approximation will match to (29) if w2 (0) = 2|u1(0)| and k0 (0) is ?nite and nonzero. It is easy to verify that the approximations (28) and (29) satisfy the time-independent Smoluchowski equation on the O (?1/2 ) lengthscale near their respective turning points. 17

The appropriate (generic) approximation to the quasistationary density ρ1 near the saddle is slightly more complicated; it is an error function approximation of the sort ?rst used by Kramers [22]. We have [29] ρ1 (x, y ) ? const × ??1/2

∞ 0 ?|λy |y exp ?(xpx + p2 x /4λx )/? dpx e

2 /? 2 /?

(30) (31)

1/2 = const × erf(λx x/?1/2 )e+λx x

e?|λy |y

2 /?

on the O (?1/2 ) lengthscale. This ‘boundary layer’ approximation agrees with the inverted Gaussian approximation (29) in the far ?eld, i.e., as x/?1/2 → +∞. So under the same conditions, the tube approximation (26) will match to it. We have now approximated ρ1 (x) at all points x in the vicinity of the line segment joining S and the saddle. We must emphasize that the validity of this procedure depends on two assumptions: ? That the physical values of W (x) and K (x) at all points x along the axis arise from integration along the on-axis instanton trajectory extending from S to x. ? That the WKB tube surrounding the axis is well behaved as the saddle is approached, so that the error function approximation to the quasistationary density is valid near the saddle. This requires that w2 → 2|u1(0)|, and that k0 tend to a ?nite, nonzero limit. The ?rst assumption breaks down when the MPEP has bifurcated, and we shall see that the second assumption breaks down at the onset of bifurcation. But if both assumptions hold, it is easy to compute the weak-noise asymptotics of the quasistationary eigenvalue λ1 and its asymptotic reciprocal, the MFPT τ . The time-dependent equation ρ ˙ = ?L? ρ may be written as ρ ˙ + ? · [?(?/2)?ρ + ρu] = 0. (32) Equation (32) is a continuity equation, and j = ?(?/2)?ρ+ ρu can be viewed as a probability current density. Since λ1 is the decay rate of the eigenmode ρ1 , it may be computed as the rate at which probability is absorbed on the separatrix [22, 33]. Necessarily λ1 =

∞ ?∞

[?jx (0, y )] dy

∞ 0

∞ ?∞

ρ1 (x, y ) dy dx

(33)

where j = (jx , jy ) is computed from ρ1 . The numerator (an absorption rate) is computed from (31), and the normalization factor in the denominator from the Gaussian approximation (28). If the constant prefactors of these two approximations are chosen to ensure consistency with the intermediate WKB tube approximation (26), the quotient will acquire a factor k0 (0) exp[?w0 (0)/?], i.e., K (0, 0) exp[?W (0, 0)/?]. This computation, if carried through, yields a so-called Eyring formula for the weak-noise asymptotics of the quasi-stationary eigenvalue, i.e., the weak-noise asymptotics of the rate

18

of noise-activated hopping [8, 15]: λ1 (?) ? K (0, 0) ?

?

λx |νx | |νy |/|λy | π

?

? exp[??W/?],

? → 0.

(34)

Here the presence of the ‘frequency factor’ K (0, 0) is attributable to the non-gradient dynamics; it will equal unity if the drift u is a gradient. The formula is otherwise familiar. 1 Since τ ? λ? 1 as ? → 0, this formula predicts a pure Arrhenius growth of the MFPT in the weak-noise limit. But as noted, this conclusion depends crucially on the validity of the Kramers-type error function approximation to the quasistationary density near the saddle. This approximation will prove not to be valid in double well models undergoing a bifurcation.

5.3

Indications of a Phase Transition

We shall now explain how the bifurcation transition displays characteristic features of a phase transition, such as power-law divergences governed by critical exponents. We begin by using the standard double well model (6), and the transport equations of the last section, to reveal the nature of the ‘nascent cusp’ singularity appearing at the saddle point, at criticality. For the standard model, the stable point S is located at x = xs = 1, and the coe?cient functions (drift velocity derivatives) in the transport equations are of the form v0 (x) = ux (x, 0) = x ? x3 2v2 (x) = ? 2 ux /?y 2 (x, 0) = ?2αx u1 (x) = ?uy /?y (x, 0) = ??(1 + x2 ). (35) (36) (37)

These may be substituted into the Riccati equation (18) for the transverse second derivative w2 (x) = ? 2 W/?y 2 (x, 0), and the equation numerically integrated. As noted, the appropriate initial condition is w2 (x = xs ) = 2|u1(xs )|, i.e., w2 (x = 1) = 4?. Consider the case ? = 1 (the subject of Fig. 2), in particular. One ?nds for all α in the range 0 < α < 4 that w2 is positive on the line segment between x = xs and the saddle at x = 0. Since the WKB tube centered on the axis, which is formed by small transverse ?uctuations about the MPEP, has variance ? ?/w2 (x), this positivity implies that the tube is everywhere well-de?ned. One also ?nds that w2 → 2, i.e., w2 → 2|u1(0)|, as the saddle is approached. Moreover, by integrating the transport equation (23) one ?nds that k0 (x) = K (x, 0) tends to a ?nite, nonzero limit as x → 0. As we explained above, these two conditions are precisely what is needed to ensure Arrhenius weak-noise asymptotics, with an MFPT prefactor proportional to K (0, 0)?1. The bifurcation transition present in the ? = 1 standard model at α = 4 is re?ected in the behavior of w2 and k0 as x → 0. When α = αc = 4, equations (18) and (23) can be solved exactly; one ?nds w2 (x) = ? 2 W/?y 2 (x, 0) = 4x2 , k0 (x) = K (x, 0) = 1/x. 19 (38) (39)

12 10 8 6 4 2 0 1 2 3 4 5

Figure 5: A plot of K (0, 0), to which the weak-noise activation rate prefactor is proportional, ? as a function of the o?-axis softening parameter α. As shown, K (0, 0) diverges as α → αc . This is for the standard double well model (6), with ? = 1, for which αc = 4. We know from the Eyring formula that the activation rate prefactor, in the limit of weak noise, is proportional to K (0, 0). The fact that k0 (x = 0) = K (0, 0) is in?nite here strongly suggests that at criticality the activation rate, i.e., the rate at which the quasistationary density is absorbed on the separatrix, is anomalously large . Equivalently, it suggests that 1 at criticality the weak-noise behavior of the MFPT (which is asymptotically equal to λ? 1 ) is non-Arrhenius , with a pre-exponential factor that tends to zero as ? → 0. There is an even stronger piece of evidence that this is the case. It is not di?cult to show, by analysing the ? transport equation (23), that K (0, 0) = k0 (x = 0) ? (αc ? α)?1/2 as α → αc . Figure 5 shows the result of a numerical computation when ? = 1. The activation rate prefactor diverges ? as α → αc . Equivalently, the MFPT prefactor tends to zero. The natural deduction is that at criticality the activation rate prefactor, and its reciprocal the MFPT prefactor, become ?-dependent . This blends nicely with the behavior above the transition, since (as shown in Fig. 3) the exponential growth rate of the MFPT (the action barrier ?W , i.e., W (0, 0)) begins to decrease as α increases beyond αc . An ?-dependent activation rate prefactor at criticality, containing a negative power of ?, would unify the exit behavior both below and above criticality. We shall show in Section 8 that in critical double well models (e.g., in the standard model with α = αc ), the weak-noise activation rate λ1 = λ1 (?) indeed has asymptotics λ1 (?) ? const × ??s exp [??W/?] , ? → 0, (40)

where s > 0 is the singularity index mentioned in Section 4. The computation of the 1 ?s singularity index is nontrivial. Since the MFPT τ satis?es τ ? λ? prefactor 1 , the ? 20

in λ1 gives rise to an ?s MFPT prefactor. In critical double well models, in the weak-noise limit the growth of the MFPT is slower than pure exponential . At criticality, the action W as well as the prefactor K displays unusual behavior at the saddle. We shall see that the behavior of (38), i.e., that w2 tends to zero quadratically as x → 0, is universal. Since the WKB tube has variance ? ?/w2 (x), this implies that the tube splays out as the saddle is approached. To leading order, it splays out to in?nite width. This is an indication that the transverse ?uctuations around the MPEP, on the O (?1/2 ) lengthscale, at criticality become very strong near the saddle. In fact, that w2 (x = 0) = 0 causes some di?culty in the interpretation, near the saddle, of the WKB approximation to ρ0 and ρ1 . One might expect that even though w2 (x = 0) = 0, the quartic tube approximation ρ0 (x, y ) ? k0 (x) exp ? w0 (x) + w2 (x) y4 y2 /? , + w4 (x) 2! 4! ?→0 (41)

would su?ce for an understanding of the behavior of the WKB approximation near the saddle. If w2 (x = 0) were zero but w4 (x = 0) were ?nite and nonzero, transverse ?uctuations around the saddle would be, by (41), of magnitude O (?1/4 ) rather than O (?1/2 ). However, explicit solution of eq. (25), the transport equation for w4 = ? 4 W/?y 4 = ? 3 py /?y 3 , shows that in the ? = 1 standard model at α = αc = 4, w4 (x) ? (4/5)x?4 + (16/5)x?2 + 8 + · · · , x → 0+ . (42)

The fact that w4 (x = 0) is in?nite , coupled with the fact that w2 (x = 0) is zero, suggests that at criticality, the transverse ?uctuations near the saddle have no natural scale . In any event, at criticality the standard matched asymptotic approximations technique of the last section breaks down. We shall need to construct an approximation to the quasistationary density ρ1 near the saddle which (i) is valid at criticality, and (ii) matches to the WKB approximation (41), despite its singular character. Critical exponents , as we de?ne them, describe the weak-noise behavior of a parametrized double well model with a singularity at some point x0 , as x → x0 and as the parameters of the model tend to the values for which the singularity appears at x0 . In particular, they characterize the behavior at and near the bifurcation transition, and at and near the saddle point, of the functions W and K appearing in the WKB approximation to ρ0 and ρ1 . At criticality, the divergence rates of the WKB prefactor K as x → 0 and y → 0 supply two such exponents; the scaling form which we shall use to approximate W near a nascent cusp (which involves fractional powers) will supply others. There are also critical exponents describing what happens as one moves o? criticality. As α is increased above αc (i.e., above 4, in the ? = 1 standard model), the MPEP bifurcates. There is a critical exponent describing the separation rate of the two resulting MPEP’s, as Fig. 3 makes clear. There is also a ? critical exponent describing the divergence rate of K (0, 0) as α → αc , which as we have already noted equals 1/2. The singularity index s can be regarded as a critical exponent too, though of a di?erent kind; to compute it, one must go beyond the WKB approximation. The ‘nascent cusp,’ as a singularity, is located in a space parametrized by x, y , and the parameter(s) of the drift ?eld u. But if one restricts oneself to a single double well 21

model, the only parameters are x and y . In this case there is a natural analogy between the Lagrangian manifold M in phase space, formed by bicharacteristics emanating from (S, ?) and (S ′ , ?), and a thermodynamic surface . The action W , as a function on the manifold, corresponds to a thermodynamic potential, in fact a Gibbs free energy. The equation p = p(x) = ?W/? x corresponds to a relation between conjugate state variables, such as pressure and volume. The singularities (points of non-di?erentiability) of the physical action W (x) therefore correspond to phase transitions . The order of such a phase transition, in the traditional sense, is the lowest order of spatial derivative (of W ) which fails to be continuous. In Section 8 we shall compute the order of the nascent cusp.

5.4

The Bifurcation Transition and Classical Mechanics

We now explain how the equations of Section 5.2 allow the bifurcation transition to be interpreted in terms of classical mechanics, and how one can predict whether or not any given double well model is at criticality. We begin by considering models in which the MPEP has already bifurcated, and the instanton trajectories emanating from the stable point S focus along the axis, as in Figs. 2(c) and 2(d). Empirically, focusing occurs in models with a su?ciently large ‘o?-axis softening’ parameter α. In such models, the action W in a mathematical sense becomes multivalued near a portion of the axis. Points (x, 0) beyond the ?rst focus are reached by multiple o?-axis instanton trajectories emanating from the stable point S , and in general these trajectories will have di?erent actions. They will also have di?erent momenta p = ?W at the time they reach (x, 0). This multivaluedness has a geometric interpretation, in terms of the shape of the twodimensional Lagrangian manifold (in the four-dimensional phase space) formed by the bicharacteristics emanating from the point (x, p) = (S, ?). As x decreases from xs toward zero, the map y → py in the vicinity of y = 0 is at ?rst single-valued; the value py = 0, and no other, (1) (1) corresponds to y = 0. Beyond the ?rst focus (xf , 0), i.e., when x < xf , the map y → py (2) becomes three-valued. At the second focus (xf , 0) it becomes ?ve-valued, etc. The generic evolution is shown in Figs. 6(a) through 6(f). Up to the ?rst focus y = y (py ) near py = 0 may be modelled as a linear function; beyond the focus, as a cubic . Beyond the second focus the global description becomes more complicated, as is clear from the whorl in Fig. 6(f). A cubic approximation is still appropriate in the immediate vicinity of (y, py ) = (0, 0), however. Since the locus of all points (y, py ) at constant x is obtained by intersecting the Lagrangian manifold with the hyperplane x = const, the manifold itself becomes increasingly ‘whorled’ with each passage through a focus. The formation of convolutions in Lagrangian manifolds was ?rst considered by Berry and Balazs [5] (in a time-dependent context), and the progression in Fig. 6 resembles the ?gures in their paper. Geometrically, the linear-tocubic transition at each successive focus corresponds to the creation of a fold [12]. One can (l) ?t the shape of the manifold near the l’th focus, i.e., near (x, y, py ) = (xf , 0, 0), by the

22

p sub y

p sub y

y

y

p sub y

p sub y

y

y

p sub y

p sub y

y

y

Figure 6: Cross-sections through the Lagrangian manifold M, revealing the ‘whorling’ that takes place as one passes through any on-axis focus. These sketches show the map y → py (1) at successively decreasing values of x, as one moves from S = (xs , 0), past two foci (xf , 0) (2) (1) and (xf , 0), toward the saddle point (0, 0). Shown are the cases (a) xs > x > xf and (2) (2) (1) (1) (1) (1) x near xs , (b) xs > x > xf and x near xf , (c) x = xf , (d) xf > x > xf , (e) x = xf , (2) and (f) xf > x > 0.

23

phenomenological formula y = y (x, py ) = ?a0 p3 y ? a1

(l) (l) (l) (l)

x ? xf

(l)

py

(43)

where a0 and a1 are certain positive constants. So each successive focus resembles a Ginzburg-Landau second-order phase transition (x corresponding to temperature, the focus (l) location xf to a critical temperature, ?y to a magnetic ?eld, and py to a magnetization). We shall say more about the ‘equation of state’ (43) (which we stress is not applicable near the ‘nascent cusp’ appearing at the saddle point of critical models) in Sections 7 and 8. For on-axis (i.e., y = 0, py = 0) trajectories, the derivative w2 (x) = ?py /?y (x, y = 0) satis?es the Riccati equation (24). So the appearance of a focus, and of multiple foci, can be investigated analytically. It is clear from Figs. 6(c) and 6(e) that passage through a focus is signalled by the tangent plane to the manifold (at y = 0, py = 0) ‘turning vertical’; equivalently, by ?y/?py passing through zero, or its reciprocal w2 (a negative magnetic susceptibility, in this context) passing through ?∞. To study this, recall that the Riccati equation 2 w ˙ 2 = ?w2 ? 2u1 w2 + 4v0 v2 (44) involves a derivative with respect to instanton transit time, and that the on-axis instanton trajectory (directed anti-parallel to the drift toward S ) satis?es x ˙ = ?v0 (x). Solutions w2 can be regarded either as a function of t, for ?∞ < t < ∞, or of x, for xs > x > 0. We see from the form of the Riccati equation that w2 can indeed be driven to ?∞ in ?nite time, i.e., at some point x = xf > 0 to the right of the saddle. In fact one sees, if tf is the time + when this occurs (the focus time), that as t → t? f , i.e., x → xf , w2 (x) = ? 2 W/?y 2 (x, 0) ? ?(tf ? t)?1 , ? const × [?(x ? xf )?1 ]. (45)

Here the constant multiplier equals 1/v0 (xf ), the reciprocal speed of the on-axis instanton trajectory when it passes through the focus (xf , 0). We note in passing that by the transport equation (23), this blowup will induce a blowup of the on-axis WKB prefactor k0 . One ?nds k0 (x) = K (x, 0) ? const × (tf ? t)?1/2 , ? const × (x ? xf )?1/2 .

(46)

Equations (45)–(46) contrast markedly with eqs. (38)–(39), which apply to the ? = 1 standard model at criticality (where, in a formal sense, xf = 0, since there is a nascent focus at the saddle). Equations (45)–(46) are not restricted to the standard model; they hold in greater generality. But they apply only when a bona ?de focus is present at some xf > 0 (i.e., before the saddle is reached), and the MPEP has already bifurcated. By examining the Riccati equation (44), we see that w2 will be driven to ?∞, and a focus will be present, only if the inhomogeneous term 4v0 v2 on the right-hand side of (44) is su?ciently negative. (This is because u1 < 0, by assumption.) But 2v2 = ? 2 ux /?y 2 (x, y = 0), 24

w2 4

2

0

-2

-4

0

0.2

0.4

0.6

0.8

1

x

Figure 7: The transverse action derivative w2 (x) ≡ ? 2 W/?y 2 (x, 0) in the ? = 1 variant of the standard double well model (6), in the vicinity of the bifurcation transition. The three curves, from top to bottom, obtain when α = 3.9, 4.0 (the critical value αc ), and 4.1. When α < αc , w2 → 2|λy |, i.e., 2, as x → 0+ . When α > αc , w2 is driven negative as x decreases, and passes through ?∞. The focal point x = xf where this occurs, in the model with α = 4.1, is indicated by a dashed line. which we are taking to be negative when 0 < x < xs , measures the extent to which the drift toward S softens as one moves o?-axis. So our empirical observation is con?rmed analytically: a su?ciently strong o?-axis softening will create a focus, and a bifurcation of the MPEP! It is best to think of w2 = ?py /?y (y = 0) as a slope , as in Fig. 6. As such, it may rotate repeatedly through the point at in?nity as t increases, i.e., as x decreases. Each such rotation results in increased whorling of the Lagrangian manifold, and also corresponds to a passage through a focus. So by counting the number of singularities of the solution curve w2 = w2 (t), one may determine the number of foci present in any given double well model. The standard model (6) will serve as an example. For the reasons discussed in Section 5.2 w2 (t = ?∞), i.e., w2 (x = xs ), in the standard model always equals 2|u1 (xs )| = 2|u1 (1)| = 4?. Suppose that ? = 1. We noted in the last section that if 0 < α < 4, w2 is well-behaved and positive at all times t between ?∞ and ∞ inclusive, i.e., at all x satisfying 0 ≤ x ≤ 1. We also explained what happens at α = αc = 4, when the nascent cusp appears at the saddle and the MPEP begins to bifurcate. At criticality, w2 → 0 as t → ∞, i.e., as x → 0+ . If 4 < α < 9, w2 is driven negative (as t increases toward ∞), and passes through ?∞ before 25

returning (through +∞) to ?nite, positive values. The change in behavior is shown in Fig. 7. When α is raised above αc , we say that the graph of w2 acquires unit winding number , since (2) it winds once through the point at in?nity. A second transition occurs at α = αc = 9. If 9 < α < 16, w2 passes through ?∞ twice , and its graph has winding number equal to 2. (j ) Except at the critical values α = αc = (j + 1)2 , w2 in this model converges to the generic value 2|u1(0)| = 2 |λy | = 2 as t → ∞, i.e., as x → 0+ . Since each passage of w2 through ?∞ gives rise to a focus, the sequence of near-axis instanton ?ow ?elds in the ? = 1 standard model, as α is increased, displays an progressively larger number of foci. In fact the progression is precisely as displayed in Figs. 2(a) through 2(d). It is worth noting that in models with one or more foci, the WKB tube centered on the axis becomes ill-de?ned when w2 goes negative, which takes place at a location on the axis somewhat before the ?rst focus is reached. (See Fig. 7.) In Section 8 we shall determine exactly what happens at the bifurcation transition of any singly parametrized symmetric double well model. But we can now pose the question: What, physically, causes the above values for α to be critical? If in general the odd function 2v2 (x) = ? 2 ux /?y 2 (x, 0) is negative between x = 0 and x = xs and is proportional to a parameter α, is there a classical mechanical technique of predicting the values of α at which the on-axis instanton trajectory will bifurcate? The answer to this question is ‘yes.’ Our technique relies on a linear stability analysis of the on-axis instanton trajectory, and identi?es the critical values of α as the values for which a transverse soft mode is present in the zero-energy Hamiltonian dynamics. This is reminiscent of Langer’s analysis of metastability in one-dimensional models [24, 36]. But because we shall consider transverse, rather than longitudinal, ?uctuations around the instanton trajectory, our stability analysis will be considerably simpli?ed. Let x = x? (t) = (x? (t), 0) be the on-axis instanton trajectory, where x = x? (t) is the solution of x ˙ = ?v0 (x). Near-axis instanton trajectories, i.e., near-axis zero-energy classical trajectories emanating from S , may to leading order be written as x = x(t) = (x(t), y (t)) ? (x? (t), 0) + δ (0, Y (t)) (47)

where δ ? 1, and where Y = Y (t) is some model-dependent function satisfying Y (t = ?∞) = 0. Y (t), ?∞ < t < ∞, is a normalized transverse deviation . Similarly, near-axis trajectories have momenta p = p(t) = (px (t), py (t)) ? (p? x (t), 0) + δ (0, Py (t)) (48)

for some unknown function Py = Py (t) satisfying Py (t = ?∞) = 0. Here p = p? (t) = (p ? x (t), 0) is the momentum of the on-axis instanton trajectory at instanton transit time t. ? ? We noted before eq. (22) that as a function of x, p? x equals ?2v0 . So px (t) equals ?2v0 (x (t)). We necessarily have w2 (t) = Py (t)/Y (t), (49) on account of w2 equalling ?py /?y (y = 0). 26

Substituting equations (47) and (48) into the Hamilton equations derived from the Wentzell-Freidlin Hamiltonian H , and separating terms proportional to δ , yields the pair of equations ˙ Y = ?2H ? ?2H (x (t), p? (t)) Py (x? (t), p? (t)) Y + ?py ?y ?p2 y (50) (51)

2 2 ˙y = ? ? H (x? (t), p? (t)) Y ? ? H (x? (t), p? (t)) Py . P ?y 2 ?py ?y

Due to the special form of the Hamiltonian (12), and the expansions (19), this pair becomes ˙ = u1 (x? (t)) Y + Py Y ˙y = 4v0 (x? (t)) v2 (x? (t)) Y ? u1 (x? (t)) Py . P (52) (53)

So w2 = w2 (t) can be represented as the quotient of two functions of instanton transit time, which satisfy a pair of coupled linear di?erential equations. We note in passing that an analogous representation is possible for solutions Z = Z (t) of the matrix Riccati equation (18). The existence of such quotient representations is well known in the theory of Riccati equations, and has a geometric interpretation [37]. Equations (50)–(51), and (52)–(53), may be viewed as Hamilton’s equations for the effective (i.e., time-dependent) transverse Wentzell-Freidlin Hamiltonian He? (Y, Py , t) = 1 ?2H 2 ?2H 1 ?2H 2 P + Y P + Y y y 2 ?p2 ?py ?y 2 ?y 2 y

2 = Py /2 + u1 (x? (t)) Y Py ? 2v0 (x? (t)) v2 (x? (t)) Y 2 .

To save space we have suppressed the arguments (x? (t), p? (t)) of the partial derivatives. This quadratic Hamiltonian governs the small transverse ?uctuations about the on-axis instanton ˙ Py ? He? , namely trajectory. Its Legendre transform Y ˙ , t) = Le? (Y, Y ?2L ˙ 1 ?2L 2 1 ?2L ˙ 2 Y + Y Y + Y 2 ?y ˙2 ? y?y ˙ 2 ?y 2 ˙ ? u1 (x? (t)) Y |2 /2 + 2v0 (x? (t)) v2 (x? (t)) Y 2 , = |Y

is an e?ective transverse Onsager-Machlup Lagrangian . Here L is the Onsager-Machlup La˙ ? (t)) of the partial derivatives. grangian (13), and we have suppressed the arguments (x? (t), x The corresponding Euler-Lagrange equation for the normalized transverse deviation Y , i.e.,

? ? ? ¨ + ? d [u1 (x? (t))] ? u2 Y 1 (x (t)) ? 4v0 (x (t)) v2 (x (t)) Y = 0, dt

(54)

is called a (transverse) Jacobi equation [36]. It may be written as an equation for Y = Y (x), 0 ≤ x ≤ xs , by changing the independent variable from t to x. One gets JY ≡ d dY u2 (x) v0 (x) + u′1 (x) ? 1 ? 4v2 (x) Y = 0 dx dx v0 (x) 27 (55)

together with the boundary condition Y (x = xs ) = 0. This Jacobi equation, which is in Sturm-Liouville form, governs the behavior of the instanton trajectories near the on-axis trajectory (the MPEP, if it has not bifurcated). So it is responsible for the various behaviors shown in Fig. 2. It is clear from our derivation that the Jacobi operator J , considered as a quadratic form, de?nes the transverse second variation of the Onsager-Machlup action functional about the on-axis trajectory. Foci are by de?nition the points (x, 0) where the o?-axis instanton trajectories converge (to leading order). Equivalently, (x, 0) is a focus only if Y (x) = 0. But since w2 = Py /Y , this implies that (unless Py (x) = 0 also) w2 (x) is in?nite. This is precisely the necessary condition for a focus that we derived earlier. If Y passes through zero more than once, then w2 will pass through the point at in?nity more than once. This is the mechanism by which, e.g., the two-focus ?ow ?eld of Fig. 2(d) engenders the increasingly ‘whorled’ Lagrangian manifold of Figs. 6(a) through 6(f). We can now give a simple criterion for determining whether or not a given double well model is at criticality. Suppose that the most probable escape path (MPEP) extends along the axis from S to the saddle, so that the symmetry is as yet unbroken. We know by the discussion in Section 4 that criticality is signalled by the appearance of a nascent cusp at the saddle. The nascent cusp itself is not a focus, as Fig. 2(b) makes clear. But if the o?-axis softening is increased, the nascent cusp becomes a genuine cusp (i.e., focus); it moves inward along the axis from the saddle toward S . This picture is consistent with the interpretation of the near-axis instanton ?ow ?eld in terms of the function Y (x), 0 ≤ x ≤ xs , only if the nascent cusp, like a conventional on-axis focus, is a zero of Y . So the signal for criticality is Y equalling zero at x = 0. We can rephrase this as follows. Critical double well models are those models with unbroken symmetry for which the Jacobi equation J Y = 0 for the transverse deviation function Y , equipped with boundary condition Y (x = xs ) = 0 and also with Y (x = 0) = 0, has a nontrivial ( i.e., nonzero) solution. The nonzero solution Y = Y1 (x), 0 ≤ x ≤ xs , when it exists, can be interpreted as a transverse soft mode of the zero-energy Hamiltonian dynamics. If the o?-axis softening is increased, the on-axis MPEP will bifurcate. Just beyond criticality, there will be two symmetrically placed o?-axis MPEP’s from S to the saddle. They will be of the form (x? (t), ±δY1 (x? (t))), for some small δ . This ‘motion in the direction of a soft mode’ is a standard bifurcation e?ect. At criticality, the transverse soft mode Y1 describes the way in which the two MPEP’s separate. Suppose that the double well model is parametrized by an o?-axis softening parameter α, i.e., that v2 = αv ?2 for some odd function v ?2 , and that v0 and u1 are independent of α. Then by rewriting the Jacobi equation, one sees that the model will be at a bifurcation point if and only if the Sturm-Liouville equation ?Y ≡ J 1 dY 1 u2 (x) d v0 (x) + u′1 (x) ? 1 Y = αY, 4? v2 (x) dx dx 4? v2 (x) v0 (x) (56)

equipped with Dirichlet boundary conditions Y (x = 0) = Y (x = xs ) = 0, has a nonzero ? may be called a normalized Jacobi operator . solution. The Sturm-Liouville operator J 28

We see that the set of critical values of α is precisely the spectrum of the normalized (1) Jacobi operator! Only the ?rst critical value (i.e., lowest eigenvalue) αc = αc will yield (j ) an actual bifurcation of the MPEP. To each higher critical value αc , j = 2, 3, . . ., there corresponds a transverse eigenmode Yj . But after the ?rst bifurcation, the on-axis instanton trajectory is no longer the physical MPEP. The higher eigenmodes Yj , j = 2, 3, . . ., which are oscillatory, govern the further bifurcations of the on-axis instanton trajectory rather than the further bifurcations (if any) of the physical MPEP’s, which have already moved o?-axis. The case of the standard model (6) is instructive. Substituting from (35)–(37) one ?nds as normalized Jacobi operator

2 2 2 ? = ? 1 d (x ? x3 ) d + ? + ? (1 + x ) . J 4x dx dx 2 4x2 (1 ? x2 )

(57)

It is easily veri?ed that on the interval from x = 0 to x = xs = 1, this operator (when equipped with Dirichlet boundary conditions) has spectrum

(j ) αc = j 2 + (3? ? 1)j + (2?2 ? ?),

j = 1, 2 , 3 , . . .

(58)

(1) So in the standard model, the bifurcation of the physical MPEP occurs at αc = αc = (j ) 2 2?(? + 1). Also, the standard model with ? = 1 has αc = (j + 1) , so the on-axis instanton trajectory bifurcates at α = 4, 9, 16, . . . We have several times mentioned this (j ) curious progression of squares. The eigenfunctions Yj corresponding to the eigenvalues αc , (j ) i.e., the transverse soft modes appearing at α = αc , turn out to be of the form

Yj (x) = (x ? x3 )? qj (x),

(59)

where qj is an even polynomial of degree 2j ? 2. Substituting into the transverse Hamilton equation (52) yields the analogous transverse momentum deviations (Py )j . One gets (Py )j (x) = ?v0 (x)Yj′ (x) ? u1 (x)Yj (x)

(j ) , so that in the standard model at α = αc ′ w2 (x) = (Py )j (x)/Yj (x) = 4?x2 ? (x ? x3 )qj (x)/qj (x).

′ (x) , = (x ? x3 )? 4?x2 qj (x) ? (x ? x3 )qj

(60)

If j = 1 then qj reduces to a constant, and the second term is absent. So at the physical bifurcation point [i.e., at α = αc = 2?(? + 1)], w2 (x) equals 4?x2 . Moreover, Y1 (x) = (x ? x3 )? . This transverse soft mode is seen clearly in Figs. 2(c) and 2(d), which show the behavior of the ? = 1 standard model beyond the bifurcation point. In those ?gures the o?-axis MPEP’s are roughly proportional to ±(x ? x3 ), i.e., to ±Y1 . As α is increased above αc , the MPEP’s move in the direction of the transverse soft mode. Recall that the pro?le of the WKB tube of probability density centered on the x-axis is asymptotically Gaussian, and that at speci?ed x this Gaussian has variance ? ?/w2 (x). But 29

in the standard model, at the ?rst (and only physical) critical value α = αc = 2?(? +1), w2 (x) equals 4?x2 . That w2 → 0 as x → 0 implies that at criticality, the WKB tube splays out as the saddle is approached. We have already seen the ? = 1 case of this in Section 5.3. The splayout is what one would expect from our picture of the bifurcation of the MPEP, which begins at the saddle, as a phase transition. It simply says that on the O (?1/2 ) transverse lengthscale, the Gaussian ?uctuations about the MPEP grow without bound as the nascent cusp is approached. It is easy to see that this behavior is universal : it occurs in any critical double well model with a bifurcating MPEP. If the (diagonal) linearization of the drift ?eld u at the saddle has eigenvalues (λx , λy ), and ? is de?ned as usual to equal |λy |/λx , then examination of the Jacobi equation shows that the soft mode Y1 has asymptotics Y1 (x) ? Cx? , x → 0+ , for some nonzero constant C . We have mentioned this ‘approach path’ property elsewhere [28]. Also, examination of the Hamilton equation for (Py )1 shows that (Py )1 (x) ? C ′ x?+2 for some nonzero C ′ . So at criticality, the quotient w2 (x) satis?es (for any ?) w2 (x) = ? 2 W/?y 2 (x, 0) ? const × x2 , (61)

as x → 0+ , and the tube splayout always occurs. Incidentally, it follows by integrating the transport equation (23) that k0 (x) = K (x, 0) ? const × x?? , (62)

as x → 0+ . Equations (61)–(62) summarize the universal behavior of the WKB tube near the saddle, in any critical double well model. They are the extension to arbitrary critical models of eqs. (38)–(39), which applied only to the critical variant (α = αc = 4) of the ? = 1 standard model. We stressed in Section 5.2 that a Kramers-type error function approximation to the quasistationary density ρ1 near the saddle is appropriate only if w2 → 2|u1 (0)| as the saddle is approached. At criticality, since w2 → 0 instead, in order to apply the method of matched asymptotic approximations we shall need to construct a di?erent boundary layer approximation. This will give rise to the universal non-Arrhenius MFPT asymptotics for models at criticality.

6

Maslov-WKB Asymptotics

By building on the previous sections, we can analyse the weak-noise behavior of double-well models with singularities. We have seen that singularities may appear in the WKB approximation K (x) exp[?W (x)/?] for the stationary density ρ0 and quasistationary density ρ1 . The possible singular behaviors are summed up in eqs. (45)–(46), which apply to models in which the MPEP has already bifurcated, and eqs. (61)–(62), which apply to models which are critical in the sense of bifurcations. Models in which the MPEP has already bifurcated have the property that the instanton trajectories emerging from S focus at a point (xf , 0) on the axis, with xf > 0. The prefactor K of the WKB approximation will diverge there. 30

In critical models, there is no actual on-axis focusing. But the prefactor will nonetheless diverge at the saddle point (0, 0). There is a standard procedure for extending the WKB approximation to such singular points, by ‘glueing in’ auxiliary, non-WKB approximations. It originated with the work of Keller and Rubinow on short-wave asymptotics [20], and has been most extensively developed by Maslov [32]. For a mathematically rigorous treatment, see Duistermaat [10]. See also Eckmann and S? en? eor [13], for a partly pedagogical one-dimensional treatment. The procedure may be applied to the (formal) asymptotic solutions of any partial di?erential equation of the form H (x, ???)ρ = 0, where H is a speci?ed Hamiltonian. Here we discuss its application to the Smoluchowski equation, in arbitrary dimensionality n. We know from Section 5.1 that mathematically, the WKB approximation to ρ0 and ρ1 is determined by (i) a Lagrangian manifold M in the 2n-dimensional phase space, formed by the bicharacteristics emanating from (S, ?) and (S ′ , ?), and (ii) functions W and K de?ned on this manifold, and computable by integration along the bicharacteristics. Of the points P (i) = (x, p(i) ) ‘over’ any point x, only the one with least action is physical. The values there of W and K are the values W (x) and K (x) appearing in the WKB approximation. This geometric interpretation motivates the introduction of a new, ‘di?raction integral’ way of formulating the WKB approximation. At any point P = (x, p) on the Lagrangian manifold M, we have W (P ) = p · dx, (63) the line integral being taken along the bicharacteristic terminating at P . We can de?ne a Legendre transform W , satisfying W = x · p ? W , by W (P ) = x · d p. (64)

It is natural to think of W as a function of momentum p, by projecting ‘sideways’ onto momentum space. Of course W (p) is potentially multivalued, like W (x). For W , it is the least of the possible values that is physical; for W , it is the most . But if one ignores the multivaluedness of W (p), one can write K (x) exp[?W (x)/?] ? ??n/2 ··· K (p) exp ?x · p + W (p) /? dp1 · · · dpn , (65)

where K (x) and K (p) are related by K (x) ∝ K (p) det ? ?2W (p) = K (p) ?pi ?pj det ? ?2W (x ) , ?xi ?xj (66)

the correspondence between p and x being given by p(x) = ?W/? x, or x(p) = ? W /? p. The asymptotic equality in (65), as ? → 0, is justi?ed by the method of steepest descent. (It may be necessary to cut o? the integral at large momentum to ensure convergence.) The method of steepest descent automatically picks out the point P = (x, p) ‘over’ x with the 31

least action W (P ). We shall call (65) a di?raction integral representation , since (if ? is pure imaginary) it resembles the di?raction integrals used in physical optics [4]. We have assumed that the Hessian matrix ? 2 W/?xi ?xj = ?pi /?xj , whose inverse is the matrix ? 2 W /?pi ?pj = ?xi /?pj , is negative de?nite. Actually it is often possible to make sense of the above formul? even when this is not the case, by analytic continuation. It is also possible to avoid the problem of positive eigenvalues by taking the Legendre transform with respect to a partial (incomplete) set of variables. In n = 2 dimensions, this means with respect to a single variable only. For example, one could use the alternative integral representation K (x, y ) exp[?W (x, y )/?] ? ??1/2 K (x) (px , y ) exp ?xpx + W (x) (px , y ) /? dpx , (67)

where W (x) = xpx ? W is regarded as a function of px and y , and K and K (x) are related by K (x, y ) ∝ K (x) (px , y ) ? ? 2 W (x) (px , y ) = K (x) (px , y ) ?p2 x ? ?2W (x, y ). ?x2 (68)

Here the correspondence between (x, y ) and (px , y ) is given by px (x, y ) = ?W/?x, or equivalently x(px , y ) = ? W (x) /?px . It is clear that the transformed prefactor K (resp. K (x) , etc.) in these integral representations, like W , K , and W (resp. W (x) , etc.), can be thought of as a function on the Lagrangian manifold M. Also, the momentum integration can be viewed as an integration over M. So introducing integral representations of this sort is really a way of replacing the position-space WKB approximation K (x) exp[?W (x)] to ρ(x) by a smeared-out equivalent one, or ones, involving integration over the manifold. As derived, these ‘momentum space’ approximations are accurate only to leading order as ? → 0, since subdominant terms in ? arising from the method of steepest descent have been neglected. But such terms could be incorporated, if desired, by adding ?-dependent corrections to the transformed prefactor. If the new formulations of the WKB approximation are equivalent to the old, why have we introduced them? The reason is that the equivalence holds only at points x at which K is ?nite. At singularities of K , the new formulations provide a means of computing the true ? → 0 asymptotics of ρ. Moreover, they reveal how at least some singularities of K can be explained as artifacts , arising from the way in which K is computed from K . It follows from (66) that if the determinant of the Hessian matrix ? 2 W/?xi ?xj = ?pi /?xj diverges at some point x, then x will be a singularity of K whenever K is nonzero at the corresponding momentum p = p(x). In other words, singularities of K may be more apparent than real: they can arise from points (x, p) on the manifold where K does not actually diverge. A similar e?ect can arise from the representation (67), or from any other di?raction integral representation. The matrix ?pi /?xj is a matrix of partial slopes, which speci?es (to ?rst order) the shape of the manifold in the vicinity of the point (x, p) = (x, p(x)). Its determinant becomes in?nite only when at least one of its elements is in?nite. Such a blowup occurs only at 32

locations on the manifold where the (n-dimensional) tangent hyperplane to the manifold ‘turns vertical,’ i.e., points along a momentum direction in the 2n-dimensional phase space. This is precisely the behavior one sees at a fold , as in Figs. 2(c), 2(d), and 6. The folding over of the Lagrangian manifold can create singularities of K . This is clearly the cause of the singular behavior of K at the on-axis foci occurring in models with a bifurcated MPEP, though not of the singular behavior at the nascent cusp occurring at criticality (which cannot be transformed away). Whether or not a singularity of the prefactor K occurring at some point x = x? is an artifact of this sort, by employing an appropriate di?raction integral representation one may compute the true weak-noise asymptotics of ρ(x? ). One usually ?nds leading-order behavior of the form const × ??σ exp[?W (x? /?)], where σ is by de?nition the singularity index of x? . In fact the ?-dependent prefactor ??σ should appear at all points x within some ?-dependent distance of x? , which shrinks to zero as ? → 0. Within this local region an asymptotically exact formula for ρ, derived from the integral representation, will be uniformly valid. This asymptotic approximation (non-WKB, at least in the traditional sense) will match in the far ?eld to the WKB approximation K (x) exp[?W (x)/?]. In Sections 7 and 8 we shall see how this ‘glueing in’ procedure works, both in models with a bifurcated MPEP and in models at criticality. For the moment we note only that the construction of a local approximation to ρ, near the singular point x? , depends crucially on the determination of the behavior of W and K near the corresponding point p? in momentum space. The case when K is well-behaved (‘slowly varying’) in a neighborhood of p? , and the singularity at x = x? is an artifact, is the simplest. Suppose that W = W (p) can be expanded in a power series around p = p? . The matrix ? 2 W /?pi ?pj must have a zero eigenvalue at p = p? , since otherwise the determinant of its inverse ? 2 W/?xi ?xj would not tend to in?nity as x → x? , the Lagrangian manifold would not turn vertical there, and the singularity in K would not appear. The term catastrophe is used to describe what happens to the manifold at x = x? . It is a standard result, due largely to Arnol’d [2], that if the manifold is smooth near (x? , p? ), the catastrophic behavior at x = x? can be captured by approximating W = W (p) by one of a handful of polynomial functions. These are the ‘structurally stable’ elementary catastrophes . A single example, illustrating the similarity to the Ginzburg-Landau theory of phase transitions, will su?ce. In n dimensions, suppose that a singularity at x = x? arises as an artifact in the above sense, and that p? = ?W/? x(x? ). In appropriate (linearly transformed) coordinates, write z = (z1 , . . . , zn ) = x ? x? g = (g1 , . . . , gn ) = p ? p? . (69) (70)

A particularly common sort of catastrophe (a ‘cuspoid’) would be described locally by a single-variable Legendre transform of the form W (zn ) (z1 , . . . , zn?1 , gn ) = ?

n+2 n 2 a0 gn a1 z1 gn an?1 zn?1 gn ? ?...? + R(z1 , . . . , zn?1 ) n+2 n 2

(71)

33

where a0 , . . . , an?1 are constants, and R(z1 , . . . , zn?1 ) is a quadratic polynomial. Since zn = ? W (zn ) /?gn , this expression implies

n+1 n?1 zn = zn (z1 , . . . , zn?1 , gn ) = ?a0 gn ? a1 z1 gn ? . . . ? an?1 zn?1 gn .

(72)

The presence of a catastrophe at (z1 , . . . , zn?1 , gn ) = (0, . . . , 0, 0) is signalled by the fact that 2 ? 2 W /?gn = ?zn /?gn equals zero there. We have already seen the n = 2 version of eq. (72) in Section 5.4, as a phenomenological description of the shape of the manifold M near an on-axis focus. Recall that we interpreted eq. (43), which is the n = 2 version, in thermodynamic terms: as the equation of state of a substance undergoing a Ginzburg-Landau second-order phase transition. (E.g., z1 is T ? Tc , z2 is a negative magnetic ?eld, and g2 is magnetization.) Equation (72) is in fact a normal form for the shape of a Lagrangian manifold near a cuspoid singularity. When n = 2, the cuspoid is a cusp . If n = 1, only the ?rst term on the right-hand side of (72) is present, and the cuspoid reduces to a quadratic fold . In general, to each possible polynomial expression (normal form) for the Legendretransformed action, there corresponds a non-WKB approximation to ρ in a local region near x = x? , computed from the appropriate di?raction integral. These integrals serve to de?ne the canonical di?raction functions ?rst explored by Maslov. The canonical di?raction functions include the classical Airy and Pearcey functions, which arise from folds and cusps respectively [4]. We shall study the cusp case further in the next section, as a warmup for the study of the nascent cusp appearing at criticality. The normal form for the action near a nascent cusp will turn out to be nonpolynomial, but the Maslov-WKB technique will still apply. We close this section by noting that di?raction integral representations are also useful for incorporating symmetry constraints and boundary conditions. As an example of this, consider behavior near the saddle point of a double well model. We emphasized in Section 5.1 that if no bifurcation of the MPEP has occurred, the WKB tube of probability density centered on the axis will be well behaved as the saddle is approached. In particular, w2 (x) = ? 2 W/?y 2 (x, 0) will tend to 2|u1 (0)| as x → 0+ . Since ?px ?2W ′ (0, 0) = (0, 0) = ?2v0 (0) 2 ?x ?x (73)

and v0 (x) = ux (x, 0) is assumed to be smooth, in the absence of bifurcations W will to leading order be locally quadratic at the saddle. If u(x, y ) ≈ (λx x, ?|λy |y ) is the linearization of the ′ drift at the saddle, we have u1 (0) = ?|λy | and v0 (0) = λx . So, near (x, y ) = (0, 0), W (x, y ) ≈ W (0, 0) ? λx x2 + |λy |y 2. And

2 W (px , py ) ≈ ?W (0, 0) ? p2 x / 4 λ x + py / 4 | λ y | (x)

(74)

(75) (76)

W

(px , y ) ≈ ?W (0, 0) ? 34

p2 x / 4 λx

? | λy | y

2

will be the leading-order approximations to the Legendre-transformed actions. Since the Hessian matrix ? 2 W /?pi ?pj is not negative de?nite, an integral representation of the type (65) is not appropriate. But a representation of the type (67) may be used. In the absence of bifurcations K and K are well behaved near the saddle, so substituting (76) into (67) yields ρ(x, y ) ? const × ??1/2

?|λy |y exp ?(xpx + p2 x /4λx )/? dpx e

2 /?

.

(77)

But when approximating the quasistationary density ρ1 (x, y ) near the saddle, one needs an approximate solution of the Smoluchowski equation that is odd rather than even. Such an approximate solution is obtained by integrating px from 0 to ∞ rather than from ?∞ to ∞. If this is in fact done, eq. (77) reduces to (30), the standard Kramers-type error function approximation to the quasistationary density! Although error function approximations originated (with Kramers) in an entirely di?erent context, they ?t naturally into the Maslov-WKB framework. We conclude that di?raction integral representations can be modi?ed to incorporate the e?ects of symmetry constraints. In Section 8.3 we shall use a similar half-range integration in our integral representation for the quasistationary density near a nascent cusp.

If px here is integrated from ?∞ to ∞, this approximation will be even in x. It will therefore serve as an approximation to the stationary density ρ0 (x, y ) near the saddle. The integral may be evaluated explicitly, and the approximation reduces to the standard inverted Gaussian approximation 2 2 ρ0 (x, y ) ? const × e+λx x /? e?|λy |y /? . (78)

7

Scaling Behavior Near a Cusp

We can apply the Maslov-WKB method of the last section to symmetric double well models in which the MPEP has bifurcated, and the instanton trajectories emerging from S = (xs , 0) focus at a point (xf , 0), with 0 < xf < xs . As we shall see, behavior near the focal point (xf , 0) is best described in the language of critical phenomena. The Maslov-WKB method was ?rst applied to focusing (cusp) singularities in twodimensional models by Dykman et al. [12]. Their analysis, which does not assume any sort of symmetry, specializes in the case of symmetry about the x-axis to the following. Assume that the Legendre-transformed action W (y) = ypy ? W , regarded as a function of x and py , may be asymptotically approximated near (x, py ) = (xf , 0) by the cuspoid (codimension n = 2) normal form a1 a0 (x ? xf )p2 (79) W (y) (x, py ) ? ? p4 y ? y ? w0 (x). 4 2 Here a0 and a1 are positive constants, and w0 (x) is simply W (x, 0), i.e., ?W (y) (x, 0). Since y (x, py ) = ? W (y) /?py (x, py ), this assumption is equivalent to y (x, py ) ? ?a0 p3 y ? a1 (x ? xf )py 35 (80)

which is the phenomenological (Ginzburg-Landau) equation of state (43), discussed at length in Section 5.4. W (y) = W (y) (x, py ) can be viewed as a Helmholtz free energy, just as W = W (x, y ) can be viewed as a Gibbs free energy. The cuspoid form for W (y) is certainly consistent with the folding of the Lagrangian manifold M, as seen (in projection) in Figs. 2(c) and 2(d). It is also consistent with the quantitative asymptotics of Section 5.4. Since py = 0 corresponds to y = y (x, py ) = 0, (80) implies w2 (x) = ?py ?y ?2W (x, 0) = (x, 0) = (x, 0) 2 ?y ?y ?py 1 ?1 ? ?a? x → x+ 1 (x ? xf ) , f.

?1

(81)

This is precisely the near-focus blowup behavior of eqs. (45)–(46), which we derived analytically from the Riccati equation (44). By comparing (81) with (44), we see that the constant a1 must equal 1/v0 (x = xf ), the reciprocal speed of the on-axis instanton trajectory as it passes through the focus. Since x ? xf is analogous to T ? Tc and w2 to a (negative) magnetic susceptibility, the blowup of (81) is analogous to the critical exponent γ of the focus, in thermodynamic language, equalling unity. Dykman et al. use a one-dimensional di?raction integral representation, resembling (67) but with x and y interchanged, to approximate the stationary probability density ρ0 near x = (xf , 0). A crucial assumption is that the transformed prefactor K (y) = K (y) (x, py ), which has no direct thermodynamic interpretation, is well behaved (locally constant, or ‘slowly varying’) near (x, py ) = (xf , 0). If this is the case, and it may be approximated by a constant, one can construct the Maslov-WKB approximation K (x, y ) exp[?W (x, y )/?] ? ??1/2

∞ ?∞

(82) ?ypy + W (y) (x, py ) /? dpy

∞ ?∞

K (y) (x, py ) exp

≈ ??1/2 K (y) (xf , 0)e?W (xf ,0)/?

exp ?

a0 4 a1 p + (x ? xf )p2 y + ypy /? dpy . 4 y 2

In terms of ‘stretched’ variables X ≡ (x ? xf )/?1/2 and Y ≡ y/?3/4 this becomes ??1/4 K (y) (xf , 0)e?W (xf ,0)/? e?[XW

′ (x f ,0)/? 1/2 +(X 2 /2)W ′′ (x f ,0)]

P (a1 a0

?1/2

X, a0

?1/4

Y ),

(83)

where the primes denote derivatives with respect to x. Here the canonical di?raction function P (u, v ) ≡

∞ ?∞

exp ?

1 4 1 2 t + ut + vt 4 2

dt

(84)

is a modi?ed (real) Pearcey function (cf. Paris [35]). The expression (83) is an asymptotically (? → 0) valid approximation to the stationary density ρ0 and quasistationary density ρ1 , on the x ? xf = O (?1/2 ), y = O (?3/4 ) lengthscale near the cusp (xf , 0). It supplements the WKB approximation, which is singular there. 36

One sees that on this lengthscale, the pre-exponential factor in ρ0 and ρ1 is actually of magnitude O (??1/4 ). The singularity index of the cusp equals 1/4, as in physical optics. The absence of an ‘i ’ from the exponent gives rise to unusual asymptotic behavior of the di?raction function. The familiar Pearcey fringes of physical optics are replaced by an exponential slope, which becomes increasingly steep as ? → 0. Beyond the cusp (i.e., at x < xf , which is analogous to T < Tc ), the WKB approximation K (x, y ) exp[?W (x, y )/?] is again valid, but W is no longer di?erentiable through the x-axis [28]. This is re?ected in the far?eld asymptotics of the Pearcey function P . One can show that in the far ?eld, i.e., as X = (x ? xf )/?1/2 → ?∞, the expression (83) matches to a WKB approximation displaying this nondi?erentiability. One can also show that the fold caustic emanating from (xf , 0), as in Fig. 2(c), is nonphysical . It arises from subdominant saddle points of the Pearcey integral, and does not contribute to the leading weak-noise asymptotics for ρ0 and ρ1 . This is closely related to the fact that “optimal paths [i.e., physical instanton trajectories] do not encounter caustics,” as Dykman et al. [12] put it. Now the preceding Maslov-WKB treatment is satisfactory so far as it goes. But it leaves unresolved the issue of the validity of the Ginzburg-Landau approximation. The quartic normal form (79) for W (y) = W (y) (x, py ), and the cubic equation of state (80) for its ?rst derivative y = y (x, py ), model a second-order phase transition with mean ?eld (i.e., classical) critical exponents . Equivalently, they model the critical behavior of a system which, though it has a phase transition, has a smooth thermodynamic surface . In the present context, assuming the local validity of the Ginzburg-Landau approximation amounts to assuming that the Lagrangian manifold M is smooth through the point (x, py ) = (xf , 0). Of course the surface turns vertical there, causing ?py /?y to diverge. The assumption is that the singularity can be transformed away by using x and py , rather than x and y , as independent variables. This assumption requires proof. One could presumably justify it by analysing the smoothness (and blowup) properties of solutions of the Hamilton-Jacobi equation. But we shall give a di?erent, more physical justi?cation. First, we shall model the local behavior of W and W (y) by a scaling law, as in the modern theory of critical phenomena. Our treatment will serve as a warmup for Section 8, where we shall analyse the much more complicated (nonclassical) singularity appearing in models where the MPEP is beginning to bifurcate. To see that a scaling law is appropriate in models with a bifurcated MPEP, consider the behavior of the on-axis transverse derivatives w2m = ? 2m W/?x2m (x, 0) as x → x+ f . We know ?1 by (45)–(46) that w2 diverges as (x ? xf ) . The Riccati equation satis?ed by w2 is only the ?rst of a hierarchy of ordinary di?erential equations, describing the evolution of the functions w2m as one moves along the on-axis instanton trajectory from S (where x = xs , and t = ?∞) to the saddle (where x = 0, and t = +∞). For example, w4 satis?es the ODE (25). w2 appears in each of the higher equations, and its blowup will induce a blowup of w4 , w6 , . . . It is not di?cult to show that w2m (x) = ? 2m W (x, 0) ? const × (x ? xf )?(3m?2) , ?x2m 37 x → x+ f, (85)

for 2m = 2, 4, 6, . . . These blowup rates motivate the scaling Ansatz W (x, y ) ? W (x, 0) + (x ? xf )2 h± y , |x ? xf |3/2 x → x± f (86)

for the behavior of W near the cusp (xf , 0). Here the exponents 2 and 3/2 are determined uniquely by the m-dependence of the blowup rates, and the functions h± (·) of the scaling variable z ≡ y/|x ? xf |3/2 are not yet determined (though they must be even). This Ansatz is assumed to be accurate to O ((x ? xf )2 ), when y = O (|x ? xf |3/2 ). We could equally well posit (x ? xf )2 ′′ W (xf , 0) + (x ? xf )2 h± (z ) 2 (87) 2 as x → x± , since we are assuming the accuracy of the scaling Ansatz only up to O (( x ? x f ) ). f The ?rst three terms in this asymptotic approximation are ‘regular’; the scaling behavior appears only in the ?nal, singular term. The exponents 2 and 3/2 are typical of a mean ?eld theory. One can show that the scaling functions h± are also those of a mean ?eld theory. They may be computed by substituting the scaling Ansatz (87) into the Hamilton-Jacobi equation H (x, ?W ) = 0. For this, one needs to rewrite the Hamilton-Jacobi equation in terms of the independent variables x and z . Using the formula (12) for H , and the expansions (19), one ?nds W x, z |x ? xf |3/2 ? W (xf , 0) + (x ? xf )W ′ (xf , 0) + H (x, p) = p2 p2 x + y + ux (x, y )px + uy (x, y )py 2 2 2 px p2 y + + v0 (x) + v2 (x)y 2 + · · · px + u1 (x)y + u3 y 3 + · · · py = 2 2 p2 p2 ′ (xf ) |x ? xf | px ≈ x + y + v0 (xf ) ± v0 2 2

(88)

up to O (|x ? xf |1 ) accuracy, since y = z |x ? xf |3/2 . It follows from the scaling form (86) that up to O (|x ? xf |1 ) accuracy, px (x, y ) = ?W (x, y ) ? ?2v0 (x) ± 2h± (z ) ? (3/2)zh′± (z ) |x ? xf | ?x ′ ≈ ?2 v0 (xf ) ± v0 (xf ) |x ? xf | ± 2h± (z ) ? (3/2)zh′± (z ) |x ? x (89) f| ?W (x, y ) ? |x ? xf |1/2 h′± (z ) ?y (90)

py (x, y ) =

′ where we have used the fact (see Section 5.2) that W ′ (x, 0) = w0 (x) = ?2v0 (x). Substituting 1 (89)–(90) into (88), and setting the coe?cient of |x ? xf | equal to zero, yields the ODE

(h′± )2 = ±v0 (xf ) 4h± ? 3zh′± . 38

(91)

It is easier to solve for z as a function of h′± , than for h± as a function of z . One ?nds z = z (h′± ) = ?C (h′± )3 ? v0 (xf )?1 h′± (92)

where C is undetermined. But z = y/|x?xf |3/2 , and by (90), h′± ? py /|x?xf |1/2 . Rewriting z in terms of y and |x ? xf |, and h′ in terms of py and |x ? xf |, yields

?1 y = y (x, py ) ? ?Cp3 |x ? xf | py y ? v0 (xf )

(93)

=

?Cp3 y

? v0 (xf ) (x ? xf )py .

?1

If one identi?es the model-dependent constant C with a0 , this is precisely eq. (80), the mean ?eld (Ginzburg-Landau) equation of state! It is valid on both sides of the on-axis focus, i.e., both when x ? xf > 0 and when x ? xf < 0. This derivation illustrates how one may go from the pattern of blowup rates of the transverse derivatives w2m (x) = ? 2m W/?y 2m (x, 0) as (xf , 0) is approached, to a scaling form for W , to an equation of state. The singular behavior of the WKB prefactor K can be analysed similarly (we only summarize the analysis). We know by (46) that K (x, 0) diverges as (x ? xf )?1/2 when x → x+ f . A scaling form K (x, y ) ? const × |x ? xf |?1/2 q± y |x ? xf |

3/2

,

x → x± f,

(94)

modelled after the scaling form (86) for W , may be used to approximate K away from the x-axis. This approximation should be accurate to O (|x ? xf |?1/2 ) as x → x± f , when 3/2 y = O (|x ? xf | ). By substituting the two scaling forms (86) and (94) into the transport equation (17) for K , and working to leading order near (xf , 0), one can determine the scaling functions q± = q± (z ). It is easily veri?ed that collecting the O (|x ? xf |?3/2 ) terms in the transport equation yields the ODE

′ 2h′± ± 3v0 (xf )z q± + h′′ ± ± v0 (xf ) q± = 0,

(95)

which q± = q± (z ) must satisfy. Using elementary calculus, and the fact that h± = h± (z ) satis?es the ODE (91), one can show that eq. (95) has solution q± (z ) = const × But since z = y/|x ? xf |3/2 , we know by (90) that h′′ ± (z ) ? |x ? xf | ?py x, y = z |x ? xf |3/2 . ?y (97) ?h′′ ± (z ) . (96)

Substituting (96) and (97) into the scaling form (94) for K reduces it to K (x, y ) ? const × 39 ? ?py (x, y ). ?y (98)

This asymptotic approximation is very simple, and has a profound consequence. We know that the transformed prefactor K (y) (x, py ) can be obtained from K (x, y ) by dividing by a ‘Van Vleck factor,’ as in (68). We therefore have that K (y) (x, py ) ∝ K (x, y ) ? const, ? ?2W (x, y ) ?y 2 (99) (100)

since ? 2 W/?y 2 = ?py /?y . This constant asymptotic approximation is accurate to leading order as x → xf , when y = O (|x ? xf |3/2 ). We have just deduced that on the appropriate lengthscale near the focus, i.e., x ? xf = o(1) and y = O (|x ? xf |3/2 ), the transformed WKB prefactor K (y) does not diverge . K (y) , unlike the prefactor K itself, is asymptotically constant near the focus. This was the crucial assumption made by Dykman et al., and we see that like the Ginzburg-Landau normal form for the Legendre-transformed action, it is justi?ed by our scaling theory of local behavior. We conclude that at least in the case of a generic (cusp) singularity, by investigating the blowup rates of the transverse action derivatives as the singularity is approached, one can derive scaling relations for W and K , and ultimately construct a Maslov-WKB approximation to the stationary probability density near the singularity. This technique is not restricted to singularities of the classical Ginzburg-Landau type.

8

Scaling Behavior Near a Nascent Cusp

Finally, we can construct a scaling theory of weak-noise behavior near the ‘nascent cusp’ singularity appearing at the saddle point of any symmetric double well model, at the onset of bifurcation. The construction will closely parallel the construction of the last section. But several novel features will appear. We shall ?nd that Legendre-transformed versions of the action are approximated, in the vicinity of a nascent cusp, by nonpolynomial normal forms. Equivalently, the nascent cusp singularity, unlike an on-axis focus, will prove to have nonclassical critical exponents. The exponents will depend continuously on the parameter ? ≡ |λy |/λx , which characterizes the linearized drift ?eld at the saddle. The universal presence at criticality of a nongeneric two-sided caustic (which, as shown in Fig. 4, extends sideways from the saddle point) will follow from the nonpolynomial normal forms for the Legendre-transformed actions. Indeed, one of the normal forms will supply a nonpolynomial unfolding of the nongeneric caustic. Moreover, the fact that the critical exponents of the nascent cusp are model-dependent and continuously varying will induce a continuously varying singularity index, and a continuously varying prefactor exponent in the non-Arrhenius weak-noise MFPT asymptotics. To see this, we shall have to go beyond the WKB approximation, by applying the Maslov-WKB method. In Section 8.1 we analyse the scaling properties of the action and the WKB prefactor, and in Section 8.2, we compare our scaling formul? with numerical data. In Section 8.3 we apply the Maslov-WKB method, 40

y x

Figure 8: A sketch of the near-axis Region N, de?ned by |y |/|x|2? ≤ const. The choice of constant is immaterial, so long as it is positive. Most of our expansions for the action W and its Legendre transforms, in double well models at criticality, are valid as (x, y ) → (0, 0) from within Region N. The region could equally well be de?ned by |py |/|x|2? ≤ const, |y |/|px|2? ≤ const, or |py |/|px |2? ≤ const. To leading order near (0, 0), the four de?nitions are equivalent. and construct weak-noise approximations to the stationary and quasistationary probability densities near the saddle.

8.1

Scaling in the WKB Approximation

Our scaling treatment of the nascent cusp begins with an investigation of the blowup rates of the transverse action derivatives w2m (x) = ? 2m W/?y 2m (x, 0), as x → 0+ . Up to now we have written down only the ODE’s satis?ed by w2 (i.e., the Riccati equation (24)) and w4 (i.e., eq. (25)). The full hierarchy of ODE’s may be derived by substituting the Taylor series ∞ 2m /(2m)! for W (x, y ) into the Hamilton-Jacobi equation H (x, ?W ) = 0, and m=0 w2m (x)y separating out the coe?cients of each power of y . One ?nds

′ w ˙ 2m = ?v0 w2 m m

(101)

= ?

k =1

m?1 2m 2m ′ ′ ?2j ]w2 v0 v ?2m , [w2 [w2k /2 + u ?2k?1 ]w2m?2k+2 ? m?2j + 2? j /2 + v 2 j 2k ? 1 j =1

where u ?2j +1 ≡ (2j + 1)! u2j +1 and v ?2j ≡ (2j )! v2j , and u2j +1 and v2j are the drift velocity derivatives de?ned in (19). As usual, the time derivative here is with respect to transit time of the on-axis instanton trajectory, which satis?es x ˙ = ?v0 (x) as it moves from S = (xs , 0) to the saddle. Since v0 (x) ≡ ux (x, 0), this trajectory t → x? (t) moves anti-parallel to the drift. And since u(x, y ) ≈ (λx x, ?|λy |y ) near (0, 0), x? (t) is approximated (as t → +∞) by const × e?λx t . 41

in any double well model at the onset of bifurcation (see eqs. (61)–(62), and Fig. 7). If the fact that w2 → 0 as x → 0+ (i.e., as t → +∞) is substituted into the general ODE (101), it is easy to show, by integrating forward in time toward t = +∞, that w2m (x) ? const × x?(4m?4)? , x → 0+ , (103) for 2m = 4, 6, 8, . . . This pattern of blowup rates, as the saddle is approached, motivates the scaling Ansatz (cf. (86)) W (x, y ) ? W (x, 0) + |x|4? h(y/|x|2?), x → 0. (104)

We showed in Section 5.4, by analysing the Jacobi equation satis?ed by the transverse soft mode, that w2 (x) ? const × x2 , x → 0+ (102)

Here h(·) is some (even) scaling function, as yet undetermined, and the exponents 4? and 2? are determined uniquely by the m-dependence of the blowup rates of (103). This Ansatz is assumed to be accurate to O (|x|4? ) as x → 0, when y = O (|x|2? ). We could equally well posit a ?nite-length asymptotic expansion for W (x, z |x|2? ), namely W (x, z |x|2? ) ? ?

?

?4?? k =0

Here z ≡ y/|x|2? is the scaling variable. Only even powers of x appear in the summation, and by convention, here and below ?4?? denotes the greatest even integer less than or equal to 4?. The expansion (105) is assumed to be accurate to O (|x|4? ), at any ?xed value of z . It can be thought of as an asymptotic development of W (x), as x → ? from within the near-axis region de?ned by the condition |z | ≤ const. This condition de?nes a notch-shaped region, which we call Region N. (See Fig. 8.) It follows by di?erentiating (104) twice with respect to y that w2 (x) ? h′′ (0) as x → 0. For consistency with the ‘splayout’ behavior w2 (x) ? const × x2 of (102), we must have h′′ (0) = 0. Notice the slight discrepancy: the fallo? rate of w2 is not fully captured by the scaling Ansatz. Actually, this is unsurprising. A term proportional to x2 y 2 in W , such as would arise from the O (x2 ) fallo? of w2 as x → 0, would (in terms of x and z = y/|x|2? ) be proportional to |x|4?+2 z 2 . It would therefore be negligible in comparison to the scaling term |x|4? h(z ), as x → 0. The scaling term captures the blowup as x → 0 of w4 , w6 , w8 , . . ., but capturing the precise fallo? rate of w2 would require a more re?ned analysis. We shall not attempt to include in our Ansatz the ‘sub-scaling’ terms that such an analysis would require. The scaling function h(·) may be computed by the technique used in Section 7. By substituting the expression (105) into the Hamilton-Jacobi equation H (x, ?W ) = 0, rewriting the Hamilton-Jacobi equation in terms of the independent variables x and z , and setting the coe?cient of |x|4? equal to zero, one obtains an ODE for h = h(z ). This ODE turns out to be (cf. (91)) (h′ )2 = 2 |λy | [4h ? zh′ ] . (106) 42

xk ? ?k W (0 , 0) + |x|4? h(z ), ?xk k!

?

x → 0.

(105)

As in Section 7, it is easier to solve for z as a function of h′ , than for h as a function of z . One ?nds (cf. (92)) z = z (h′ ) = h′ /2|λy | + c(h′ )1/3 , (107)

where c is undetermined. But z = y/|x|2? and py = ?W/?y ? |x|2? h′ (z ). Rewriting z in terms of x and y , and h′ in terms of x and py , yields an asymptotically accurate equation of state (cf. (93)) /3 (108) y = y (x, py ) ? py /2|λy | + cy;x,py |x|4?/3 p1 y , where cy;x,py = c. In practice the model-dependent constant c would be computed numerically, by ?tting (108) to the ?ow ?eld of instanton trajectories in the vicinity of the saddle. We must have c > 0, since the map py → y is necessarily monotone increasing near the saddle. This is because ‘whorling,’ as in Fig. 6, occurs only in models with a bifurcated MPEP. Whorling is absent at criticality, i.e., at the onset of bifurcation. The equation of state (108) is certainly not of the classical Ginzburg-Landau form. By anti-di?erentiating it, we can obtain an equally unusual approximation to the Legendretransformed action W (y) = ypy ? W , where py = ?W/?y . Since y = ? W (y) /?py , we necessarily have

4?/3 4/3 py W (y) (x, py ) ? ?W (x, 0) + p2 y /4|λy | + Cx,py |x|

where Cx,py = 3c/4. This asymptotic approximation should be accurate to O (|x|4? ), when py = O (|x|2?) [i.e., when y = O (|x|2? ), or when x → ? from within Region N]. The formula (110) can be called a nonpolynomial normal form for the transformed action W (y) near the nascent cusp. Notice that as x → ?, the ?nal, nonpolynomial term 4/3 is signi?cant in a relative sense only within Region N. In the far ?eld of the Cx,py |x|4?/3 py 2? py = O (|x| ) lengthscale, as x → 0 it is increasingly dominated by the p2 y term, and the 4?/3 4/3 py term plays a much more impornormal form reduces to a polynomial. The Cx,py |x| tant role in the near ?eld. One can think of (110) as providing an interpolation between the non-polynomial asymptotic development that is valid as x → ? from within Region N, and the polynomial development that is valid as x → ? from within its far ?eld. The scaling behavior is visible only within Region N. It is worth noting that despite its asymptotic validity, the nonpolynomial normal form (110) does not fully capture the py → 0 behavior of W (y) (x, py ) at ?xed , nonzero x. If the 4/3 nonanalytic py fallo? were exact, it would follow by di?erentiating twice with respect to py that ? 2 W (y) /?p2 y , i.e., ?y/?py , would diverge as py → 0. This would imply that w2 = ?py /?y (y = 0) would be identically zero at any nonzero x. But we know that w2 (x) ? const × x2 , x → 0. The discrepancy is due to the fact that the nonzero w2 near x = 0 arises from ‘sub-scaling’ behavior that we are not attempting to model. It is not di?cult to see that at ?xed nonzero x, the apparent nonanalyticity at py = 0 must be ‘rounded’ at a lengthscale py = O (|x|2?+3 ), or equivalently at y = O (|x|2?+1 ), to yield consistency with the 43

? ??

?

?4?? k =0

xk ? ?k W 4?/3 4/3 py , (0 , 0) + p2 y /4|λy | + Cx,py |x| ?xk k!

?

(109) (110)

w2 (x) ? const × x2 asymptotics. However, on the py = O (|x|2? ) lengthscale the rounding becomes invisible as x → 0. With its continuously varying (in general, irrational) exponent 4?/3, the normal form for W (y) in Region N looks quite di?erent from the normal forms of catastrophe theory [2, 4]. Its most striking feature is the non-analyticity at (x, py ) = (0, 0), which can be interpreted thermodynamically. Recall that W (y) (like W (x) , W , and W ) can be viewed as a thermodynamic potential on the thermodynamic surface (i.e., Lagrangian manifold) M. In fact, through its derivatives it determines the shape of M. So at (x, py ) = (0, 0), or equivalently at (x, y ) = (0, 0), the surface M will itself be non-analytic. However, as ? increases, W (y) becomes increasingly di?erentiable (with respect to x, at least) at x = 0. The order of the ‘phase transition’ appearing at the saddle at criticality is, therefore, an increasing function of ?. We can Legendre-transform the normal form for W (y) to obtain a normal form for the double Legendre transform W = x · p ? W = xpx + W (y) , as a function of px and py . A further Legendre transform will yield a normal form for the remaining thermodynamic potential, W (x) . We sketch only the ?rst of these two computations. Di?erentiating (109)– (110) with respect to x, and using px = ?W/?x = ?? W (y) /?x, yields

/3 px = px (x, py ) ? px (x, 0) + cpx ;x,py |x|4?/3?1 sgn x p4 y

(111)

?

/3 ?2λx x + · · · + const × x?4???1 + cpx ;x,py |x|4?/3?1 sgn x p4 y ,(112)

where cpx ;x,py = ?(4?/3)Cx,py = ??c. Here we have used the fact that W (x, 0) = ?2v0 (x), so that W ′′ (0, 0) = ?2λx , etc. This approximation is accurate to O (|x|4??1 ) as x → 0, when py = O (|x|2? ) [i.e., when x → ? from within Region N]. If ? > 1/2, it is easy to invert the series (112) to approximate x = x(px , py ). The ?2λx x term is dominant, and inversion yields

/3 x = x(px , py ) ? x(px , 0) + cx;px ,py |px |4?/3?1 sgn px p4 y

(113)

?

/3 ?4???1 , ?px /2λx + · · · + const × px + cx;px ,py |px |4?/3?1 sgn px p4 y (114)

where cx;px ,py = ?(2λx )?4?/3 cpx ;x,py = (2λx )?4?/3 ?c. Since x(px , py ) = ? W /?px (px , py ), we must have

4?/3 4/3 py W (px , py ) ? W (px , 0) + p2 y /4|λy | + Cpx ,py |px |

(115) + p2 y / 4 | λy | +

4/3 (116) Cpx ,py |px |4?/3 p y

?

4?/3 4/3 2 py , ≈ ?W (0, 0) ? p2 x /4λx + py /4|λy | + Cpx ,py |px |

?W (0, 0) ?

p2 x / 4 λx

+ · · · + const ×

?4?? px

(117)

where Cpx ,py = (3/4?)cx;px,py = (3c/4)(2λx )?4?/3 . The momentum-space normal forms (115) and (116) should be accurate to O (|px|4? ) as px → 0, when py = O (|px |2? ). This is simply a momentum-space version of the condition that x → ? from within Region N. It is useful to compare the truncated normal form (117) with (75), the quadratic approximation to W that is valid near the saddle point in the absence of focusing. We see that 44

4?/3 4/3 py W (px , py ) ? W (px , 0) + p2 y /4|λy | + Cpx ,py |px | 4?/3 4/3 2 ?4?? py + p2 ? ?W (0, 0) ? px /4λx + · · · + const × px y /4|λy | + Cpx ,py |px |

W (x) (px , y ) ? W (x) (px , 0) ? |λy |y 2 + Cpx ,y |px |4?/3 y 4/3 ?4?? ? ?W (0, 0) ? p2 ? |λy |y 2 + Cpx ,y |px |4?/3 y 4/3 x /4λx + · · · + const × px

4?/3 4/3 py W (y) (x, py ) ? W (y) (x, 0) + p2 y /4|λy | + Cx,py |x| 4?/3 4/3 2 py ? ?W (0, 0) + λx x + · · · + const × x?4?? + p2 y /4|λy | + Cx,py |x|

W (x, y ) ? W (x, 0) + |x|4? h(y/|x|2? ) ? W (0, 0) ? λx x2 + · · · + const × x?4?? + |x|4? h(y/|x|2? ) Figure 9: Normal forms for the thermodynamic potentials (the Legendre transforms of the action W , and W itself) in the vicinity of a nascent cusp. In terms of the model-dependent constant c, Cpx ,py = (3c/4)(2λx )?4?/3 , Cpx ,y = (3c/4)(2λx)?4?/3 (2|λy |)4/3 , and Cx,py = 3c/4. These asymptotic expansions are valid in Region N, in critical double well models with ? > 1/2. They are accurate to O (|x|4? ), or equivalently to O (|px |4? ). the fact that a double well model is ‘critical’ modi?es the double Legendre transform W near the saddle in a very simple way: it adds the ?nal, nonpolynomial term. In a sense, the coe?cient Cpx ,py measures the strength of the nascent cusp singularity at the saddle. The computation of the remaining thermodynamic potential, W (x) , is left to the reader. In Figure 9 we list the normal forms for W , W (x) , and W (y) , as well as the scaling form for W . The expressions listed there are accurate to O (|x|2? ), i.e., to O (|px |2? ), as the nascent cusp is approached from within Region N. In Figure 10 we list the four possible equations of state for x and y . They are accurate to O (|x|2??1 ), i.e., to O (|px |2??1 ), in the same limit. We emphasize that the normal form (116) for W , the normal form for W (x) , and the equations of state that follow from them, are valid only for critical models with ? > 1/2. The reason is that when ? ≤ 1/2, the ?nal term in (112), which when py = O (|x|2? ) is of magnitude O (x4??1 ), is at least as large as the ?2λx x term as x → 0. In fact when ? < 1/2, in Region N (except on the x-axis) the leading asymptotics of px = px (x, py ) are not linear in x. This makes di?cult the computation of asymptotic approximations to x = x(px , py ) and W = W (px , py ). For this reason we shall assume ? > 1/2 henceforth. It is a reasonable conjecture that in critical models where the symmetrical approximation (117) to W = W (px , py ) is valid, it is valid not merely near the px -axis (i.e., in Region N), but uniformly as p → ?. One would like to substitute it into the Maslov-WKB di?raction integral (65), so as to obtain boundary layer approximations to the stationary and quasistationary probability densities near the saddle point (0, 0). The approximation to the quasistationary density would be a replacement for the usual Kramers-type error function approximation, (31). From it, one could derive an Eyring formula for the MFPT 45

x(px , y ) ? x(px , 0) + cx;px ,y |px |4?/3?1 sgn px y 4/3 ?

?4???1 ?px /2λx + . . . + const × px + cx;px,y |px |4?/3?1 sgn px y 4/3

/3 x(px , py ) ? x(px , 0) + cx;px ,py |px |4?/3?1 sgn px p4 y

?

/3 ?4???1 ?px /2λx + . . . + const × px + cx;px,py |px |4?/3?1 sgn px p4 y

1/3 y (x, py ) ? py /2 |λy | + cy;x,py |x|4?/3 py /3 y (px , py ) ? py /2 |λy | + cy;px ,py |px |4?/3 p1 y

Figure 10: The four asymptotic equations of state, which describe the shape of the Lagrangian manifold M in the vicinity of a nascent cusp. It follows by di?erentiating the normal forms listed in Fig. 9 that cx;I,J = (4?/3)CI,J and cy;I,J = (4/3)CI,J . So cx;px ,py = ?c/(2λx )4?/3 , cx;px ,y = ?c(2|λy |)4/3 /(2λx )4?/3 , cy;px ,py = c/(2λx )4?/3 , and cy;x,py = c. These asymptotic expansions are valid in Region N, in critical double well models with ? > 1/2. asymptotics, as in Section 5.2. Unfortunately there is a problem. If Cpx ,py = 0 and (117) becomes quadratic, the Hessian matrix ? 2 W /?pi ?pj is clearly not negative de?nite. As we noted in Section 6, this precludes the use of the two-dimensional di?raction integral (65). The situation does not improve much if Cpx ,py is positive, so it is preferable to use an alternative integral representation. The asymptotic approximation to the Legendre transform W (x) = xpx ? W = ?ypy + W , as a function of px and y , is listed in the table in Figure 9. A truncated version of it would be

4?/3 4/3 2 y , W (x) (px , y ) ≈ ?W (0, 0) ? p2 x /4λx ? |λy | y + Cy,px |px |

(118)

which is a nonpolynomial modi?cation of the Gaussian approximation (76). This approximation is precisely what is needed in the one-dimensional di?raction integral (67), which is what we shall use instead of (65). The reader may wonder about the domain of validity of the approximation (118) to (x) W = W (x) (px , y ). Is it valid outside Region N? In Section 8.2 we present numerical evidence that it is, in fact, a useful asymptotic approximation near the y -axis, even at ?xed, nonzero y . Indeed, it explains the mysterious ‘sideways’ caustic of Fig. 4! To see this, di?erentiate (118) with respect to px to get x = x(px , y ) ≈ ?px /2λx + cx;px ,y y 4/3 |px |(4?/3)?1 sgn px , (119)

which is a truncated version of the asymptotic expansion of x = x(px , y ) listed in Figure 10. If ? < 3/2, the formula (119) predicts that at any nonzero y , the map px → x will not be monotone . This is because the coe?cient cx;px ,y is positive. By examination, if

?1 < |x| ? const × |y |(3/2??) ,

y → 0,

(120)

46

then the inverse map px = px (x, y ) (and hence W = W (x, y )) will be multivalued. This inequality de?nes a two-sided nongeneric caustic , which emanates from (0, 0) along the positive and negative y -axes. In the language of catastrophe theory, the formula (119) is a (non-smooth) unfolding of the nongeneric caustic emanating from the nascent cusp. It also resembles a thermodynamic equation of state in the vicinity of a phase transition. However, the thermodynamic interpretation of the variables di?ers from the case of an on-axis focus, as analysed in the last section. Here |y | is analogous to Tc ? T , for example. And by examination, the thermodynamic critical exponent γ is nonzero whenever ? < 3/2; it equals (3/2 ? ?)?1 . In any event, the critical exponents of the nascent cusp are clearly nonclassical: they depend continuously on the parameter ?. We caution the reader that in arbitrary double well models at criticality, the nonpolynomial approximation (118) and the nonclassical equation of state (119) may not necessarily describe the px → 0 behavior of x = x(px , y ) at ?xed, nonzero y . If ? < 3/4, the y 4/3 |px |(4?/3)?1 sgn px term in (119) would cause x = x(px , y ) to diverge as px → 0, at any nonzero y or py . Such a divergence would greatly distort the shape of the Lagrangian manifold M. So we shall assume ? > 3/4 henceforth. In Section 8.2 we present numerical evidence of the need for the ? > 3/4 restriction, and also verify that the nongeneric caustic is present if and only if ? < 3/2. Incidentally, our numerical results indicate that at ?xed nonzero y , the apparent non-analyticity at px = 0 is ‘rounded’ at a su?ciently small (|y |-dependent) lengthscale, as px → 0. This is analogous to the abovementioned rounding of W (y) (x, py ), and its derivative y = y (x, py ), as py → 0 at ?xed nonzero x. There is also a problem with the nonpolynomial approximation (118) and the nonclassical equation of state (119) when ? ≥ 3. To see this, note that a more complete asymptotic expansion of x = x(px , y ) near px = 0 would presumably be of the form

?4???1 x = x(px , y ) ≈ ?px /2λx + . . . + const × px + cx;px ,y y 4/3 |px |4?/3?1 sgn px .

(121)

Such an asymptotic expansion is listed in Fig. 10, and is certainly valid as x → ? from within Region N. If a similar expansion is valid near the y -axis, we see that there will be a crossover at ? = 3 between two regimes. When ? = 3, the nonpolynomial cx;px ,y y 4/3 |px |4?/3?1 sgn px 3 term in (121) becomes cx,px ,y y 4/3 p3 x . This is increasingly dominated by the px term in (121) as y → 0. In fact when ? is raised above 3, at small |y | the leading corrections to the naive x ? ?px /2λx behavior are no longer given by the nonpolynomial term, but rather by the p3 x term. For this reason we shall assume for the remainder of our analysis that ? < 3 as well as ? > 3/4. To use the one-dimensional Maslov-WKB di?raction integral (67) as promised, we need to approximate in the vicinity of the nascent cusp at (px , y ) = (0, 0) not only the Legendretransformed action W (x) (px , y ), but also the transformed prefactor K (x) (px , y ). It may be approximated in a very similar way, which we only summarize. By (61)–(62), K (x, 0) diverges in any critical model as |x|?? , when x → 0. A scaling form K (x, y ) ? |x|?? q (y/|x|2?), 47 (122)

modelled after (94) and (104), may be used to approximate K away from the x-axis. The approximation (122) should be accurate to O (|x|??) as x → 0, when y = O (|x|2?) [i.e., as x → ? from within Region N]. By substituting (122) and the scaling form (104) for W into the amplitude transport equation (17), and working to leading order near (x, y ) = (0, 0), one can determine the scaling function q = q (z ). The procedure closely resembles the procedure used in Section 7. It is easily veri?ed that collecting the O (|x|?? ) terms in the transport equation yields the ODE 2 h′ + |λy | z q ′ + h′′ q = 0, (123) which q = q (z ) must satisfy. (Cf. (95).) Here h = h(z ) is the scaling function for W , which satis?es the ODE (106). Using elementary calculus, one can show that eq. (123) has solution q (z ) = const × |h′ (z )|

?1/3

?h′′ (z ).

(124)

(Cf. (96).) But since py = ?W/?y ? |x|2? h′ , one may write |x|?2? py for h′ , and ?py /?y for h′′ . Substituting (124) into the scaling form (122), and performing the indicated rewriting, yields K (x, y ) ? const × |x|??/3 |py |?1/3 ? ?py (x, y ). ?y (125)

(Cf. (98).) This asymptotic approximation is exact to leading order as x → ? from within Region N. The formula (125) facilitates the computation of the transformed prefactor K (y) = (y ) K (x, py ). It may be computed from K as in (99), by dividing by the appropriate ‘Van Vleck factor.’ We immediately ?nd K (y) (x, py ) ∝ K (x, y ) ? ?2W (x, y ) ?y 2 (126) (127)

? const × |x|??/3 |py |?1/3 ,

since ? 2 W/?y 2 = ?py /?y . (Cf. (99)–(100).) The uncomplicated asymptotic approximation (127) should be accurate to O (|x|??) as x → 0, when py = O (|x|2? ). It simply says that K (y) (x, a|x|2? ) ? const × |a|?1/3 |x|?? as x → 0, for any nonzero a. The two remaining transformed prefactors, K and K (x) , may be computed from K (y) by dividing (or multiplying) by the appropriate Van Vleck factors. (Cf. (66) and (68).) For example, K (px , py ) ∝ K (y) (x, py ) ? The details are left to the reader. One ?nds K (px , py ) ? const × |px |??/3 |py |?1/3 ,

(x)

?2W (p x , p y ). ?p2 x

(128)

(129) (130)

K

(px , y ) ? const × |px | 48

??/3

|y |

?1/3

.

The transformed prefactor K (x) is the one we need for the Maslov-WKB di?raction integral. In Section 8.2, we examine the numerical evidence for the validity of this asymptotic approximation to K (x) = K (x) (px , y ), when px → 0 at ?xed nonzero y . Remarkably, the formula (130) predicts that K (x) (px , y ) diverges at the location (px , y ) = (0, 0) of the nascent cusp. This is di?erent from the case of a generic (cusp) singularity, treated in Section 7. It is also di?erent from the geometrical optics limit of physical optics, where the transformed amplitude function near a singularity is normally a ‘slowly varying’ (i.e., non-singular) function [4]. The fact that the transformed prefactor K (x) diverges at the nascent cusp is at least as important to the weak-noise behavior of critical double well models as the fact that the normal form for the transformed action W (x) is nonpolynomial.

8.2

Comparison with Numerics

We now summarize the numerical evidence for the validity, in double well models at criticality, of our nonpolynomial normal form for the Legendre-transformed action W (x) = W (x) (px , y ), and our approximation to the transformed WKB prefactor K (x) = K (x) (px , y ). We shall see that both are valid approximations near the y -axis separatrix; in particular, near the saddle point. This justi?es their use in the Maslov-WKB method, which we shall employ in Section 8.3 to construct boundary layer approximations to the stationary and quasistationary probability distributions of double well models at criticality. We begin by examining the evidence for the nonpolynomial normal form (118) for W (x) = (x) W (px , y ). Actually we shall study the related nonpolynomial approximation (119) to its ?rst derivative x = x(px , y ), i.e., x = x(px , y ) ≈ ?px /2λx + cx;px ,y y 4/3 |px |(4?/3)?1 sgn px . (131)

As explained above, we expect on theoretical grounds that this approximation is generically valid near the saddle point, in critical models in which the quotient ? ≡ |λy | /λx satis?es 3/4 < ? < 3. The formula (131) predicts that at nonzero y , the correspondence px → x is monotone if ? ≥ 3/2, but non-monotone at nonzero y if ? < 3/2. When ? < 3/2, the correspondence px → x is analogous to the correspondence m → ?h, in a ferromagnet, between magnetization and (negative) magnetic ?eld. It is easily checked that when

?1 < |x| ? const × |y |(3/2??) ,

y → 0,

(132)

the inverse map x → px is three-valued rather than single-valued. In this region the three possible values for px are by examination of the same magnitude as x, i.e., px = O |y |(3/2??)

?1

,

y → 0.

(133)

We interpret the inequality (132) as de?ning a two-sided nongeneric caustic centered on the y -axis, in the interior of which the action W , and its gradient p = ?W , are three-valued. 49

1

1

0.5

0.5

0

0

-0.5

-0.5

-1 -1.5

-1

-0.5

0

0.5

1

1.5

-1 -1.5

-1

-0.5

0

0.5

1

1.5

1

1

0.5

0.5

0

0

-0.5

-0.5

-1 -1.5

-1

-0.5

0

0.5

1

1.5

-1 -1.5

-1

-0.5

0

0.5

1

1.5

Figure 11: The ?ow ?eld of instanton trajectories emanating from both stable ?xed points S and S ′ , in critical versions of the standard double well model (6). Parts (a), (b), (c), (d) correspond to models with ? = 0.725, 0.85, 1.15, and 1.6. In all cases the parameter α is set equal to the critical value αc = 2?(? + 1), at which the MPEP bifurcates. The two-sided nongeneric caustic of Fig. 4 is visible in parts (b) and (c), but it has separated into two generic caustics in part (d). 50

Recall that the Lagrangian manifold M, which is traced out by WKB bicharacteristics, comprises all points in phase space of the form (x, p(x)). The three-valuedness of W and p within the caustic accordingly implies that there are three points on M ‘above’ any point x in the interior of the caustic. We have already seen in Fig. 4 that a nongeneric caustic qualitatively agreeing with this prediction does indeed appear in the ? = 1 standard double well model (6), at criticality. Figures 11(a)–11(d) show the ?ow ?eld of instanton trajectories, i.e., projected bicharacteristics, in several more critical variants of the standard double well model. (At criticality α = αc = 2?(? + 1), by eq. (58).) Figures 11(b) and 11(c), with ? = 0.85 and ? = 1.15, illustrate the fact that the two-sided caustic of Fig. 4 appears at criticality in any double well model whose parameter ? ≡ |λy |/λx satis?es 3/4 < ? < 3/2. The caustic disappears, as expected, in critical models with ? ≥ 3/2. Figure 11(d) shows what happens. As ? in the standard model is raised above 3/2 (with α set equal to αc = αc (?)), the two-sided caustic separates into two one-sided generic caustics, whose cusps move out along the positive and negative y -axes, away from the saddle. In any critical model with ? > 3/2, there is a portion of the separatrix near the saddle that is not crossed by any instanton trajectory. Figure 11(a) illustrates the bizarre behavior that occurs in critical models with ? ≤ 3/4. At ?rst glance it seems that the now familiar two-sided caustic is present, but closer study reveals that points in its interior are reached by only two instanton trajectories, rather than three. Apparently, in the ? ≤ 3/4 regime the approximation (131) breaks down near the separatrix. Empirically, when ? ≤ 3/4 the |px |(4?/3)?1 factor in (131) must be replaced by unity. The ? ≤ 3/4 regime is still under investigation, and we shall not consider it further in this paper. We can now compare the predictions of our scaling theory with numerical data. Figure 12(a) is a section through the caustic of Fig. 11(b), i.e., a cross-section through the corresponding Lagrangian manifold. It shows the correspondence px → x, at y = 0.05, in the ? = 0.85 standard model at criticality. The qualitative shape of the curve certainly resembles the prediction of formula (131). But before making a quantitative comparison, we need to discuss the interpretation of (131). It was derived from an asymptotic development of the action about the saddle point. To what extent does it describe the small-|px | asymptotics of x = x(px , y ) at ?xed , nonzero y ? That is what is plotted in Fig. 12. From a rigorous point of view, when ? < 3/2 the formula (131) provides a two-term ?1 asymptotic expansion of x = x(px , y ) as y → 0, on the px = O (|y |(3/2??) ) lengthscale on which the nongeneric caustic is visible. This is strongly reminiscent of the ‘Region N’ constraint of the last section. There we began by approximating W = W (x, y ) in the notch-shaped region of Fig. 8. Here we are approximating x = x(px , y ), and by extension its anti-derivative W (x) = W (x) (px , y ), in a region that is similarly notch-shaped, but is centered on the y -axis rather than the x-axis. We shall not attempt to expand x and W (x) systematically, but the basic procedure is plain. If we de?ne a new scaling variable z ? ≡ px / |y |(3/2??)

?1

,

(134)

then formula (131) can be interpreted as comprising the ?rst two terms in an asymptotic 51

x 0.03 0.02 0.01 0 -0.01 -0.02 -0.03 -0.08-0.06-0.04-0.02 0 0.02 0.04 0.06 0.08 px

0.01

sca

lin

g

0.001 1e-06 1e-05 0.0001 0.001

sub-sca

ling

0.01

0.1

Figure 12: A section through the nongeneric caustic appearing at criticality in the ? = 0.85 standard double well model, as shown in Fig. 11(b). This section is taken at y = 0.05. Parts (a) and (b) show a linear and a logarithmic plot of x = x(px , y = 0.05). The scaling behavior (x ∝ |px |(4?/3)?1 ) and the sub-scaling behavior (x ∝ px ) are both visible. Crossover between ?1 . the two regimes occurs at px = O (|y |(5/2)(3/2??) ), i.e., at px = 10?5 . development of x(? z |y |(3/2??) , y ) as y → 0. The development should be valid at any ?xed z ?. Though this restatement is a bit pedantic, it suggests that on smaller lengthscales than ?1 px = O (|y |(3/2??) ), the formula (131) might be invalid. Actually there are strong reasons for believing that the nonanalytic |px |(4?/3)?1 sgn px behavior as px → 0 does not appear at ?xed, nonzero y . If it did, ?x/?px (px = 0) would diverge at criticality, at y = 0, in any model with 3/4 < ? < 3/2. Equivalently, ?px /?x(x = 0), i.e., ? 2 W/?x2 (x = 0), would be identically zero , irrespective of the choice of nonzero y . But this prediction is too simple: it ignores the presence of ‘sub-scaling’ terms. We noted in the last section that at criticality, W = W (x, y ) should contain an x2 y 2 term, for consistency with the sub-scaling w2 (x) ? const × x2 behavior. In other words, ? 2 W/?x2 (x = 0) near the saddle should be nonzero and proportional to y 2 . For consistency with this prediction, the nonanalytic |px |(4?/3)?1 sgn px behavior of x = x(px , y ) must be ‘rounded’ at su?ciently small |px |. It is easy to check ?1 that px = O (|y |(5/2)(3/2??) ) is the correct lengthscale. On that lengthscale, one should ?nd W (x, y ) ≈ const × x2 y 2 , i.e., px ≈ const × xy 2 , or x ≈ const × y ?2 px . What we conclude from this discussion is that at ?xed, nonzero y , the nonpolynomial ?1 formula (131) for x = x(px , y ) should be valid on the ‘caustic lengthscale’ px = O (|y |(3/2??) ), ?1 but that it will break down when px is decreased to O (|y |(5/2)(3/2??) ). On that smaller lengthscale, one expects a crossover to a linear regime, where x is proportional to px . We can now proceed to our comparison with numerics. In Fig. 12(b) we plot the correspondence px → x(px , y = 0.05) of Fig. 12(a) on a logarithmic scale. We also ?t two trendlines to it: 52

?1

x ∝ |px |(4?/3)?1 and x ∝ px . As the two trendlines reveal, our theoretical analysis is perfectly con?rmed. There is indeed a crossover to a linear, sub-scaling regime when px is decreased ?1 ?1 to O (|y |(5/2)(3/2??) ). But at larger lengthscales, e.g., px = O (|y |(3/2??) ), the fractional power |px |(4?/3)?1 sgn px of the scaling formula (131) is clearly visible. Similar crossover plots can be obtained for other critical double well models, whose parameter ? = |λy |/λx lies in the range 3/4 < ? < 3/2. Our asymptotic approximation K (x) (px , y ) ? const × |px |??/3 |y |?1/3 to the transformed prefactor K (x) at criticality, derived in Section 8.1, can also be numerically tested. The approximation should be valid in the vicinity of the y -axis separatrix, i.e., as px → 0 at ?xed, nonzero y . There are two separate cases: 3/4 < ? < 3/2, when a caustic is present, and ? ≥ 3/2, when one is not. For simplicity we consider only the latter. When ? ≥ 3/2, it follows from (131) that x(px , y ) ≈ ?px /2λx as px → 0; the nonpolynomial term is subdominant. So our asymptotic approximation to K (x) implies that (cf. (68)) K (x, y ) ∝ K (x) (px , y ) ? ?2W (p x , y ) ?x2 ?px = K (x) (px , y ) ? (p x , y ) ?x ≈ const × |x|??/3 |y |?1/3 , (135) (136) (137)

as x → 0 at ?xed, nonzero y . This comparatively slow power-law divergence as the separatrix is approached at (small) nonzero y is to be contrasted with the K (x, 0) ? const × |x|?? divergence that occurs when the saddle point is approached along the x-axis. (See (62).) It is susceptible to numerical test. In Fig. 13 we graph K = K (x, y ) as a function of x, at y = 0.05, for the critical version of the standard double well model with ? = 1.6. (This is the same model whose instanton trajectories are shown in Fig. 11(d).) The curve is ?tted to high accuracy by a power-law const × x?0.533 , i.e., const × x??/3 . This con?rms the prediction of our scaling theory. No subscaling regime is evident at small |x|. Similar plots can be obtained for the near-separatrix behavior of K in critical models with other values of ?. We conclude that the asymptotic approximations to W (x) = W (x) (px , y ) and K (x) = (x) K (px , y ) derived from our scaling theory have a wide domain of validity, and may be employed in the Maslov-WKB method.

8.3

Scaling Beyond the WKB Approximation

We can now compute the Maslov-WKB boundary layer approximations to the stationary density ρ0 and the quasistationary density ρ1 near the nascent cusp at the saddle point, where the conventional WKB approximation breaks down. The boundary layer approximations are determined by (118) and (130), the asymptotic approximations to the Legendre-transformed action W (x) and the transformed prefactor K (x) respectively. As discussed above, these approximations should be valid when the eigenvalue ratio ? ≡ |λy |/λx satis?es 3/4 < ? < 3. 53

K 1000

100

10 1e-05

0.0001

0.001

0.01

0.1

x

Figure 13: Behavior of the WKB prefactor near the y -axis separatrix, in a critical version of the ? = 1.6 standard double well model. K = K (x, y = 0.05) is plotted on a logarithmic scale, revealing the scaling behavior (K ∝ x??/3 at nonzero y ). Substituting (118) and (130) into the one-dimensional di?raction integral (67) yields the rather complicated expression K (x, y ) exp[?W (x, y )/?] ? ??1/2 (138) ?xpx + W (x) (px , y ) /? dpx

2

K (x) (px , y ) exp

≈ const × ??1/2 e?W (0,0)/? |y |?1/3 e?|λy |y /? p2 x × |px |??/3 exp ? ? Cpx ,y y 4/3 |px |4?/3 + xpx 4 λx

?

dpx

which requires a bit of explanation. The ?rst problem to be resolved is the lengthscale near (x, y ) = (0, 0) on which this di?raction integral de?nes a valid Maslov-WKB approxima2 tion. The p2 x and xpx terms in the argument of exp(·) are O (1) when px and xpx are O (?); i.e., when x and px are O (?1/2 ). The term Cpx ,y y 4/3 |px |4?/3 /? is O (1) when, also, y = O (?3/4??/2 ). This is the case when ? < 3/2, at least. If ? ≥ 3/2 then the Cpx ,y y 4/3 |px |4?/3 /? term is negligible whenever y = o(1). We conclude that in the weak-noise (? → 0) limit of models with ? < 3/2, the di?raction integral (138) de?nes a valid Maslov-WKB approximation on the x = O (?1/2 ), y = O (?3/4??/2 ) lengthscale near the saddle point. This is precisely the caus?1 ?1 tic lengthscale x = O (|y |(3/2??) ), i.e., px = O (|y |(3/2??) ), of the last section. If ? ≥ 3/2, so that no caustic is present, then y = O (?3/4??/2 ) must be replaced by y = o(1). On the appropriate lengthscale, the di?raction integral de?nes a noncanonical di?raction function. 54

The di?raction integral, being one-dimensional, cannot resolve the singularity at y = 0. This is because it does not include an integration over py . So one cannot expect the MaslovWKB approximation to be valid at arbitrarily small y . The stationary and quasistationary densities in critical models are expected to be tightly concentrated on the y = O (?1/2 ) transverse lengthscale near the saddle point, as we saw in (29)–(31) (which apply in the absence of focusing). At most, the Maslov-WKB approximation will be valid in the far ?eld 2 of the y = O (?1/2 ) lengthscale, where the factor e?|λy |y /? is exponentially small. So it will not be directly comparable to (29)–(31). But it proves to be very useful nonetheless. De?ne ‘stretched’ variables X ≡ x/?1/2 and Y ≡ y/?1/2 ; also, change the integration variable to Px = px /?1/2 . The approximation becomes const × ??(?+1)/6 e?W (0,0)/? |Y |?1/3 e?|λy |Y

2

|Px |??/3 e?XPx e?Px /4λx dPx ,

2

(139)

since when px , y = O (?1/2 ), the term Cpx ,y y 4/3 |px |4?/3 /? is negligible. At this point we must explain how to interpret the integration over Px , or px . The stationary density ρ0 (x, y ) = ρ0 (X?1/2 , Y ?1/2 ) must be even in X , and the quasistationary density ρ1 (x, y ) = ρ1 (X?1/2 , Y ?1/2 ) must be odd. To get approximations with these symmetry properties, we may integrate Px from ?∞ to ∞ and 0 to ∞ respectively. We remarked at the end of Section 6 that performing a half-range integration is one way of incorporating an antisymmetry constraint, and that is the technique we shall use. As summarized in Abramowitz and Stegun (Ref. [1], §19.5), de?nite integrals resembling (139) de?ne parabolic cylinder functions . Evaluating the integral (139), with the two possible choices of the range of integration (full range and half range), yields the MaslovWKB approximations

/2 1/2 +λx x ρ0 (x, y ) ? const × ??(?+1)/6 F0 (λ1 )e x x/? 1/2 ρ1 (x, y ) ? const × ??(?+1)/6 F1 (λx x/?1/2 )e+λx x

2 /?

y/?1/2 y/?1/2

?1/3 ?1/3

e?|λy |y e?|λy |y

2 /?

,

(140) (141)

2 /?

2 /?

to the stationary and quasistationary probability densities, on the O (?1/2 ) lengthscale near the saddle point. Here the so-called boundary layer functions Fi = Fi (Z ), where Z ≡ /2 1/2 1/2 λ1 , are de?ned by x X = λx x/? Fi (Z ) ≡ yi+1 (1/2 ? ?/3, 21/2 Z )e?Z

2 /2

,

(142)

in the notation of Abramowitz and Stegun. y1 (1/2 ? ?/3, ?) and y2 (1/2 ? ?/3, ?) are even and odd parabolic cylinder functions, respectively. We could equally well de?ne the boundary layer functions Fi in terms of an Hermite function of non-integer index , by F0 (Z ) ≡ H(?/3)?1 (Z )e?Z

2

even

,

F1 (Z ) ≡ H(?/3)?1 (Z )e?Z

2

odd

.

(143)

Here [·]even and [·]odd signify even and odd parts, under the re?ection Z → ?Z . The de?nitions (143) are meaningful whenever the index n ≡ (?/3) ? 1 is not an integer, in which 55

case the Hermite function Hn (Z ) is not a conventional Hermite polynomial, and is neither even nor odd. But since we are assuming 3/4 < ? < 3, this is always the case. Irrespective of the choice of de?nitions, F0 (Z ) ? const × |Z |??/3 , F1 (Z ) ? const × |Z |??/3 sgn Z,

(144)

as Z → ±∞. (Cf. Ref. [1], §19.8.) We stress that the Maslov-WKB approximations to ρ0 and ρ1 are strictly valid only in the transverse far ?eld , i.e., as Y = y/?1/2 → ±∞. But they make very clear how critical double well models di?er from non-critical double well models. By comparing (140)–(141) with (29)–(31), we see that at criticality, the boundary layer functions F0 (·) and F1 (·) replace the boundary layer functions 1 and erf(·) respectively. The approximations (140)–(141) are guaranteed to match to the standard WKB approximation K (x) exp[?W (x)], as one moves in a transverse direction away from the saddle point. For example, the |Z |??/3 fallo? of eqs. (144) will match to the |x|??/3 prefactor fallo? of eqs. (125) and (127), which is seen, e.g., in Figure 13. The Maslov-WKB approximations to ρ0 and ρ1 , and the nonpolynomial normal form for the Legendre-transformed action that engendered them, have several striking consequences for double well models at criticality. ? A nongeneric caustic, emerging sideways from the nascent cusp at the saddle. In Section 8.1 we predicted from the nonpolynomial normal form for W that when 3/4 < ?1 < ? < 3/2, a caustic is located at |x| ? const × |y |(3/2??) . Our prediction was con?rmed by Figure 11. This caustic has an unusual (continuously varying) exponent. It is nongeneric , in the sense of singularity theory. ? An unusual (continuously varying) singularity index. As ? → 0, the fallo? of the stationary density ρ0 at the saddle point (0, 0) is not pure exponential, on account of the ??(?+1)/6 prefactor in the Maslov-WKB approximation (140). This is interpreted as a statement that the nascent cusp has singularity index s = s(?) = (? + 1)/6, as mentioned in Section 4. It too is nongeneric, in the sense of singularity theory. ? Non-Arrhenius MFPT asymptotics. If one computes the rate at which the quasistationary density ρ1 is absorbed on the separatrix near the saddle, the ??(?+1)/6 prefactor in the Maslov-WKB approximation (141) will appear in the ? → 0 asymptotics. Equivalently, the exponentially decaying eigenvalue λ1 = λ1 (?) of the Smoluchowski operator will have an asymptotic ??(?+1)/6 prefactor, as well as the usual Arrhenius factor [i.e., exp(??W/?)]. And the MFPT will be asymptotic to const × ?+(?+1)/6 exp(+?W/?), as ? → 0. At criticality, the weak-noise growth of the MFPT is slower than pure exponential . ? A non-Gaussian limiting exit location distribution. In the absence of MPEP bifurcation, for a symmetric double well model the location of the point of exit from either 56

of the two wells would have an asymptotic Gaussian distribution, on the transverse 2 O (?1/2 ) lengthscale near the saddle. In fact, its density would fall o? as e?|λy |y /? . We see from the Maslov-WKB approximation to ρ1 that at criticality, the exit location density on the separatrix includes scaling corrections . In the transverse far ?eld it falls 2 2 o? as |y |?1/3e?|λy |y /? , rather than e?|λy |y /? . It has a non-Gaussian tail . These phenomena look natural from the point of view of the theory of critical phenomena, though the stochastic escape problem has not previously been considered from that point of view.

9

Discussion

We can now step back and review our results. We began with a WKB treatment of the weak-noise asymptotics of stationary (and quasistationary) solutions of the Smoluchowski equation. The WKB analysis led to instanton trajectories, which have a physical interpretation as most probable weak-noise ?uctuational paths. The instanton trajectories turned out to be zero-energy trajectories of an associated Hamiltonian dynamical system. This is because the phase space versions of the instanton trajectories (i.e., WKB bicharacteristics) trace out a Lagrangian manifold in phase space. In double well models the onset of bifurcation is associated with the ?eeting appearance of an unusual singularity (a nascent cusp) in the shape of this manifold, as the parameters of the model are varied. There is a formal analogy between the Lagrangian manifold of a dynamical system perturbed by weak noise, and the thermodynamic surface of a condensed matter system. This analogy led us to construct a scaling theory of the shape of the Lagrangian manifold near the nascent cusp. To date, most work on Lagrangian manifolds has assumed that they are smooth , and that any apparent singularities in their shape can be transformed away by a change of coordinates. This is analogous to assuming that thermodynamic surfaces are real analytic, and that non-analyticities in thermodynamic behavior (i.e., phase transitions) can be transformed away by working in terms of the appropriate thermodynamic potential. Equivalently, it is analogous to assuming that all phase transitions have classical critical exponents. Our scaling theory makes it clear that the nascent cusp singularity is a genuine point of non-smoothness of the Lagrangian manifold. In thermodynamic terms, it has nonclassical, indeed continuously varying, critical exponents. Applying the Maslov-WKB method to the nascent cusp yielded several interesting predictions, which we summed up in the four bulleted items at the end of the last section. One normally expects that in a double well system perturbed by weak noise of strength ?, the rate of inter-well hopping λ1 = λ1 (?) will be asymptotic to a constant multiple of the Arrhenius factor exp(??W/?), where ?W is an e?ective barrier height. Also, one expects that the distribution of exit locations (from either well) will asymptotically become a Gaussian of O (?1/2 ) standard deviation, centered on the saddle point between the two wells. The Maslov-WKB method predicts that at criticality, both these phenomena are strongly altered. In particular, the factor exp(??W/?) must be replaced by ??s exp(??W/?), where s = (? + 1)/6 is the 57

log λ1

1/ε

Figure 14: A sketch, on a logarithmic scale, of the rate of noise-activated inter-well hopping λ1 as a function of the reciprocal noise strength 1/?. O? criticality (solid curve), λ1 displays a pure exponential fallo? as ? → 0, i.e., λ1 ? const × exp(??W/?). At the onset of bifurcation (dashed curve), λ1 ? const × ??s exp(??W/?) instead, where s is the singularity index of the nascent cusp appearing at the saddle point. singularity index of the nascent cusp. (As we noted in Section 5.3, the singularity index is a sort of critical exponent.) In Fig. 14 we sketch an Arrhenius plot, showing this anomalous (non-Arrhenius) behavior. In mathematical terms, the nascent cusp appearing at the onset of bifurcation is a nongeneric singularity, i.e., a singularity di?erent from any of the now classical singularities of catastrophe theory. As shown in Fig. 4, in many double well models it induces an unusual caustic in the ?ow ?eld of instanton trajectories. This caustic is itself nongeneric, in that its exponent is not equal to 3/2. As we have seen (see, e.g., Fig. 12), its presence quantitatively con?rms the validity of our scaling theory. It is remarkable that such nongeneric phenomena are a generic feature of singly parametrized symmetric double well models. At least as developed in this paper, our scaling theory is a scaling theory of weaknoise behavior near the nascent cusp, precisely at criticality. It would be useful to treat as well models that are nearly critical, but not exactly so. Such models should display a crossover from non-Arrhenius behavior to Arrhenius behavior at su?ciently weak noise strength. By developing a joint scaling theory , one of the variables in which measures the distance from criticality, it should be possible to analyse this phenomenon. We expect that it is possible to derive a ‘Ginzburg criterion’ [27], expressing how close to criticality any given double well model should be, for the non-Arrhenius behavior of Fig. 14 to be visible. Work 58

on this is under way. We brie?y mention two geometric features of models ‘o? criticality’ that cry out for a theoretical explanation. A nongeneric caustic appears near the separatrix not only at criticality, but in many non-critical models as well [12]. Also, as criticality is approached (e.g., as ? α → αc in the standard model of (6)), it frequently happens that the nascent cusp is formed by a collision of two generic cusps, which move along the separatrix toward the saddle point. These phenomena can presumably be explained by an appropriate joint unfolding , but that is for the future. We close by mentioning a possible extension of a more theoretical sort. In this paper we have focused exclusively on the asymptotic solutions of the time-independent weak-noise Smoluchowski equation. There is reason to believe that nongeneric singularities resembling the nascent cusp can occur, and are perhaps even widespread, in the asymptotic solutions of other singularly perturbed elliptic partial di?erential equations. Most WKB treatments of singularly perturbed elliptic PDE’s (see, e.g., Duistermaat [10]) assume that each WKB characteristic (i.e., instanton trajectory) eventually leaves any bounded region of space. This assumption is violated in the Smoluchowski equation for any double well model, since the MPEP(s) terminate on the saddle point, rather than extending to in?nity. We expect that when it is violated in other PDE’s, analogous nongeneric singularities in formal asymptotic solutions can occur. The nongeneric singular phenomena that we have seen in this paper may simply be representatives of a larger class.

References

[1] M. Abramowitz and I. A. Stegun, editors. Handbook of Mathematical Functions . Dover, New York, 1965. [2] V. I. Arnol’d. Critical Points of Smooth Functions and Their Normal Forms. Russian Math. Surveys 30:5 (1975), 1–75. [3] A. Auerbach and S. Kivelson. The Path Decomposition and Multidimensional Tunneling. Nuclear Phys. B 257 [FS 14] (1985), 799–858. [4] M. V. Berry. Waves and Thom’s Theorem. Adv. in Phys. 25 (1976), 1–26. [5] M. V. Berry and N. L. Balazs. Evolution of Semiclassical Quantum States in Phase Space. J. Phys. A 12 (1979), 625–642. [6] A. J. Bray and A. J. McKane. Instanton Calculation of the Escape Rate for Activation over a Potential Barrier Driven by Colored Noise. Phys. Rev. Lett. 62 (1989), 493–496. [7] B. Carmeli, V. Mujica, and A. Nitzan. Dynamics of Multidimensional Barrier Crossing in the Overdamped Limit. Berichte der Bunsen-Gesellschaft 95 (1991), 319–326. [8] B. Caroli, C. Caroli, B. Roulet, and J.-F. Gouyet. A WKB Treatment of Di?usion in a Multidimensional Bistable Potential. J. Statist. Phys. 22 (1980), 515–536.

59

[9] M. V. Day. Recent Progress on the Small Parameter Exit Problem. Stochastics 20 (1987), 121–150. [10] J. J. Duistermaat. Oscillatory Integrals, Lagrange Immersions, and Unfolding of Singularities. Comm. Pure Appl. Math. 27 (1974), 207–281. [11] M. I. Dykman, P. V. E. McClintock, V. N. Smelyanskiy, N. D. Stein, and N. G. Stocks. Optimal Paths and the Prehistory Problem for Large Fluctuations in Noise-Driven Systems. Phys. Rev. Lett. 68 (1992), 2718–2721. [12] M. I. Dykman, M. M. Millonas, and V. N. Smelyanskiy. Observable and Hidden Singular Features of Large Fluctuations in Nonequilibrium Systems. Phys. Lett. A 195 (1994), 53–58. Available as cond-mat/9410056. [13] J.-P. Eckmann and R. S? en? eor. The Maslov-WKB Method for the (an-)Harmonic Oscillator. Arch. Rational Mech. Anal. 61 (1976), 153–173. [14] M. I. Freidlin and A. D. Wentzell. Random Perturbations of Dynamical Systems . SpringerVerlag, New York/Berlin, 1984. [15] S. Glasstone, K. J. Laidler, and H. Eyring. The Theory of Rate Processes . McGraw-Hill, New York, 1941. [16] R. Graham. Macroscopic Potentials, Bifurcations and Noise in Dissipative Systems. In Theory of Continuous Fokker-Planck Systems , edited by F. Moss and P. V. E. McClintock, volume 1 of Noise in Nonlinear Dynamical Systems , chapter 7, pp. 225–278. Cambridge University Press, Cambridge, 1989. [17] M. C. Gutzwiller. Periodic Orbits and Classical Quantization Conditions. J. Math. Phys. 12 (1971), 343–358. [18] M. C. Gutzwiller. Chaos in Classical and Quantum Mechanics , volume 1 of Interdisciplinary Applied Mathematics . Springer-Verlag, New York/Berlin, 1990. [19] H. R. Jauslin. Nondi?erentiable Potentials for Nonequilibrium Steady States. Physica A 144 (1987), 179–191. [20] J. B. Keller and S. I. Rubinow. Asymptotic Solution of Eigenvalue Problems. Ann. Physics 9 (1960), 24–75. [21] M. M. Klosek-Dygas, B. M. Ho?man, B. J. Matkowsky, A. Nitzan, M. A. Ratner, and Z. Schuss. Di?usion Theory of Multidimensional Activated Rate Processes: The Role of Anisotropy. J. Chem. Phys. 90 (1989), 1141–1148. [22] H. A. Kramers. Brownian Motion on a Field of Force and the Di?usion Model of Chemical Reactions. Physica 7 (1940), 284–304. [23] R. Landauer. Motion out of Noisy States. J. Statist. Phys. 53 (1988), 233–248. [24] J. S. Langer. Statistical Theory of the Decay of Metastable States. Ann. Physics 54 (1969), 258–275.

60

[25] R. G. Littlejohn. The Van Vleck Formula, Maslov Theory, and Phase Space Geometry. J. Statist. Phys. 68 (1992), 7–50. [26] D. Ludwig. Persistence of Dynamical Systems Under Random Perturbations. SIAM Rev. 17 (1975), 605–640. [27] S.-K. Ma. Modern Theory of Critical Phenomena . W. A. Benjamin, Reading, Massachusetts, 1976. [28] R. S. Maier and D. L. Stein. The E?ect of Focusing and Caustics on Exit Phenomena in Systems Lacking Detailed Balance. Phys. Rev. Lett. 71 (1993), 1783–1786. Available as chaodyn/9305010. [29] R. S. Maier and D. L. Stein. The Escape Problem for Irreversible Systems. Phys. Rev. E 48 (1993), 931–938. Available as chao-dyn/9303017. [30] R. S. Maier and D. L. Stein. Asymptotic Exit Location Distributions in the Stochastic Exit Problem. Submitted to SIAM J. Appl. Math., 1994. Available as adap-org/9407003. [31] M. Mangel. Barrier Transitions Driven by Fluctuations, with Applications to Evolution and Ecology. Theoret. Population Biol. 45 (1994), 16–40. [32] V. P. Maslov and M. V. Fedoriuk. Semi-Classical Approximation in Quantum Mechanics. Reidel, Boston/Dordrecht, 1981. [33] T. Naeh, M. M. Klosek, B. J. Matkowsky, and Z. Schuss. A Direct Approach to the Exit Problem. SIAM J. Appl. Math. 50 (1990), 595–627. [34] L. Onsager and S. Machlup. Fluctuations and Irreversible Processes. Phys. Rev. 91 (1953), 1505–1512. [35] R. B. Paris. The Asymptotic Behaviour of Pearcey’s Integral for Complex Variables. Proc. Roy. Soc. London Ser. A 432 (1991), 391–426. [36] L. S. Schulman. Techniques and Applications of Path Integration . Wiley, New York, 1981. [37] M. A. Shayman. A Geometric View of the Matrix Riccati Equation. In The Riccati Equation , edited by S. Bittanti, A. J. Laub, and J. C. Willems, pp. 89–112. Springer-Verlag, New York/Berlin, 1991. [38] P. Talkner. Mean First Passage Times and the Lifetime of a Metastable State. Z. Phys. B 68 (1987), 201–207.

61

赞助商链接

- Gaussian scaling for the critical spread-out contact process above the upper critical dimen
- Developments for Reference--State One--Particle Density--Matrix Theory
- Excited states from time-dependent density functional theory
- Measuring the elements of the optical density matrix
- Kinetic Theory of Plasmas Translational Energy

更多相关文章：
更多相关标签：