Diffusion in a class of exactly solvable non-harmonic potentials. Intrinsic effects induced by non-linearities

Physica D 41 (1990) 79-88 North-Holland DIFFUSION IN A CLASS OF EXACTLY SOLVABLE NON-HARMONIC POTENTIALS. INTRINSIC EFFECTS INDUCED BY NON-LINEARITIES J.L. ROMERO and J. RAMIREZ Departamento de Matem~ticas, Unioersidad de Cddiz, P.O. Box 40, Puerto Real (Cadiz), Spain F. ROMERO and J.F.R. ARCHILLA Departamento de Fisica Terrica, Facultad de Fisica, P.O. Box 1065, Seoilla, Spain Received 1 March 1988 Revised manuscript received 11 December 1988 Accepted 17 June 1989 Communicated by A.V. Holden This paper deals with the problem of a particle that diffuses in a potential with a reflecting barrier and has a point of stable equilibrium and a point of unstable equilibrium. Based on the exact solutions obtained earlier for the Fokker-Planck equation of a class of these models, we analyze the behavior of the probability density, the mean path and the onset time which determines the transition from unimodal to bimodal probability densities. The study is made over different initial positions, two of them very close to the unstable point, which permits a clear comparison among the subsequent evolutions, and the observation of some intrinsic effects induced by non-linearities. 1. Introduction The problem of diffusion in one-dimensional potentials, and particularly in non-harmonic potentials, has been the subject of many recent studies. One of the motives for the importar~ce of this problem lies in its connection with Langevin's description of non-equilibrium phase transitions. When the white noise constitutes a good idealization of the fluctuations concerning a macrovariable, the probability density of this macrovariable obeys a Fokker-Planck equation (FPE), which must be solved under appropriate boundary conditions. In the theory of reaction kinetics, many studies use a model of a particle moving in a one-dimensional potential under the effect of thermal noise and damping. Following Kramers [1], many use a model of a potential with a single well and one barrier. Initially the particle is placed in the potential well and can only escape by passing over the potential barrier. The objective is to calculate the rate of escape out of the well, and the dependence of this rate on temperature, friction, and the parameters of the potential. In the highfriction regime one is led to analyze a FPE of the Smoluchowski type [2, 3]. The model of Kramers has played a central role in many areas, other than chemical reactions in condensed phases [4], such as surface desorption [5] and surface catalysis [6]. Since the original analysis by Kramers a number of authors have considered parametric forms of the potential in an attempt both to describe the evolution of different processes and to check the validity of Kramers' results. One of the first investigations is due to van 0167-2789/90/$03.50 © Elsevier Science Publishers B.V. (North-Holland) 80 J.L. Romero et aL / Diffusion in a class of exactly solvable non-harmonic potentials Kampen [7], who showed that Kramers' calculation leads to a correct value for the case of a solvable double-well symmetric potential. Owing to its importance there has been considerable research both in the development of approximation schemes [8-10] and in obtaining exactly solvable models [4, 11] which cover a rich variety of behaviours. Although the mathematical expressions for the potentials obtained for these exactly solvable models are relatively complicated, their qualitative forms are very close to some interesting potentials for which no exact solutions exist. We think that the study of these models is important for two reasons. First, to learn about intrinsic effects induced by non-linearities, such as time-controlled transitions from unimodal to bimodal probability densities, branching times and boomerang behaviour of the mean path. Second, exact solutions can serve to test the validity of the different approximation methods. In ref. [12] we showed the existence of a large class of non-linear stochastic processes, with constant diffusion, which admits exact solutions. This class is defined by a two-parameter family of non-harmonic potentials, and includes potentials with either reflecting or absorbing barriers. This paper concerns with a particular model of that family: the model of the potential has a reflecting barrier at the origin of coordinates, a stable equilibrium point and an unstable equilibrium point. It is analogous to the potential used by Kramers [1] in the context of reaction kinetics, although in our case we can only treat the high-friction regime, as we only know the exact solutions corresponding to the FPE of the Smoluchowski type. This paper is organized as follows. First, we expose a summary of some previous results, showing a particular class of models which has exact solutions, its explicit solutions, and the three qualitatively different kinds of potentials that are included. Second, we analyze a model of the type mentioned above by adequately selecting two values for the parameters that define the family of potentials. Having centered on a specific model, we expose details on the evolution of the probability density, the behaviour of the mean path and the splitting of the probability density, for different types of initial conditions. We conclude with some comments about the observed phenomena and some further studies that are underway. 2. Smmnary of previous results We consider the stochastic differential equation (Langevin equation) for a driven variable x(t): dx _ D(x) + 71(t) = dU(x) dt d------x~ + T/(t), (1) where U(x) is a potential field from which the drift D(x) is derived, and 7/(t) is a Gaussian white noise, characterized by (,/(t)) =0, (2) (~/(t) ~i(t')) = 28(t - t'), (3) with the initial condition being x(0) =x o. (4) Under these conditions, the corresponding FPE for the distribution function P(x, t) is OP(x,t)~t =-~[[~]a [[~U(X) Ip¢ ,t)] a2p(x,t) + ~X 2 with the properties P(x,t)>_O, Vxen, teR+=[0, oo), (Sa) frP(x,t)dx=l, Vt~a +, (5b) P(x,o) -- a(x - go). (5c) For some cases it may be required that x ~ R +. J. L Romero et al./ Diffusion in a class of exactly solvable non-harmonic potentials 81 There exist exact solutions for this FPE for the special choice of the potential U(x): U(x) = -21n V(x), (6) with V(x) = (½x 2)"+ 1/4 exp (- ¼x 2) X ,Fx(A, B, ½x2), (7) where ,FI(A, B, z) denotes Kummer's hypergeometric function of the first kind [13], and A=½+a-~fl, B=l+Za, (8) with lit,) o] A=O B:6 cl I : .18/5 I~ : 4.5 Fig. 1. Shapes of U(x) with a = 0.25 and (a) fl = 6, (b) /8 = 5.99, (c) fl = 4.5. a>-}, A>0, (9) with When a > - ¼, the potential may adopt some physically interesting behaviours. The three qualitatively different forms of U(x) when a > - ¼ are represented in fig. 1. All the curves tend asymptotically to + oo when x tends to zero, i.e. the point x = 0 represents a reflecting barrier. When A = 0 the potentials tend to + m for x tending to + o0, i.e. the potentials are attractive Vx ~ [i +, like that represented in fig. la. When A q. 0, the potentials tend to - c~ for x tending to + oo. For a given value of a >- ¼, there always exists a value A(a) > 0, such that if A ~ (0, A(a)) the potential adopts a form like that represented in fig. lb, i.e. with one stable point and one unstable point. The value A(a) may be found by solving a transcendental equation• In other cases, with A 4: 0, the potentials have no depth, i.e. they are like fig. lc. For the particular case a = ¼, U(x) adopts the form O(x) = - 2 ln[(x/¢~)exp(- ¼x 2) XxFl( A,a2, ½x2)], (10) a=l-i#. With the prescription that x ~ R +, the solution for P(x, t), with P(x,O) = 8(x - Xo), is P(x, t) = Pl(X, t; Xo) + Pl(X, t; - go) , x ~ R +, (11) where Px( x, t; Xo) = (4~) - '/2(sinh t ) -'/2(X/Xo) iFI(A,3v_; ½x 2) • 1 2 X ,Fx(A, 3, _~Xo) X exp(¼flt)exp(-¼(X-Xoe-t)2e') sinh t (12) 3. Analysis of the model We shall now show that our class of models covers a great variety of behaviours induced by non-linearity. Nevertheless, it is not possible to obtain an exact explicit calculation of all the rele- 82 J.L. Romero et aL / Diffusion in a class of exactly solvable non-harmonic potentials vant physical quantities. The principal reason is that there is no possibility of choosing an initial value x 0 that makes (11) a symmetrical function. Thus, for example, for the calculation of the positions and the nature of the extrema of P(x, t), which are determined by means of the equation P'(x, t)---0, we are led to the resolution of a transcendental equation, which cannot generally be solved analytically. To overcome this difficulty partially we have made exhaustive computer determinations, using (11) and (12), of the most relevant quantities and of the time evolution of P(x, t) for several different values of the initial position x 0. In what follows we will expose details of some results concerning these computations. In our study we have fixed the parameters a and fl to the values a = 0.25 and fl --- 5.99. These give a potential U(x) like that of fig. lb, with one stable point x s and one unstable point x~. The computation of x s and of x u gives the values x s = 1.416 and x~ = 4.537. These parameters correspond to a potential with enough depth that one can expect a rich variety of phenomena to appear easily. For this case U(xu) - U(xs) --- 5.39. 4. Time evolution of the probability density For the study of the behaviour of the probability density P(x, t), and of some relevant physical quantities, we may select different initial positions x0, which correspond to qualitatively different initial physical situations. We shall sketch and comment on the time evolution of the probability density P(x, t) for cases which clearly show phenomena induced by nonlinearity. P(xA) 0.0 1.6 3.~ / x,= 4.50 4.B 6.4 8.0 9.6 II.2 12.8 1( Fig. 2. Probability density P(x, t); x o = 4.50, t = 0.1 ~ (0.1) 0.6. evolution, for t = 0.1 --* (0.1) ~ 0.6. At t = 0.1 the curve is practically Gaussian, the abscissa of the maximum Xu is such that x M < xu. Nevertheless by t = 0.2 one may appreciate the existence of two peaks for the probability density. There exists a critical time t a, defined as the time for this splitring of the probability density such that t a _< 0.2. We observe that the ordinate of the first peak is greater than the ordinate of the second one, and that the abscissa of this second peak is greater than xu. The formation of this second peak in the region (x~, oo) may be viewed as a "tunneling process" through the potential barrier drawn in fig. lb. Initially the abscissa of the first peak moves towards x s, and the abscissa of the second peak moves towards oo. 4.1. Case x 0 = 4.50 As the unstable point corresponds to xu= 4.5379 , the initial situation is near this unstable point. One may expect phenomena such as the formation of two maxima for the probability density, i.e. the occurrence of a transition from unimodal to bimodal states. Fig. 2 shows us this 4.2. Case x o = 4.60 Now x 0 >__ x u. There is also a transition to bimodal states after a time greater than 0.1, as we can see in fig. 3. The "tunneling process" leads to the formation of a maximum in the region [x s, xu], i.e. it is produced in backward direction. If we J.L. Romero et aL / Diffusion in a class of exactly solvable non-harmonic potentials 83 Ptx,L) 6o 0 2 V 6 8 lO 12 14 16 x Fig. 3. Probability density P(x, t); x o = 4.60, t = 0.1 --* (0.1) --, 1.1. compare the forms of P(x, t) for t--0.3 corresponding to the cases x 0 = 4.5 and x 0 = 4.6, we can see that for the second case the curve is more symmetrical, which is a consequence of the asymmetrical form of the potential in the region around x,. Later on we shall comment on this asymmetry. For the times considered in fig. 3, i.e. t = 0.1 (0.1) ---, 1.1, there exists a monotonic movement of the abscissa of the first peak towards x,. At time t = 1.1 an increase in the ordinate of the first peak with respect to the ordinate for t = 0.9 is observed. The existence of this increase is much better appreciated in fig. 4, where for t = 60 the ordinate is greater than 0.2. At successive instants of time the ordinates decrease monotonically and abscissas are stabilized near x s. The existence of a reflecting wall at x = 0 prevents the probability density from passing through it. Reflection is responsible of the transient growing of the first peak and of the existence of a second "tunneling process" which causes a further decrease in the ordinates. P(x,LI 6 8- ,e~_ _ 0.0 x, : LGO I.fi 3.2 4.8 6.4 8.0 x Fig. 4. Probability density P(x, t); x o = 4.60, t = 60 --, (60) --, 110. 4. 3. Case Xo = 5 This is the last case we consider. There exists also a transition to bimodal states as a consequence of a "tunneling process" through the potential barrier. Now the form of the curves of P(x,t) after this transition are much more asymmetrical (see fig. 5). 5. Behaviour of the mean path Let us analyze the behaviour of the mean path (x(t)> corresponding to the process (1) when the probability density is given by (11) and (12). This magnitude is determined by calculating the integral (x(t)> = fo°°Xe(x, t) dx. (13) The analytical operations required to calculate this integral are quite tedious, so we have preferred to determine this integral via computer calculation only for some special initial values x 0. No phenomena of interest appear for the cases 84 J.L. Romero et al./ Diffusion in a class of exactly solvable non-harmonic potentials PIx,[I Xt= 5 0 4 8 12 16 20 24 X Fig. 5. Probability density P(x, t); x 0 = 5, t = 0.1 ~ (0.1) --* 0.6. x 0 < Xs, or x 0 > x u. We have selected two cases corresponding to x s < x 0 < x u. In fig. 6 we represent (x(t)) versus t for x 0 = 1.5. We observe that initially the values of (x(t)) decrease until a time t* ---0.1, and then the mean path increases monotonically, i.e. the mean path is stabilized for t < t* and suffers a destabilization for t > t*. The minimum value of (x(t)) does not attain x s, (X>n~n ~ 1.491. This situation is called boomerang behaviour of the mean path [14]. Fig. 7 corresponds to the case x o = 3. We observe also the existence of an initial period of stabilization followed by a destabilization at t* 1.5. As before the minimum value of (x(t)> does not attain xs, for this case (x)~n --- 1.81. For the case x0=4.50 we detect also a boomerang behaviour, but in this case t* is of the order of 0.06, i.e. the initial period of stabilization is negligible. The numerical studies that we have made demonstrate the existence of an interval of initial values x o for which the process presents a period <x>~ I0-' ,d-, 0.0 0.2 0.4 0.6 0.8 1.0 1.2 1.4 t Fig. 6. Boomerang behaviour of the mean path for x o = 1.5. of stabilization followed by a destabilization of the mean path. Although not demonstrated mathematically, it could be inferred that the condition for the existence of the so-called "boomerang effect" is that x s < x 0 < x u. 6. Transition to bimodal states. Times for bifurcation As we have pointed out in section 4, for a given value x 0' the probability density P(x, t) presents two peaks after a critical time t m called the time for bifurcation or for transition to bimodal states, which depends obviously on the initial position x 0. Bifurcation such that the two peaks have significant values, occurs only for a limited interval of values x o around x u. We can consider the space (x o, t) as having the role of a parameter space C, and (x) the role of a behaviour space X. In this way the probability density is represented by a smooth map: P: C x X ~ R. (14) J.L. Romero et aL / Diffusion in a class of exactly solvable non-harmonic potentials 85 o 8 ,-J. ...... , , , , , , 0.0 0.4 0.8 1.2 1.6 2.0 t Fig. 7. Boomerang behaviour for the mean path for x 0 = 3. , , , , , , , , , , , , , , , , 1.2 2.0 2.8 3.6 4.4 5.2 6.0 6.8 7.6 x, The time-dependent positions of the extrema of P(x, t) are determined by the equation 3P(x, t) _ 0. (15) P'(x, t) - Ox The set M = ((x 0, t, x)/P'(x, t) = 0} constitutes a manifold, of dimension two, that can be considered as a catastrophe manifold. Each point of M corresponds to a minimum or a maximum of P(x, t). In plates I and II we represent two perspective drawings of the manifold M. Axes are represented in blue, red and black, and correspond to x o, t, x, respectively. Plate I corresponds to a two-dimensional projection defined by the vector ( -0.5, 5,4) of the space (x 0, t, x), and plate II is a projection defined by the vector (0.6, -0.6,5). In this way, it is possible to give the same impression that one would obtain from looking through the actual three-dimensional surface, along the direction of projection. In both figures, the color black, of the projected surfaces, corresponds to the second maxima of P(x, t), red corresponds to the minima of P(x, t), and blue corresponds well to the first maxima, or well to ~=, ,=; 8 Fig. 8. Graphical representation of the behaviour of the time for bifurcation, t~, for different values of x 0. the unique maxima for times prior to bifurcation. This manifold is like that corresponding to the so-called "dual cusp catastrophe" of catastrophe theory. In both pictures the fold line is relatively well distinguished, this being the smooth curve: Mr3 ((Xo, t,x)/P"(x,t ) =0}. (16) Projection of the fold line, onto the plane (x0, t), produces a curve, represented in fig. 8, that constitutes the bifurcation set, i.e. it gives us for each initial value x0, the corresponding time for bifurcation t B. This bifurcation set has a singularity for X0m ~ 4.577, tBm = 0.137, and we observe that this value of Xom is slightly greater than the abscissa of the unstable point x u - 4.537. Our potential U(x) has an asymmetrical form for a neigbourhood of x u, i.e. the unstable region has not a symmetrical form. This can be seen by studying the behaviour of the second derivate of U(x) with respect to x, U"(x), in that unstable region. The curve, represented in fig. 9, has a minimum at x m --- 4.620, U"(Xm) -~ -7.481, and 86 J.L. Romero et al./ Diffusion in a class of exactly solvable non-harmonic potentials Plate I. Two-dimensional projection of the catastrophe manifold M, defined by the vector ( - 0.5, 5, 4). these values permit us to determinate the coefficients of the following approximate form of U'(x) for a neigbourhood of Xu: U'( x ) ~ a( x - xu) 3 + b( x - Xu) 2 +c(x-Xu), (17) resulting in a -~ 7.756, b -~ - 1.931, c ~- - 7.295. The ~act that b #: 0 indicates a non-symmetrical form of U(x) in the unstable region. The minimum time for bifurcation tBm corresponds to an initial position x 0 = X0m, such that x u < XOm < Xnr We observe a relation of tBm with U"(Xom ) given by 1 /Bm ~ IU,,(XOm)[. (18) J.L. Romero et al. / Diffusion in a class of exactly solvable non-harmonic potentials 87 Plate II. Two-dimensional projection of the manifold, defined by the vector (0.6, - 0.6, 5). By studying the same question with an exactly solvable model with a symmetrical unstable region, as one of the models proposed by Hongler and Zheng [11], one obtains that the minimum time for bifurcation corresponds exactly to an initial position located at the unstable point. For an asymmetrical unstable region there is a displacement of x0m towards the point where the curvature is greater, i.e. the point where the absolute value of the second derivative of U(x) is greater. 7. Final comment In this paper we have limited ourselves to some aspects of the behaviour of the solutions for an exactly solvable model of the Kramers type. By considering different initial conditions it has been possible to observe phenomena induced by nonlinearity, such as boomerang behaviour of the mean path and splitting of the probability density. We have quantitative estimates both for the onset time t B of this splitting and the time t* which