Nonadiabatic escape and stochastic resonance

We analyze the fluctuation-driven escape of particles from a metastable state under the influence of a weak periodic force. We develop an asymptotic method to solve the appropriate Fokker–Planck equation with mixed natural and absorbing boundary conditions. The approach uses two boundary layers flanking an interior region; most of the probability is concentrated within the boundary layer near the metastable point of the potential and particles transit the interior region before exiting the domain through the other boundary layer, which is near the unstable maximal point of the potential. The dominant processes in each region are given by approximate time-dependent solutions matched to construct the approximate composite solution, which gives the rate of escape with weak periodic forcing. Using reflection we extend the method to a double well potential influenced by white noise and weak periodic forcing, and thereby derive a two-state stochastic model—the simplest treatment of stochastic resonance theory—in the nonadiabatic limit.


Introduction
The escape of particles from a metastable state under the influence of noise is a classical problem in non-equilibrium statistical mechanics [1]. Calculating the rate of escape can be approached in a variety of ways [see e.g. [2][3][4][5]and Refs. therein]. An important extension of the original escape problem includes the influence of periodic forcing, with a phase that impacts the escape rate [6,7]. Of particular relevance here is the problem of stochastic resonance, wherein the combined effect of background noise and weak periodic forcing control the state of the system [8][9][10]. Indeed, because of the compelling consequences of such resonances, there are many methods that have been developed to calculate the escape rate, ranging from eigenfunction expansions [11] to path-integrals [12,13]. However, the simplest solution used to study the principle characteristics of stochastic resonance appeals to the approx imation of the adiabatic limit [14].
Recently, concepts of stochastic resonance have been utilized in numerous fields including sensory biology [e.g. 15,16], image processing [e.g. 17,18], signal detection and processing [e.g. 19,20], and energy harvesting [e.g. 21,22]. The broad impact of stochastic resonance is often viewed as counter-intuitive because rather than background noise obscuring the detection of a weak signal, it leads instead to an enhancement of that signal.
The canonical configuration of stochastic resonance focuses on particles in a double-well potential influenced by white noise and weak periodic forcing. Although there are multiple time-scales involved in the dynamics, the principal interest concerns the time it takes for a particle to transition from one stable point in the potential to the other. Hence, it is common to consider a two-state model using a master equation that describes the time-evolution of the probability density of two discrete states and their exchange rates [14]. Thus one uses the classical escape rate from one metastable point; when the rate is independent of the slowly varying phase of the periodic forcing, this is called the adiabatic limit. However, considering the extent of the fields in which stochastic resonance plays a role, from climate to engineering to biology [23], it is of interest to go beyond the adiabatic limit. Here we address the nonadiabatic situation, in the two-state framework, to determine the escape rate when the phase of the periodic forcing does not vary slowly. To achieve this we introduce an asymptotic method to obtain an explicit expression for the escape rate, and in so doing we show how to transform the original double-well potential problem into the two-state model.

Escape rate under periodic forcing
First, we consider the escape rate from the metastable region of a potential U(x) under the influence of weak periodic forcing Acos(ωt) and noise induced fluctuations, as shown in figure 1. The Fokker-Planck equation for the probability density in this situation is with boundary conditions P(x = −∞,t) = P(x = x max + χ,t) = 0, wherein the tilde's denote dimensional variables. We assume that the magnitude of the noise, σ, and the amplitude of the periodic forcing, A, are both small (in the precise sense outlined below), and that The potential U has the characteristic 'diffusivity' scale, ∆U = U max − U min , and length scale, ∆x = x max − x min , which leads to the three parameters Using the dimensionless variables with boundary conditions P(x = −∞, t) = P(x = χ/∆x, t) = 0, and the local minimum and maximum are at x = −1 and x = 0 respectively. The underlying scaling assumptions are that 1 and r = O(1). The small magnitude of is associated with the kinetic energy of a particle near the minimum of the potential being much less than the potential energy, ∆U, necessary to escape from it. The assumption that the external forcing, Acos(ωt), is weak relative to the thermal noise is embodied by r ∼ O( 2 ). The dimensionless frequency, Ω, is the ratio of the time-scale for a particle to reach a quasistationary state near the potential minimum, (∆x) 2 /∆U , to the oscillation time-scale ω −1 of the potential. The assumption that Ω 1 implies that the period of the external forcing is much longer than the time required for a particle to reach a quasi-stationary state near the minimum. This assumption facilitates the asymptotic matching procedure near the maximum.  Near x = −1, the potential U(x) can be approximated as U −1 + 1 2 a(x + 1) 2 , where a ≡ |U min |. Thus, we rewrite equation (3) in terms of a state variable P B1 (η, t) that depends on the stretched coordinate η = (x + 1)/ : The leading-order solution, written in terms of the original position variable, is where n 1 , which will be determined as part of the matching procedure, is a slowly-varying function of time satisfying 1 n1 dn1 dt 1. Implicit in this solution is therefore that the probability density reaches a quasi-steady state around the potential minimum, which arises because the weak noise in system drives only a small leakage of probability across the barrier at the maximum.

Potential maximum-boundary Layer 2 (BL2): |x| ∼ O( )
Near the maximum ( x = x max ) we let x = ζ and the potential U can be approximated as (3) is rewritten in terms of the state vari- Equation (6) must be solved subject to P B2 (ζ = ∞, t) = 0 and that the solution match to that for the interior region between the extrema of the potential, as outlined presently. The match, however, implies that the solution for the interior delivers a probability flux to the boundary layer around the maximum that varies periodically in time with frequency Ω. This precludes a straightforward solution of equation (6) if Ω = O(1). Instead, we avoid solving the boundary-layer problem for BL2 for general Ω, and adopt the convenient approximation that the oscillation frequency is small, Ω 1. This allows us to neglect the left-hand side of equation (6) and write the quasi-stationary approximation,

Interior region
Within the interior region between the two extrema in the potential, the rapid exponential decline of the probability P suggests that we adopt a WKBJ-type ansatz, viz., from which the leading-order Fokker-Planck equation is The characteristic curves of equation (9) are given by where Note that, in the limit that Ω 1, an integration by parts furnishes the simpler approximation, This approximation exposes an issue with general solution in equation (11): the transit time function T(x) in equation (12) diverges logarithmically for x → x min or x → x max . This does not present a problem at the potential minimum in view of the integration limits in either equation (11) or (13), but it does obscure the limit to the maximum. In fact, for |x| = O(δ) with 1 δ , the interior solution should actually be matched to the boundarylayer solution, leaving T ∼ −b −1 log δ + O(1). Thus we write equation (11) as where T max ∼ t 0 − b −1 log δ and t 0 is an (undetermined) order-one constant time shift that should, in principle, be fixed by a matching argument.

Asymptotic matching & uniform approximation
Asymptotic matching near the local minimum leads to and near the local maximum it is required that The preceding results suggest an approximation that is valid throughout the two boundary layers and the interior region: In the small frequency (Ω 1) approximation, equation (17) becomes Near x = ±1, this approximation reduces to the two boundary layer solutions in equations (5) and (7), with D 0 given by the relevant approximation of equation (16), whereas in the interior it reduces to the solution implied by equation (13). Note that the limit of the last integral in equation (18) introduces the approximation δ ≈ , thereby furnishing a solution that depends on only a single small parameter and avoids any exercise in matching.

Exit rate
Now, we construct the exit rate by a suitable integration of the Fokker-Planck equation (3) as follows. Note that and because the probability is principally concentrated near the minimum, The flux at the origin is given by Hence, by combining equations (16) and (17) we obtain thereby giving the escape rate R ≡ 1 n1 dn1 dt as In the small frequency (Ω 1) approximation, by substituting equation (18) into

Adiabatic limit
In the adiabatic limit, we discard the terms containing factors of Ω in equation (18), to arrive at with the associated escape rate

The cubic potential
For the cubic potential we have Hence we obtain where so that if we again set δ = − , for which Υ(Ω) = then we can capture the suppression of the periodic adiabatic variation of the escape rate by non-adiabatic effects.
In figure 2(a) we compare a numerical solution of the Fokker-Planck equation, (3), with the non-adiabatic (equations (17) and (18)) and adiabatic (equation (25)) analytical solutions. Our numerical method for equation (3) is based on the implicit finite difference scheme introduced by Chang and Cooper [24]. Both of the non-adiabatic analytical solutions match the numerical solution at the percentage level of accuracy, save for the transition region from the interior to the boundary layer near the maximum (x = 0). However, the adiabatic solution differs substantially from both of the others, as is particularly evident when Ωt = π/2 where the non-adiabatic contribution, sin(Ωt), is maximal.
In figure 2(b) we show the time evolution of the escape rates for the four solutions, defined as J x=0 / 0 −∞ Pdx for the numerical solution. To calculate the analytical solutions, we use the approximation n 1 = 1, which is accurate to machine precision for the parameter settings and times used in the figure. (Likewise, for the double-well potential below, we use the approximation n 1 = n 2 = 0.5 in comparing the asymptotic predictions with numerics in figure 6.) In both amplitude and phase, the non-adiabatic solution in (31) compares more favourably with the numerical solution than the adiabatic result in (24).
Finally, in figure 3 we show the deviations of the various approximations as a function of frequency. Note the substantial differences in phase and amplitude between figures 3(a) and (c). In particular, while the non-adiabatic analytic solution (31) compares well with the numerical solution, there is a pronounced deviation of the adiabatic solution in both phase and amplitude of the maximum escape rates.

Double-well potential and stochastic resonance
We now treat Brownian particles in a double-well potential under the influence of weak periodic forcing, which is the original configuration of stochastic resonance [8,9]. By reflection of figure 1 we extend the approach described above to construct the approximate solutions in the five regions shown in figure 4. The potential Ũ is scaled as before, so that U = (Ũ − U max )/∆U, where ∆U is now a measure of the height of the barrier, and we define ∆x as half the distance between the two minima. As the potential may not be symmetrical, this translates to a scaled potential that vanishes at x = 0 and takes the values U 1 and U 2 at the two minima x 1 and x 2 , respectively. We replace the absorbing boundary near the local maximum with the usual boundary condition, P(±∞, t) = 0, insuring that the probability is conserved throughout the entire domain as particles move between the two minima.
The asymptotic solution is 2} ), and the 'constants' of integration, n 1 , n 2 , S 1 , S 2 , D 0 and D 1 must be connected by matching the five solutions together. In particular, we find

cos[Ω(t + T(x) − T(− ))]dx and
if we again make the approximation that the match can be accomplished at x = ± . More compactly we have (35) In the small frequency (Ω 1) approximation, equation (35) becomes (36) Global conservation of probability implies that ∞ −∞ P n 1 + n 2 = 1, where dn1 dt = J| x=xmax , which leads to or equivalently where (40) Equation (38) has the same form as the two-state Master equation used in other approaches to stochastic resonance theory [e.g. 14]. The factors Υ {1,2} (Ω) again represent the suppression of the adiabatic variation of the escape rates by non-adiabatic effects; Θ {1,2} (Ω) represent additional phase shifts. The functions Υ {1,2} (Ω) are equal for a symmetrical base potential, and are illustrated in figure 5 for the quartic potential with U = 1 4 x 4 − 1 2 x 2 . In the theory of stochastic resonance presented in [14], the adiabatic variation of the escape rates (given by the replacements Υ {1,2} (Ω) → 1 and Θ {1,2} (Ω) → 0 in (39)) is responsible for the characteristic improvement in the signal-to-noise ratio. Thus, the suppression factors Υ {1,2} (Ω) determine the destructive effects of non-adiabaticity in a generalization of that analysis. In view of our small-frequency approximation, the exponential decline of Υ 1 = Υ 2 for large Ω seen in figure 5 is again only approximate.
In figure 6 we compare the numerical and analytic solutions from equations (35) and (36), for the double well potential U = 1 4 x 4 − 1 2 x 2 . Because our approach is a reflection of the   escape calculation described in the previous section, we expect the match to be good, save for the same transition region, here between the two interior regions and the central maximum.

Conclusion
We have developed an asymptotic method of calculating the probability density function and the associated escape rate of Brownian particles from a metastable state under weak periodic forcing. The approach uses boundary layers near the two extremes, where the potential U(x) is approximately quadratic and the time-dependent linear Fokker-Planck equations can be solved. In the interior layer separating these, an advection-dominated solution is constructed and the three approximate solutions are matched. Because the evolution of the total probability is equal to the probability flux at the absorbing boundary, we can integrate the Fokker-Planck equation over the complete domain and determine the escape rate in the non-adiabatic limit. Finally, by reflection we extended this asymptotic approach to the problem of a double-well potential with weak periodic forcing to find a solution to the problem of stochastic resonance in the non-adiabatic case. In particular, the ease with which equation (39) can be used, and its limits understood through figure 5, provide substantial applicability. Given the ubiquity of stochastic resonance, this result is likely of the broadest relevance. Finally, the approach we take here is complimentary to other general approaches, which focus on the universality of fast-slow systems in stochastic resonance and two state systems [25][26][27].