- Email: [email protected]
Diffusion of Brownian particles: dependence on the structure of the periodic potentials
Solid State Ionics 159 (2003) 331 – 343 www.elsevier.com/locate/ssi
Diffusion of Brownian particles: dependence on the structure of the periodic potentials A. Asaklil, Y. Boughaleb, M. Mazroui *, M. Chhib, L. El Arroum Laboratoire de Physique de la Matie`re Condense´e, Faculte´ des Sciences Ben M’sik, BP 7955, Casablanca, Morocco Received 18 July 2002; received in revised form 13 November 2002; accepted 14 November 2002
Abstract We present here a study of the diffusive motion of a particle submerged in bistable and metastable periodic potentials within the framework of the Brownian motion theory. This is done through an investigation of the quasi-elastic peak in the dynamic structure factor Ss( q,x). Its width is found to contain valuable information on the mechanism for the diffusion process. For this study, we use the Fokker-Planck equation, which is solved numerically by the matrix-continued fraction method (MCFM) for wide system parameters’ range and for the both types of potentials. It is the purpose of the present work to study how the transport properties of the system are modified by going from bistable to metastable periodic potential. Our finding results indicate large difference between transport properties in bistable and metastable potentials essentially at low temperature. In fact, in the former case, the mechanism process results from a combination of inter-cell liquid-like and intra-cells hopping motion of the particle. While for the second case, the diffusive process consists only of hopping motion, with different jump lengths, inside and between the cells. So, in metastable potential, a simple jump model describes the diffusive motion quite well. Further, a direct comparison between the numerical diffusion coefficient D and the analytical approximation for both potential shapes in the intermediate friction limit is presented and discussed. D 2002 Elsevier Science B.V. All rights reserved. PACS: 05.40.a Fluctuation phenomena, random processes, noise and Brownian motion; 66.30.Dn Theory of diffusion and Ionic conduction in solids; 82.20.Fd Stochastic and trajectory models, other theories and models; 05.60.k Transport processes Keywords: Diffusion; Brownian particles; Periodic potentials
1. Introduction The one-dimensional diffusion problem of a Brownian particle in a periodic potential [1] is ubiquitous in almost all scientific areas. It represents a model that can be applied to numerous systems, ranging from superionic conductors and intercalation compounds to submonolayer films adsorbed on surfaces of crystal* Corresponding author. E-mail address: [email protected] (M. Mazroui).
line substrates. The common features of these systems are that they consist of two species of particles; one species is fixed around certain equilibrium sites and forms a regular lattice while the other species is mobile and moves through this lattice. These systems are characterized by diffusion constants in the range of liquids, but the periodicity, due to the crystalline substrate, is not completely lost. If the diffusing particle is sufficiently massive, the diffusion problem can be treated within a classical approach. The theoretical method to understand this important problem
0167-2738/02/$ - see front matter D 2002 Elsevier Science B.V. All rights reserved. doi:10.1016/S0167-2738(02)00890-1
332
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
has been developed mainly in two directions. The first is numerical molecular dynamics or the Monte Carlo simulation method [2]. The second method is an analytical study based on the Langevin or the Fokker-Planck equation (also called a forward Kolmogorov equation) describing the Brownian motion in a periodic potential [3,4]. The largest part of the results on the Fokker-Planck dynamics problem has been obtained for a simple cosine potential [5,6]. While the diffusion problem in symmetric and asymmetric double-well potentials has been much less thoroughly investigated [7 –9]. Experimental studies performed on different superionic conductors materials show that a number of these ionic compounds are known in which diffusion of ions occurs in double-well potentials. A classical examples include diffusion in h aluminas (h-Al2O3) and RbAg4I5. The major treatments with these potential forms are regarded only for a Kramers problem, in order to elucidate some points in reaction-rate theory. Due to his important contribution, the escape-rate problem is now commonly known as the Kramers problems [10,11] and the corresponding Fokker-Planck equation is called the Kramers equation [3]. Here, we shall be mainly concerned with the Fokker-Planck dynamics of a classical particle moving in both symmetric bistable and metastable potentials, which are of primary interest in the study of ionic diffusion in solids and of chemical kinetics [10]. This is conveniently studied through the quasi-elastic peak in the dynamical structure factor Ss( q,x), which is proportional to the quasi-elastic scattering intensity both in neutron [12,13] and in atom scattering experiments [14]. Its width is found to contain valuable information on the mechanism for the diffusion process. The purpose of this paper is twofold. Firstly, we want to elucidate the mechanism for the diffusion process and we shall do that by investigating the quasi-elastic peak in Ss( q,x). Secondly, we will compare and discuss the difference between our results obtained both in bistable and metastable potentials. The paper is structured as follows. In Section 2, we describe briefly the matrix-continued fraction method for the solution of the Fokker-Planck equation. The numerical method is applied in the case of bistable and metastable periodic potential. Section 3 contains the description of two potentials used in the numerical calculations. In Sections 4 and 5, we establish a
comparison between the diffusion process in the symmetric bistable and in the metastable potential by computing two important physical quantities for the two forms of periodic potential: the half-width k( q) of quasi-elastic peak of Ss( q,x) and the diffusion coefficient D. The numerical results for the diffusion coefficient are also compared with those obtained by analytical approximation. The conclusions are outlined in Section 6.
2. Model and method: the matrix-continued fraction method (MCFM) The Fokker-Planck equation is of course not the only equation to describe the evolution of the distribution function, but it is one of the simplest equations for continuous macroscopic variables. The diffusing particle is subjected to three forces: a periodic deterministic force, derived from the potential V(x), a frictional force cmv (c is the friction coefficient which is often treated as a phenomenological parameter), and a Gaussian white noise, related to the friction via the fluctuation – dissipation theorem [15]. In these conditions, the phase space probability density f(x,v,t) reads: Bf ðx; v; tÞ ¼ LFP f ðx; v; tÞ Bt
ð1Þ
where LFP is the Fokker-Planck operator: LFP
B 1 BV ðxÞ B B kB T B þ þc vþ ¼ v Bx m Bx Bv Bv m Bv ð2Þ
f(x,v,t)dxdv is the probability of finding the particle in the phase space element between (x,v) and (x + dx,v + dv), T is the temperature and kB is the Boltzmann constant. Various methods of solution of this equation have been used such as transformation of a Fokker-Planck equation to a Schro¨dinger equation, numerical integration methods and analytical solutions for certain model potentials. However, the most efficient numerical method, adopted here, to solve this equation is the matrix-continued fraction method (MCFM), developed by Risken (see Ref. [3], where a complete
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
description of the method with many applications is given). In the periodic case, the MCFM is based on the expansion of the solution into a basis set of plane waves for the position variable and of Hermite functions for the velocity variable. In this way, the Green function (i.e. of the conditional probability distribution) of the FPE and the dynamical structure factor Ss can be obtained. From Ss( q,x), many correlation functions(the velocity correlation function, the meansquare displacement, jump rate, and the diffusion coefficient) may be computed via Kubo relations [15 – 17]. The method used here is quite general; it can be applied to different lattices, position-dependent friction [18] and symmetric and asymmetric doublewell potentials [7– 9]. In this paper, we will analyze in detail the results concerning the diffusion process, in the case of symmetric bistable and metastable periodic potentials. Their analytical expression is given by considering the first two terms of the Fourier expansion of the periodic potential V ðxÞ ¼ A cosðq0 xÞ þ B cosð2q0 xÞ
ð3Þ
where q0=(2p/a) denotes the reciprocal lattice vector and a is the lattice constant. The details about the choice of this form of potential can be found in the following section. The details of the numerical method are given in Refs. [3,7]. Here we give only the essential lines. First of all, it is convenient to define some dimensionless variable: rffiffiffiffiffiffiffiffiffi rffiffiffiffiffiffiffiffiffi 2p 2p kB T a m x¯ ¼ x; ¯t ¼ t; c¯ ¼ c a a m 2p kB T rffiffiffiffiffiffiffiffiffi m V ðxÞ ð4Þ v¯ ¼ v; V¯ ð¯xÞ ¼ kB T kB T with this choice for x, the unit cell goes from p to p and from 1/2 to 1/2 in the real and reciprocal axes, respectively. In the following, we will use these dimensionless quantities, if not otherwise specified. By analogy to the notation adopted by the groups of Carrati et al. [18] and Dieterich et al. [19], the dimensionless variables will be rewritten without overline and the Fokker-Planck operator becomes LFP
B BV ðxÞ B B B þ þc vþ ¼ v : Bx Bx Bv Bv Bv
ð5Þ
333
The stochastic theory based on the FPE is able to describe different diffusion mechanisms. Quasi-continuous diffusion and hopping mechanism (by single or multiple jumps) correspond to different ranges of the friction and of the temperature or, more physically, to different ratio between some typical times scales. The studies performed by means of FPE in periodic potential showed that, in general, at fixed q, the dynamical structure factor Ss( q,x) presents a quasielastic peak centered around x = 0, which properly describes the diffusive motion of the particle. At very high friction, Ss( q,x) is a monotonically decreasing function of x, but for low value of friction and at high barriers, the dynamical structure factor presents secondary peaks connected to the oscillatory dynamics around the potential minima [19]. It can be defined as the time Fourier transform of the characteristic function F( q,t); 1 S ðq; xÞ ¼ 2p s
Z
þl
eixt Fðq; tÞdt;
ð6Þ
l
where Fðq; tÞ ¼ hexpðiqðxðtÞ xð0ÞÞÞi:
ð7Þ
The angle brackets refer to the thermal average F( q,t); ¯hq and ¯hx are the momentum transfer and the energy loss, respectively. If the FPE is solved, the latter Eq. (7) can be developed as:
Fðq; tÞ ¼
Z
Z
þp
p
Z
l
dx0 l
Z
l
dv0
l
dx l
dvPst ðx0 ; v0 Þ
l
Pðx; v; t=x0 ; v0 ; 0Þexpðiqðx x0 ÞÞ
ð8Þ
where Pst is the stationary probability density given by the Boltzmann distribution and P(x,v,t/x0,v0,0) is the conditional probability of having the particle in x,v at time t, if it was in x0,v0 at time zero. P(x,v,t/x0,v0,0) is then the Green function of the probability density f(x,v,t) (the solution of the FPE with initial d-condition in both position and velocity). The dynamical structure factor can be obtained from the characteristic function after expanding the probability density into a basis set of plane waves for the position variable, and of hermit functions for the velocity variable which form together a complete
334
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
system. The final expression for S( q,x) can then be written as: ( ) l X p;r ˜ Sðq; xÞ ¼ N R * ð9Þ G ðk; ixÞMpl Mrl 0;0
p;r¼l
where q= (2p/a)(k + l) with l integer, (1/2) < k V (1/2), and Z p 1 V ðxÞ þ irx dx ð10Þ Mr ¼ exp 2p p 2
˜ 0;0 ðk; ixÞ ¼ G
N 1 ¼
Z
dxexp½V ðxÞ :
H p;r ðkÞ ¼ ðp þ kÞdp;r þ
A ðdp;r1 dp;rþ1 Þ 4kB T
iterations may be very large. The matrix-continued fractions seem to converge even for very small damping constants.
In this section, we address in greater detail the physical interests of the structured potentials used here. Their common analytical expression is given by
A ðdp;r1 dp;rþ1 Þ 4kB T
B ðdp;r2 dp;rþ2 Þ 2kB T
ð12Þ
3. Applications of bi- and metastable potential models
B ðdp;r2 dp;rþ2 Þ þ 2kB T
H˜ p;r ðkÞ ¼ ðp þ kÞdp;r
ð11Þ
The integrals Mr and N are completely defined by the potential V(x). The Green function ˜ 0,0(k,ix) is expressed by a matrix-continued fracG tion:
I rffiffiffiffiffiffiffiffiffi : m I a rffiffiffiffiffiffiffiffiffi i 2p H˜ xI þ H m I kB T a rffiffiffiffiffiffiffiffiffi H˜ 2p kB T ðix þ cÞI þ 2H a m ðix þ 2cÞI þ . . . 2p kB T
˜ are In Eq. (12), I is the identity matrix and H and H found to have the following expressions [7]
is the modified Bessel function. The normalization factor N is given by
V ðxÞ ¼ A cosðxÞ þ B cosð2xÞ: ð13Þ
where dp,r represents the usual Kronecker symbol. ˜ are equal to For free Brownian motion both H and H the differential operator B/Bx. The continued fraction ˜ 0,0 is explicitly written above. needed for obtaining G As for the Brownian motion problem in simple cosine potential, the matrix-continued fractions may be calculated numerically. The speed of convergence of the method depends on the friction constant C (C = 2pc/x0 with x0=(2p/a)(V1/2m)1/2). For large friction, only a few terms need be taken into account whereas for very small friction, the number of
ð14Þ
Dimensionless units are employed and the potential is normalized. These models represent a large class of potentials of interest in physical problems. The constants A and B determine the amplitudes of both potential barriers (V1 and V2) and are chosen such that the barrier energy V1 of the periodic potential is equal to 0.1 eV (i.e. of the order of the one found in superionic conductors). The computations were performed using this potential for various values of D which is defined as the ratio of the two different barriers of potential V(x): D = V2/V1. In order to take different values of D, we vary only the second barrier V2. So, the ratio can vary between zero and one. By these form potentials, it is
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
possible to give a reasonably accurate description of the systems. In fact, it can describe a real situation of many superionic conductors. Our investigations were performed using two different forms of potential. In the first one, the choice of value of two constants A and B in Eq. (14) is positive. In this case, we speak about a symmetric bistable potential (see Fig. 1). When a particle moves in this type of periodic potential, two different mechanisms are possible, depending on the importance of the barrier height V2 compared to the thermal energy kBT. For the low value of barrier V2, the particle performs small-amplitude oscillations around the well bottoms and jumps from well to another, by overcoming the second barrier height V1. If the barrier V2 is sufficiently high, the particle makes a jump from a well to another (jump from a to b and from b to a + 1). In a previous paper [20], we have shown that for strongly interacting Brownian particles in two-dimensional periodic potential, the effective potential computed along the direction where the system is incommensurate presents the same shape as the one presented in the Fig. 1. Indeed, in the low temperature regime, the effective potential exhibits a complicated structure characterized by the appearance of new equilibrium sites. Most ions are displaced away from the lattice sites. Such structure was observed by the means of X-ray diffraction along the tunnel axis deduced by Weber and Schultz [21]. Also, the experiment performed on RbAg4I5 has shown effectively that the potential felt by a silver ion located at a tetrahedral site has two different barriers [22].
335
Fig. 2. Inequivalent site model with deep a-type and shallow b-type wells of symmetric metastable potential V(x).
In the second one, by inverting the symmetric bistable potential, we get the metastable periodic potential, which we present in Fig. 2. In this form of potential, we notice two inequivalent sites a and b. Once the particle has escaped from a-type well, it may be retrapped in a nearest-neighbor one (b-type) or it may continue its flight towards another a-type. The particle jumps rapidly from a b-type to an unoccupied a-type well, but the jump from a-type to b-type well is much slower. Very recently, Montalenti and Ferrando [23] showed that for Au atoms, an important contribution to diffusion (especially for long-jump events) comes from metastable walks, in which an adatom starts from a minimum in the channel reaching another in channel minimum by an out-of-channel trajectory, passing through a metastable potential energy region. For all these reasons, we choose the expression of potential given by Eq. (14), in order to describe a real situation in wide variety of systems and specially the superionic conductors which have attracted considerable interest for the last several years [24,25].
4. Half-width of the quasi-elastic peak: diffusion mechanisms
Fig. 1. Structure of symmetric bistable potential V(x). Different types of jump of a Brownian particle are shown.
As mentioned in Section 2, information on the microscopic diffusion mechanism can be obtained from quasi-elastic line of the dynamical structure
336
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
factor Ss( q,x). Generally, the latter quantity consist of a quasi-elastic line and of oscillatory side peaks, which may become overdamped for large friction or high temperature. Let us first discuss the oscillatory part of the spectrum. The weight of the resonance near the harmonic frequency x0 increases strongly with q. At still larger q ( q2kBT/mx0g1), the intensity ceases to increase and finally decreases in favor of higher order harmonics. The peak moves towards lower frequencies as T increases and finally merges in a broad background as kBT becomes comparable to V1. This is due essentially to the anharmonic motion of the particle. Let us turn to the quasi-elastic line of Ss( q,x) whose half-width contains valuable information about the interaction of the mobile ions with the rigid framework and with each other. In effect, the behavior of the width of the quasi-elastic line of the dynamical structure factor with the wavescattering vector q and temperature should provide valuable indication on the diffusion process. In particular, it should clarify the question whether a continuous diffusion or a jump diffusion description is appropriate. This question is interesting in connection with superionic conductors where the potential barriers height V1 and kBT may become of the same order. The study of the dynamical structure factor in the reciprocal space and in particular of its
half-width at half maximum (HWHM) k( q) shows more interesting features, depending on the form of potential. Let us analyze the behavior of the HWHM k( q) in two cases of potential (bistable and metastable potential) and for different values of the ratio of two potential barriers D. In Fig. 3, the HWHM of the quasi-elastic peak as a function of the scattering wave-vector q is plotted in the case of bistable potential for different values of D and at low temperature (kBT/V1b1). The calculation of k( q) shows that the FPE at low temperature describes a thermally activated jump diffusion. Thus, for D = 0, the half-width k( q) is found to be approximately periodic in q with a period equal to q0 = 2p/a. For reciprocal lattice vectors q/q0 = 1,2,3,. . ., the halfwidth vanishes and reaches its maxima at q/q0 = 1/ 2,3/2. . ., etc. In this case, we recover the results obtained by Ferrando et al. [6] for an usual cosine potential. The system is described by a jump diffusion process and the jump length is close to the lattice constant a. This behavior can be derived exactly from the hopping model where the only jumps considered are those connecting nearest neighbor sites. The behavior of k( q) for D = 3/4 presents a complicated structure which is far from being simply a periodic function of q, indicating that the dynamic of mobile ions cannot be described by a simple
Fig. 3. Half-width of the quasi-elastic line of dynamical structure factor as a function of wave vector, for different shapes of symmetric bistable potential and at high friction limit. The parameters for this figure are: V1 = 0.1 eV and kBT = V1/6.
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
jump. Furthermore by increasing still more the value of D(D = 9/10), we note clearly that the width of the minima at q0 and 3q0 become narrower; this implies that the jump probabilities of length a decrease more and more. Consequently, the jump-length probabilities are not equivalent for all values of D. However, the diffusion mechanism becomes dominated by jumps of length a/2, inside and between the unit cells. While for D = 1(V1 = V2), the periodic potential becomes simple with one barrier and with a spatial period equal to a/2. The half-width k( q) recovers then its cosine shape. The diffusion mechanism is entirely represented by instantaneous jumps from an equilibrium site to another one with jump length a/2. In this case, k( q) describes a thermally activated jump diffusion. In Fig. 4, we report k( q) for a symmetric metastable potential; we see that the behavior of k( q) is practically the same as in Fig. 3 except for D = 3/4. For this last value of D, the behavior of k( q) is close to a simple cosine and the dynamic of mobile ions is described by a simple jump diffusion process with jump length close to lattice constant a. For high value of D, the diffusion process of Brownian particle
337
moving in a symmetric metastable potential is entirely represented by a jump motion with different value of jump length (a and a/2). In Figs. 5 and 6, the HWHM k( q) is represented as a function of the wave vector respectively in the case of bistable and metastable periodic potential and for low value of D (1/10 V D V 1/2). A comparison between this two figures shows that the behavior of k( q) in a symmetric bistable potential differs strongly from the one in a metastable potential. In fact, in Fig. 5, we see that at very low value of D (D = 1/10), the shape of k( q) is very close to a simple cosine but by increasing slightly the value of D, we notice the appearance of the pointed peaks which become more important for the high value of q and tend to the hydrodynamic limit. This reflect the fact that the diffusion mechanism is described by jump diffusion model with jump length close to the constant a and also by the liquid-like motion inside the unit cell. While in Fig. 6, contrary to Fig. 5, the behavior of k( q) shows that for low value of D (0 < D V 1/2), the system is described only by a jumps of length a. The general features of the behavior of the HWHM k( q) in metastable potential with increas-
Fig. 4. The q-dependance of the half-width of the quasi-elastic line k( q) of Ss( q,x) associated with different form of symmetric metastable potential. The parameters are the same as for Fig. 3.
338
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
Fig. 5. Half-width of the quasi-elastic line of Ss( q,x) as a function of wave vector, for different high value of D (D is the ratio of barriers bistable potential). The parameters are the same as for Fig. 3.
ing the temperature (V1 = 2kBT) are represented in Fig. 7. We observe that the form of k( q) is the same as that presented in Ref. [7] where we report the q-dependence of the half-width of the quasielastic line k( q) of Ss( q,x) associated with differ-
ent forms of the bistable potential and with the same parameters as for Fig. 7. The shape of k( q) is no more periodic reflecting that the diffusion process is far from being simply jump diffusion: the liquid-like motion of the particle dominates as
Fig. 6. The same as for Fig. 5 but for a symmetric metastable potential.
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
339
Fig. 7. Half-width of the quasi-elastic line of Ss( q,x) as a function of wave vector, for different form of symmetric metastable potential and at high temperature kBT = V1/2.
we increase the temperature independently of the value of D.
5. Diffusion coefficient: dependance on the structure periodic potentials In this section, we will present the results for the diffusion coefficient D of a Brownian particle moving in both symmetric bistable and metastable potentials. We are able to calculate this quantity which can describe the intrinsic properties of system, in three different limits: the high friction regime, in which the escape particle is dominated by spatial diffusion, the intermediate friction and the underdamped regime, characterized by energy diffusion. A complete study of diffusion in all damping regimes can be performed only by numerical methods developed in Section 2. So, the diffusion coefficient D is obtained as: D¼
a 4p
rffiffiffiffiffiffiffiffiffi m kðqÞ lim 2 kB T q!0 q
ð15Þ
where a is the lattice constant and k( q) is the halfwidth of the quasi-elastic peak of Ss( q,x). Eq. (15) means that the small q behavior of k( q), which is
usually investigated in neutron and atom scattering experiments, gives the diffusion coefficient. The results for the diffusion coefficient will be presented and discussed separately in the different damping regimes. The exact numerical results will be also compared to those of different approximations. Firstly, we represent the general features of the behavior of the diffusion coefficient D in a symmetric bistable potential with different values of friction C. In Fig. 8, we report D as a function of the ratio of two potential barriers D (D = V1/V2) at low temperature. Several conclusions can be drawn directly from a visual inspection of this figure. The behavior of D is practically the same for different values of C. We can also see in this figure a fastest decay of D from its initial value D(0) within a short interval of D (0 < D < 1/10). This is due essentially to the form of periodic potential, which seems to be rectangular as the width of its well increases; the vibration frequency tends rapidly towards zero and explains the steep decrease of D for small values of D. Three different diffusion process are shown. For low value of D (0 V D < 1/4), the system can be described by jump diffusion process with jump length close to lattice constant a. In the range of intermediate value of D (1/4 < D < 3/4), the diffusion process consist of a superposition of both hopping and
340
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
Fig. 8. Diffusion coefficient D of Brownian particle moving in a bistable periodic potential as a function of the ratio D and for different regime of friction. The parameters for this figure are: V1 = 0.1 eV and kBT = V1/6.
liquid-like motions inside the unit cell. In the last region (3/4 < D V 1), the diffusion process is also described by hopping motions but with jump length equal to a/2. Secondly, we analyze the results concerning the diffusion coefficient D of Brownian particle submerged in a symmetric metastable periodic potential. The results of our calculations are presented in Fig. 9
where the quantity D is reported down to different values of friction. Contrary to Fig. 8, the behavior of D depends on the range of friction, specially for the low value of D. In order to make a quantitative comparison, we have plotted in one figure (Figs. 10 and 11) D versus D for both potentials at high and low friction limits. The results corresponding to bistable and metastable potentials are represented by solid lines and
Fig. 9. Diffusion coefficient D in a metastable periodic potential as a function of the ratio D for different values of friction. The parameters are the same as Fig. 8.
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
341
and essentially in the intermediate damping regime, which is between the case of low friction and the case of high friction. Their analytical approximation for diffusion coefficient is derived in the framework of the Linear Response Theory Dð0Þ ¼ aðdÞDHF ð0Þ þ bðdÞDLF ð0Þ
Fig. 10. Diffusion coefficient D of the Brownian particle in both bistable and metastable periodic potential as function of D. The parameters are: V1 = 0.1 eV, kBT = V1/6 and C = 36.
dashed lines, respectively. From Fig. 10, we observe that we have coincidence between the two curves in the regime of high friction. This agreement can be well explained in terms of vibration frequencies at the bottom xb and the top xs point of the potential V(x). In fact, the product xbxs has the same value in both bistable and metastable periodic potential for every value of D. This is due essentially to the fact that we can get the metastable potential only by inversion of the bistable potential and vice versa. We recall that, in this regime of friction, the jump probability / (/ = D/ a2) is approximately defined by Kramers [26]: /~
xb xs ; 2pc
while for low friction limit, we observe a great deviation between the two results of the diffusion coefficient (Fig. 11). As above, this discrepancy can be easily understood from vibration frequency. The jump probability / is proportional to xb(/~xb). Hence, by increasing slightly the value of D, the vibration frequency at the bottom increases in the case of metastable potential and decreases in bistable potential. We note clearly that this interpretation is valid only for very values of D (D < 1/10).
ð16Þ
where DHF(0) and DLF(0) are analytical expressions of diffusion coefficient obtained for high and low friction limits, respectively, and are given elsewhere [29,30]. The weight parameters a(d) and b(d) are given in terms of the energy lost d between the crossing of a barrier and are given by the flowing expressions [27,28] aðdÞ ¼ 1 exp bðdÞ ¼ exp
d 2kB T
d : 2kB T
ð17Þ
The mean energy dissipated on a lattice spacing a can be calculated approximately from Langevin equation Z s Z a mcvdx ¼ mcv2 dt; s ¼ 2p=x0 : d¼ 0
0
In the case of symmetric double-well potential, the energy dissipated on a cell depends on the ratio
5.1. Comparison with analytical approximations It is interesting also to compare our numerical results with those of an analytical approximation. Very recently, Mazroui et al. [27,28] have given a detailed study of diffusion coefficient in wide damping range
Fig. 11. Same as Fig. 10 but for low value of friction C = 1.5.
342
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
Fig. 12. Comparison between analytical and numerical results of the diffusion coefficient in both bistable and metastable potential and for low value of friction. The parameters are: V1 = 0.1 eV, kBT = V1/6 and C = 1.5. ( – ): analytical approximation results. (- -): numerical (MCFM) results.
of two potential barriers and the dimensionless friction C. d ¼ dðDÞ ¼
CV2 CDV1 ¼ : 2 2
In Fig. 12, D is reported as a function of the ratio D in symmetric bistable and metastable potential for moderate values of friction (C = 1.5). From these two figures, we observe clearly that the analytical approximation result (Eq. (16)), obtained earlier, is qualitatively in good agreement with the numerical results (dashed line) for both shapes of potential.
Numerical results have been obtained both for biand metastable potentials, which are ubiquitous in almost all scientific areas. The results presented here indicate large difference between the transport properties in these two types of periodic potential. In effect, for metastable potential, the process of migration of a particle is very well described by a jump diffusion motion with different jump lengths. While for bistable potential, the diffusion process is more complicated and it is not possible to reproduce the jump diffusion model results. The corresponding results of bistable periodic potential demonstrate that the diffusion process combines liquid-like behavior with lattice-like properties. Further, we have calculated the diffusion coefficient D in both potentials in wide range of theory parameters: the strength of friction and the amplitude of the potential barriers. At high friction, the diffusion coefficient obtained in bistable potential close exactly to the one obtained for metastable potential independently of the ratio D. While at low friction except for D = 0,1 the diffusion coefficient behaves differently in these two types of potential. The excellent coincidence and the greater deviation of D is well explained in terms of vibration frequency at the bottom and at the saddle point of the periodic potential. Also, in the framework of this work, we have made a direct comparison between the numerical diffusion coefficient D and the analytical approximation for two potentials shape in intermediate friction limit. The agreement between the numerical results and the analytical approximation is qualitatively good.
6. Summary and conclusion In this paper, we have discussed the diffusive motion of a particle subject to symmetric double-well potential (both in bi- and metastable potentials). The investigation is essentially based on the Brownian motion theory (Fokker-Planck model) for these processes, which can be applied equally well to diffusion in solids. The Fokker-Planck equation has been numerically solved by using the matrix-continued fraction method. In this way, the Green function of the Fokker-Planck equation and the dynamical structure factor can be obtained. The diffusion coefficient as well as diffusion mechanisms can be extracted from calculations of the dynamical structure factor Ss( q,x).
References [1] A. Bunde, Z. Phys., B Condens. Matter 44 (1981) 225 – 232. [2] S. Paik, S. DasSarma, Surf. Sci. Lett., L 53 (1989) 208; J.D. Doll, H.K. McDowell, J. Chem. Phys. 77 (1982) 479. [3] H. Risken, The Fokker-Planck Equation, Springer, Berlin, 1989. [4] G. Wahnsto¨m, Surf. Sci. 159 (1985) 311. [5] R. Ferrando, R. Spadacini, G.E. Tommei, Phys. Rev., E 51 (1995) 1. [6] R. Ferrando, R. Spadacini, G.E. Tommei, Surf. Sci. 311 (1994) 411 – 421. [7] A. Asaklil, M. Mazroui, Y. Boughaleb, Eur. Phys. J., B 10 (1999) 91. [8] M. Hairibi, A. Asaklil, M. Mazroui, Y. Boughaleb, M.J. Condens. Matter l1 (1999) 1.
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343 [9] H. Risken, K. Voigtla¨nder, J. Stat. Phys. 41 (1985) 825; K. Voigtla¨nder, H. Risken, Chem. Phys. Lett. 105 (1984) 506. [10] P. Ha¨nngi, P. Talkner, M. Brokovec, Rev. Mod. Phys. 62 (1990) 251. [11] V.I. Mel’nikov, Phys. Rep. 209 (1991) 1. [12] J.W.M. Frenken, B.J. Hinch, J.P. Toennies, Ch. Wo¨ll, Phys. Rev., B 41 (1990) 938. [13] G. Radons, T. Geisel, J. Kell, Solid State Ionics 13 (75) (1984) 80 (North-Holland, Amsterdam). [14] A.C. Levi, R. Spadacini, G.E. Tommei, Surf. Sci. 121 (1982) 357. [15] R. Kubo, Rep. Prog. Phys. 29 (1966) 255. [16] J.P. Hansen, I.R. McDonald, Theory of Simple Liquids, Academic Press, London, 1986. [17] R. Ferrando, R. Spadacini, G.E. Tommei, G. Caratti, Physica, A 195 (1993) 506. [18] G. Caratti, R. Ferrando, R. Spadacini, G.E. Tommei, Phys. Rev., E 55 (1997) 4. [19] W. Dieterich, P. Fulde, I. Peschel, Adv. Phys. 29 (3) (1980) 527 – 605.
343
[20] M. Mazroui, A. Asaklil, Y. Boughaleb, J. Phys., I Fr. 7 (1997) 675 – 690. [21] H.P. Weber, H. Schulz, J. Chem. Phys. 85 (1986) 475. [22] K. Funke, in: A.L. Laskar, S. Chandra (Eds.), Superionic Solids and Solid Electrolytes, Academic Press, New York, 1989. [23] F. Montalenti, R. Ferrando, Phys. Rev., B 59 (1999) 8. [24] M. Aniya, Solid State Ionics 70 (1994) 673. [25] M. Aniya, Physica, B 239 (1997) 144. [26] H.A. Kramers, Physica 7 (1940) 284. [27] M. Mazroui, A. Asaklil, Y. Boughaleb, Surf. Sci. 409 (1998) 528 – 540; M. Mazroui, Y. Boughaleb, Int. J. Mod. Phys., B 15 (2000) 2193. [28] M. Mazroui, M. Hairibi, Y. Boughaleb, Phys. Lett., A 253 (1999) 227. [29] H. Risken, H.D. Vollmer, Phys. Lett. 6 (1979) 69. [30] A.R. Bishop, W. Dieterich, I. Peschel, Z. Phys., B Condens. Matter 33 (1979) 187.
Diffusion of Brownian particles: dependence on the structure of the periodic potentials A. Asaklil, Y. Boughaleb, M. Mazroui *, M. Chhib, L. El Arroum Laboratoire de Physique de la Matie`re Condense´e, Faculte´ des Sciences Ben M’sik, BP 7955, Casablanca, Morocco Received 18 July 2002; received in revised form 13 November 2002; accepted 14 November 2002
Abstract We present here a study of the diffusive motion of a particle submerged in bistable and metastable periodic potentials within the framework of the Brownian motion theory. This is done through an investigation of the quasi-elastic peak in the dynamic structure factor Ss( q,x). Its width is found to contain valuable information on the mechanism for the diffusion process. For this study, we use the Fokker-Planck equation, which is solved numerically by the matrix-continued fraction method (MCFM) for wide system parameters’ range and for the both types of potentials. It is the purpose of the present work to study how the transport properties of the system are modified by going from bistable to metastable periodic potential. Our finding results indicate large difference between transport properties in bistable and metastable potentials essentially at low temperature. In fact, in the former case, the mechanism process results from a combination of inter-cell liquid-like and intra-cells hopping motion of the particle. While for the second case, the diffusive process consists only of hopping motion, with different jump lengths, inside and between the cells. So, in metastable potential, a simple jump model describes the diffusive motion quite well. Further, a direct comparison between the numerical diffusion coefficient D and the analytical approximation for both potential shapes in the intermediate friction limit is presented and discussed. D 2002 Elsevier Science B.V. All rights reserved. PACS: 05.40.a Fluctuation phenomena, random processes, noise and Brownian motion; 66.30.Dn Theory of diffusion and Ionic conduction in solids; 82.20.Fd Stochastic and trajectory models, other theories and models; 05.60.k Transport processes Keywords: Diffusion; Brownian particles; Periodic potentials
1. Introduction The one-dimensional diffusion problem of a Brownian particle in a periodic potential [1] is ubiquitous in almost all scientific areas. It represents a model that can be applied to numerous systems, ranging from superionic conductors and intercalation compounds to submonolayer films adsorbed on surfaces of crystal* Corresponding author. E-mail address: [email protected] (M. Mazroui).
line substrates. The common features of these systems are that they consist of two species of particles; one species is fixed around certain equilibrium sites and forms a regular lattice while the other species is mobile and moves through this lattice. These systems are characterized by diffusion constants in the range of liquids, but the periodicity, due to the crystalline substrate, is not completely lost. If the diffusing particle is sufficiently massive, the diffusion problem can be treated within a classical approach. The theoretical method to understand this important problem
0167-2738/02/$ - see front matter D 2002 Elsevier Science B.V. All rights reserved. doi:10.1016/S0167-2738(02)00890-1
332
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
has been developed mainly in two directions. The first is numerical molecular dynamics or the Monte Carlo simulation method [2]. The second method is an analytical study based on the Langevin or the Fokker-Planck equation (also called a forward Kolmogorov equation) describing the Brownian motion in a periodic potential [3,4]. The largest part of the results on the Fokker-Planck dynamics problem has been obtained for a simple cosine potential [5,6]. While the diffusion problem in symmetric and asymmetric double-well potentials has been much less thoroughly investigated [7 –9]. Experimental studies performed on different superionic conductors materials show that a number of these ionic compounds are known in which diffusion of ions occurs in double-well potentials. A classical examples include diffusion in h aluminas (h-Al2O3) and RbAg4I5. The major treatments with these potential forms are regarded only for a Kramers problem, in order to elucidate some points in reaction-rate theory. Due to his important contribution, the escape-rate problem is now commonly known as the Kramers problems [10,11] and the corresponding Fokker-Planck equation is called the Kramers equation [3]. Here, we shall be mainly concerned with the Fokker-Planck dynamics of a classical particle moving in both symmetric bistable and metastable potentials, which are of primary interest in the study of ionic diffusion in solids and of chemical kinetics [10]. This is conveniently studied through the quasi-elastic peak in the dynamical structure factor Ss( q,x), which is proportional to the quasi-elastic scattering intensity both in neutron [12,13] and in atom scattering experiments [14]. Its width is found to contain valuable information on the mechanism for the diffusion process. The purpose of this paper is twofold. Firstly, we want to elucidate the mechanism for the diffusion process and we shall do that by investigating the quasi-elastic peak in Ss( q,x). Secondly, we will compare and discuss the difference between our results obtained both in bistable and metastable potentials. The paper is structured as follows. In Section 2, we describe briefly the matrix-continued fraction method for the solution of the Fokker-Planck equation. The numerical method is applied in the case of bistable and metastable periodic potential. Section 3 contains the description of two potentials used in the numerical calculations. In Sections 4 and 5, we establish a
comparison between the diffusion process in the symmetric bistable and in the metastable potential by computing two important physical quantities for the two forms of periodic potential: the half-width k( q) of quasi-elastic peak of Ss( q,x) and the diffusion coefficient D. The numerical results for the diffusion coefficient are also compared with those obtained by analytical approximation. The conclusions are outlined in Section 6.
2. Model and method: the matrix-continued fraction method (MCFM) The Fokker-Planck equation is of course not the only equation to describe the evolution of the distribution function, but it is one of the simplest equations for continuous macroscopic variables. The diffusing particle is subjected to three forces: a periodic deterministic force, derived from the potential V(x), a frictional force cmv (c is the friction coefficient which is often treated as a phenomenological parameter), and a Gaussian white noise, related to the friction via the fluctuation – dissipation theorem [15]. In these conditions, the phase space probability density f(x,v,t) reads: Bf ðx; v; tÞ ¼ LFP f ðx; v; tÞ Bt
ð1Þ
where LFP is the Fokker-Planck operator: LFP
B 1 BV ðxÞ B B kB T B þ þc vþ ¼ v Bx m Bx Bv Bv m Bv ð2Þ
f(x,v,t)dxdv is the probability of finding the particle in the phase space element between (x,v) and (x + dx,v + dv), T is the temperature and kB is the Boltzmann constant. Various methods of solution of this equation have been used such as transformation of a Fokker-Planck equation to a Schro¨dinger equation, numerical integration methods and analytical solutions for certain model potentials. However, the most efficient numerical method, adopted here, to solve this equation is the matrix-continued fraction method (MCFM), developed by Risken (see Ref. [3], where a complete
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
description of the method with many applications is given). In the periodic case, the MCFM is based on the expansion of the solution into a basis set of plane waves for the position variable and of Hermite functions for the velocity variable. In this way, the Green function (i.e. of the conditional probability distribution) of the FPE and the dynamical structure factor Ss can be obtained. From Ss( q,x), many correlation functions(the velocity correlation function, the meansquare displacement, jump rate, and the diffusion coefficient) may be computed via Kubo relations [15 – 17]. The method used here is quite general; it can be applied to different lattices, position-dependent friction [18] and symmetric and asymmetric doublewell potentials [7– 9]. In this paper, we will analyze in detail the results concerning the diffusion process, in the case of symmetric bistable and metastable periodic potentials. Their analytical expression is given by considering the first two terms of the Fourier expansion of the periodic potential V ðxÞ ¼ A cosðq0 xÞ þ B cosð2q0 xÞ
ð3Þ
where q0=(2p/a) denotes the reciprocal lattice vector and a is the lattice constant. The details about the choice of this form of potential can be found in the following section. The details of the numerical method are given in Refs. [3,7]. Here we give only the essential lines. First of all, it is convenient to define some dimensionless variable: rffiffiffiffiffiffiffiffiffi rffiffiffiffiffiffiffiffiffi 2p 2p kB T a m x¯ ¼ x; ¯t ¼ t; c¯ ¼ c a a m 2p kB T rffiffiffiffiffiffiffiffiffi m V ðxÞ ð4Þ v¯ ¼ v; V¯ ð¯xÞ ¼ kB T kB T with this choice for x, the unit cell goes from p to p and from 1/2 to 1/2 in the real and reciprocal axes, respectively. In the following, we will use these dimensionless quantities, if not otherwise specified. By analogy to the notation adopted by the groups of Carrati et al. [18] and Dieterich et al. [19], the dimensionless variables will be rewritten without overline and the Fokker-Planck operator becomes LFP
B BV ðxÞ B B B þ þc vþ ¼ v : Bx Bx Bv Bv Bv
ð5Þ
333
The stochastic theory based on the FPE is able to describe different diffusion mechanisms. Quasi-continuous diffusion and hopping mechanism (by single or multiple jumps) correspond to different ranges of the friction and of the temperature or, more physically, to different ratio between some typical times scales. The studies performed by means of FPE in periodic potential showed that, in general, at fixed q, the dynamical structure factor Ss( q,x) presents a quasielastic peak centered around x = 0, which properly describes the diffusive motion of the particle. At very high friction, Ss( q,x) is a monotonically decreasing function of x, but for low value of friction and at high barriers, the dynamical structure factor presents secondary peaks connected to the oscillatory dynamics around the potential minima [19]. It can be defined as the time Fourier transform of the characteristic function F( q,t); 1 S ðq; xÞ ¼ 2p s
Z
þl
eixt Fðq; tÞdt;
ð6Þ
l
where Fðq; tÞ ¼ hexpðiqðxðtÞ xð0ÞÞÞi:
ð7Þ
The angle brackets refer to the thermal average F( q,t); ¯hq and ¯hx are the momentum transfer and the energy loss, respectively. If the FPE is solved, the latter Eq. (7) can be developed as:
Fðq; tÞ ¼
Z
Z
þp
p
Z
l
dx0 l
Z
l
dv0
l
dx l
dvPst ðx0 ; v0 Þ
l
Pðx; v; t=x0 ; v0 ; 0Þexpðiqðx x0 ÞÞ
ð8Þ
where Pst is the stationary probability density given by the Boltzmann distribution and P(x,v,t/x0,v0,0) is the conditional probability of having the particle in x,v at time t, if it was in x0,v0 at time zero. P(x,v,t/x0,v0,0) is then the Green function of the probability density f(x,v,t) (the solution of the FPE with initial d-condition in both position and velocity). The dynamical structure factor can be obtained from the characteristic function after expanding the probability density into a basis set of plane waves for the position variable, and of hermit functions for the velocity variable which form together a complete
334
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
system. The final expression for S( q,x) can then be written as: ( ) l X p;r ˜ Sðq; xÞ ¼ N R * ð9Þ G ðk; ixÞMpl Mrl 0;0
p;r¼l
where q= (2p/a)(k + l) with l integer, (1/2) < k V (1/2), and Z p 1 V ðxÞ þ irx dx ð10Þ Mr ¼ exp 2p p 2
˜ 0;0 ðk; ixÞ ¼ G
N 1 ¼
Z
dxexp½V ðxÞ :
H p;r ðkÞ ¼ ðp þ kÞdp;r þ
A ðdp;r1 dp;rþ1 Þ 4kB T
iterations may be very large. The matrix-continued fractions seem to converge even for very small damping constants.
In this section, we address in greater detail the physical interests of the structured potentials used here. Their common analytical expression is given by
A ðdp;r1 dp;rþ1 Þ 4kB T
B ðdp;r2 dp;rþ2 Þ 2kB T
ð12Þ
3. Applications of bi- and metastable potential models
B ðdp;r2 dp;rþ2 Þ þ 2kB T
H˜ p;r ðkÞ ¼ ðp þ kÞdp;r
ð11Þ
The integrals Mr and N are completely defined by the potential V(x). The Green function ˜ 0,0(k,ix) is expressed by a matrix-continued fracG tion:
I rffiffiffiffiffiffiffiffiffi : m I a rffiffiffiffiffiffiffiffiffi i 2p H˜ xI þ H m I kB T a rffiffiffiffiffiffiffiffiffi H˜ 2p kB T ðix þ cÞI þ 2H a m ðix þ 2cÞI þ . . . 2p kB T
˜ are In Eq. (12), I is the identity matrix and H and H found to have the following expressions [7]
is the modified Bessel function. The normalization factor N is given by
V ðxÞ ¼ A cosðxÞ þ B cosð2xÞ: ð13Þ
where dp,r represents the usual Kronecker symbol. ˜ are equal to For free Brownian motion both H and H the differential operator B/Bx. The continued fraction ˜ 0,0 is explicitly written above. needed for obtaining G As for the Brownian motion problem in simple cosine potential, the matrix-continued fractions may be calculated numerically. The speed of convergence of the method depends on the friction constant C (C = 2pc/x0 with x0=(2p/a)(V1/2m)1/2). For large friction, only a few terms need be taken into account whereas for very small friction, the number of
ð14Þ
Dimensionless units are employed and the potential is normalized. These models represent a large class of potentials of interest in physical problems. The constants A and B determine the amplitudes of both potential barriers (V1 and V2) and are chosen such that the barrier energy V1 of the periodic potential is equal to 0.1 eV (i.e. of the order of the one found in superionic conductors). The computations were performed using this potential for various values of D which is defined as the ratio of the two different barriers of potential V(x): D = V2/V1. In order to take different values of D, we vary only the second barrier V2. So, the ratio can vary between zero and one. By these form potentials, it is
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
possible to give a reasonably accurate description of the systems. In fact, it can describe a real situation of many superionic conductors. Our investigations were performed using two different forms of potential. In the first one, the choice of value of two constants A and B in Eq. (14) is positive. In this case, we speak about a symmetric bistable potential (see Fig. 1). When a particle moves in this type of periodic potential, two different mechanisms are possible, depending on the importance of the barrier height V2 compared to the thermal energy kBT. For the low value of barrier V2, the particle performs small-amplitude oscillations around the well bottoms and jumps from well to another, by overcoming the second barrier height V1. If the barrier V2 is sufficiently high, the particle makes a jump from a well to another (jump from a to b and from b to a + 1). In a previous paper [20], we have shown that for strongly interacting Brownian particles in two-dimensional periodic potential, the effective potential computed along the direction where the system is incommensurate presents the same shape as the one presented in the Fig. 1. Indeed, in the low temperature regime, the effective potential exhibits a complicated structure characterized by the appearance of new equilibrium sites. Most ions are displaced away from the lattice sites. Such structure was observed by the means of X-ray diffraction along the tunnel axis deduced by Weber and Schultz [21]. Also, the experiment performed on RbAg4I5 has shown effectively that the potential felt by a silver ion located at a tetrahedral site has two different barriers [22].
335
Fig. 2. Inequivalent site model with deep a-type and shallow b-type wells of symmetric metastable potential V(x).
In the second one, by inverting the symmetric bistable potential, we get the metastable periodic potential, which we present in Fig. 2. In this form of potential, we notice two inequivalent sites a and b. Once the particle has escaped from a-type well, it may be retrapped in a nearest-neighbor one (b-type) or it may continue its flight towards another a-type. The particle jumps rapidly from a b-type to an unoccupied a-type well, but the jump from a-type to b-type well is much slower. Very recently, Montalenti and Ferrando [23] showed that for Au atoms, an important contribution to diffusion (especially for long-jump events) comes from metastable walks, in which an adatom starts from a minimum in the channel reaching another in channel minimum by an out-of-channel trajectory, passing through a metastable potential energy region. For all these reasons, we choose the expression of potential given by Eq. (14), in order to describe a real situation in wide variety of systems and specially the superionic conductors which have attracted considerable interest for the last several years [24,25].
4. Half-width of the quasi-elastic peak: diffusion mechanisms
Fig. 1. Structure of symmetric bistable potential V(x). Different types of jump of a Brownian particle are shown.
As mentioned in Section 2, information on the microscopic diffusion mechanism can be obtained from quasi-elastic line of the dynamical structure
336
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
factor Ss( q,x). Generally, the latter quantity consist of a quasi-elastic line and of oscillatory side peaks, which may become overdamped for large friction or high temperature. Let us first discuss the oscillatory part of the spectrum. The weight of the resonance near the harmonic frequency x0 increases strongly with q. At still larger q ( q2kBT/mx0g1), the intensity ceases to increase and finally decreases in favor of higher order harmonics. The peak moves towards lower frequencies as T increases and finally merges in a broad background as kBT becomes comparable to V1. This is due essentially to the anharmonic motion of the particle. Let us turn to the quasi-elastic line of Ss( q,x) whose half-width contains valuable information about the interaction of the mobile ions with the rigid framework and with each other. In effect, the behavior of the width of the quasi-elastic line of the dynamical structure factor with the wavescattering vector q and temperature should provide valuable indication on the diffusion process. In particular, it should clarify the question whether a continuous diffusion or a jump diffusion description is appropriate. This question is interesting in connection with superionic conductors where the potential barriers height V1 and kBT may become of the same order. The study of the dynamical structure factor in the reciprocal space and in particular of its
half-width at half maximum (HWHM) k( q) shows more interesting features, depending on the form of potential. Let us analyze the behavior of the HWHM k( q) in two cases of potential (bistable and metastable potential) and for different values of the ratio of two potential barriers D. In Fig. 3, the HWHM of the quasi-elastic peak as a function of the scattering wave-vector q is plotted in the case of bistable potential for different values of D and at low temperature (kBT/V1b1). The calculation of k( q) shows that the FPE at low temperature describes a thermally activated jump diffusion. Thus, for D = 0, the half-width k( q) is found to be approximately periodic in q with a period equal to q0 = 2p/a. For reciprocal lattice vectors q/q0 = 1,2,3,. . ., the halfwidth vanishes and reaches its maxima at q/q0 = 1/ 2,3/2. . ., etc. In this case, we recover the results obtained by Ferrando et al. [6] for an usual cosine potential. The system is described by a jump diffusion process and the jump length is close to the lattice constant a. This behavior can be derived exactly from the hopping model where the only jumps considered are those connecting nearest neighbor sites. The behavior of k( q) for D = 3/4 presents a complicated structure which is far from being simply a periodic function of q, indicating that the dynamic of mobile ions cannot be described by a simple
Fig. 3. Half-width of the quasi-elastic line of dynamical structure factor as a function of wave vector, for different shapes of symmetric bistable potential and at high friction limit. The parameters for this figure are: V1 = 0.1 eV and kBT = V1/6.
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
jump. Furthermore by increasing still more the value of D(D = 9/10), we note clearly that the width of the minima at q0 and 3q0 become narrower; this implies that the jump probabilities of length a decrease more and more. Consequently, the jump-length probabilities are not equivalent for all values of D. However, the diffusion mechanism becomes dominated by jumps of length a/2, inside and between the unit cells. While for D = 1(V1 = V2), the periodic potential becomes simple with one barrier and with a spatial period equal to a/2. The half-width k( q) recovers then its cosine shape. The diffusion mechanism is entirely represented by instantaneous jumps from an equilibrium site to another one with jump length a/2. In this case, k( q) describes a thermally activated jump diffusion. In Fig. 4, we report k( q) for a symmetric metastable potential; we see that the behavior of k( q) is practically the same as in Fig. 3 except for D = 3/4. For this last value of D, the behavior of k( q) is close to a simple cosine and the dynamic of mobile ions is described by a simple jump diffusion process with jump length close to lattice constant a. For high value of D, the diffusion process of Brownian particle
337
moving in a symmetric metastable potential is entirely represented by a jump motion with different value of jump length (a and a/2). In Figs. 5 and 6, the HWHM k( q) is represented as a function of the wave vector respectively in the case of bistable and metastable periodic potential and for low value of D (1/10 V D V 1/2). A comparison between this two figures shows that the behavior of k( q) in a symmetric bistable potential differs strongly from the one in a metastable potential. In fact, in Fig. 5, we see that at very low value of D (D = 1/10), the shape of k( q) is very close to a simple cosine but by increasing slightly the value of D, we notice the appearance of the pointed peaks which become more important for the high value of q and tend to the hydrodynamic limit. This reflect the fact that the diffusion mechanism is described by jump diffusion model with jump length close to the constant a and also by the liquid-like motion inside the unit cell. While in Fig. 6, contrary to Fig. 5, the behavior of k( q) shows that for low value of D (0 < D V 1/2), the system is described only by a jumps of length a. The general features of the behavior of the HWHM k( q) in metastable potential with increas-
Fig. 4. The q-dependance of the half-width of the quasi-elastic line k( q) of Ss( q,x) associated with different form of symmetric metastable potential. The parameters are the same as for Fig. 3.
338
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
Fig. 5. Half-width of the quasi-elastic line of Ss( q,x) as a function of wave vector, for different high value of D (D is the ratio of barriers bistable potential). The parameters are the same as for Fig. 3.
ing the temperature (V1 = 2kBT) are represented in Fig. 7. We observe that the form of k( q) is the same as that presented in Ref. [7] where we report the q-dependence of the half-width of the quasielastic line k( q) of Ss( q,x) associated with differ-
ent forms of the bistable potential and with the same parameters as for Fig. 7. The shape of k( q) is no more periodic reflecting that the diffusion process is far from being simply jump diffusion: the liquid-like motion of the particle dominates as
Fig. 6. The same as for Fig. 5 but for a symmetric metastable potential.
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
339
Fig. 7. Half-width of the quasi-elastic line of Ss( q,x) as a function of wave vector, for different form of symmetric metastable potential and at high temperature kBT = V1/2.
we increase the temperature independently of the value of D.
5. Diffusion coefficient: dependance on the structure periodic potentials In this section, we will present the results for the diffusion coefficient D of a Brownian particle moving in both symmetric bistable and metastable potentials. We are able to calculate this quantity which can describe the intrinsic properties of system, in three different limits: the high friction regime, in which the escape particle is dominated by spatial diffusion, the intermediate friction and the underdamped regime, characterized by energy diffusion. A complete study of diffusion in all damping regimes can be performed only by numerical methods developed in Section 2. So, the diffusion coefficient D is obtained as: D¼
a 4p
rffiffiffiffiffiffiffiffiffi m kðqÞ lim 2 kB T q!0 q
ð15Þ
where a is the lattice constant and k( q) is the halfwidth of the quasi-elastic peak of Ss( q,x). Eq. (15) means that the small q behavior of k( q), which is
usually investigated in neutron and atom scattering experiments, gives the diffusion coefficient. The results for the diffusion coefficient will be presented and discussed separately in the different damping regimes. The exact numerical results will be also compared to those of different approximations. Firstly, we represent the general features of the behavior of the diffusion coefficient D in a symmetric bistable potential with different values of friction C. In Fig. 8, we report D as a function of the ratio of two potential barriers D (D = V1/V2) at low temperature. Several conclusions can be drawn directly from a visual inspection of this figure. The behavior of D is practically the same for different values of C. We can also see in this figure a fastest decay of D from its initial value D(0) within a short interval of D (0 < D < 1/10). This is due essentially to the form of periodic potential, which seems to be rectangular as the width of its well increases; the vibration frequency tends rapidly towards zero and explains the steep decrease of D for small values of D. Three different diffusion process are shown. For low value of D (0 V D < 1/4), the system can be described by jump diffusion process with jump length close to lattice constant a. In the range of intermediate value of D (1/4 < D < 3/4), the diffusion process consist of a superposition of both hopping and
340
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
Fig. 8. Diffusion coefficient D of Brownian particle moving in a bistable periodic potential as a function of the ratio D and for different regime of friction. The parameters for this figure are: V1 = 0.1 eV and kBT = V1/6.
liquid-like motions inside the unit cell. In the last region (3/4 < D V 1), the diffusion process is also described by hopping motions but with jump length equal to a/2. Secondly, we analyze the results concerning the diffusion coefficient D of Brownian particle submerged in a symmetric metastable periodic potential. The results of our calculations are presented in Fig. 9
where the quantity D is reported down to different values of friction. Contrary to Fig. 8, the behavior of D depends on the range of friction, specially for the low value of D. In order to make a quantitative comparison, we have plotted in one figure (Figs. 10 and 11) D versus D for both potentials at high and low friction limits. The results corresponding to bistable and metastable potentials are represented by solid lines and
Fig. 9. Diffusion coefficient D in a metastable periodic potential as a function of the ratio D for different values of friction. The parameters are the same as Fig. 8.
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
341
and essentially in the intermediate damping regime, which is between the case of low friction and the case of high friction. Their analytical approximation for diffusion coefficient is derived in the framework of the Linear Response Theory Dð0Þ ¼ aðdÞDHF ð0Þ þ bðdÞDLF ð0Þ
Fig. 10. Diffusion coefficient D of the Brownian particle in both bistable and metastable periodic potential as function of D. The parameters are: V1 = 0.1 eV, kBT = V1/6 and C = 36.
dashed lines, respectively. From Fig. 10, we observe that we have coincidence between the two curves in the regime of high friction. This agreement can be well explained in terms of vibration frequencies at the bottom xb and the top xs point of the potential V(x). In fact, the product xbxs has the same value in both bistable and metastable periodic potential for every value of D. This is due essentially to the fact that we can get the metastable potential only by inversion of the bistable potential and vice versa. We recall that, in this regime of friction, the jump probability / (/ = D/ a2) is approximately defined by Kramers [26]: /~
xb xs ; 2pc
while for low friction limit, we observe a great deviation between the two results of the diffusion coefficient (Fig. 11). As above, this discrepancy can be easily understood from vibration frequency. The jump probability / is proportional to xb(/~xb). Hence, by increasing slightly the value of D, the vibration frequency at the bottom increases in the case of metastable potential and decreases in bistable potential. We note clearly that this interpretation is valid only for very values of D (D < 1/10).
ð16Þ
where DHF(0) and DLF(0) are analytical expressions of diffusion coefficient obtained for high and low friction limits, respectively, and are given elsewhere [29,30]. The weight parameters a(d) and b(d) are given in terms of the energy lost d between the crossing of a barrier and are given by the flowing expressions [27,28] aðdÞ ¼ 1 exp bðdÞ ¼ exp
d 2kB T
d : 2kB T
ð17Þ
The mean energy dissipated on a lattice spacing a can be calculated approximately from Langevin equation Z s Z a mcvdx ¼ mcv2 dt; s ¼ 2p=x0 : d¼ 0
0
In the case of symmetric double-well potential, the energy dissipated on a cell depends on the ratio
5.1. Comparison with analytical approximations It is interesting also to compare our numerical results with those of an analytical approximation. Very recently, Mazroui et al. [27,28] have given a detailed study of diffusion coefficient in wide damping range
Fig. 11. Same as Fig. 10 but for low value of friction C = 1.5.
342
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343
Fig. 12. Comparison between analytical and numerical results of the diffusion coefficient in both bistable and metastable potential and for low value of friction. The parameters are: V1 = 0.1 eV, kBT = V1/6 and C = 1.5. ( – ): analytical approximation results. (- -): numerical (MCFM) results.
of two potential barriers and the dimensionless friction C. d ¼ dðDÞ ¼
CV2 CDV1 ¼ : 2 2
In Fig. 12, D is reported as a function of the ratio D in symmetric bistable and metastable potential for moderate values of friction (C = 1.5). From these two figures, we observe clearly that the analytical approximation result (Eq. (16)), obtained earlier, is qualitatively in good agreement with the numerical results (dashed line) for both shapes of potential.
Numerical results have been obtained both for biand metastable potentials, which are ubiquitous in almost all scientific areas. The results presented here indicate large difference between the transport properties in these two types of periodic potential. In effect, for metastable potential, the process of migration of a particle is very well described by a jump diffusion motion with different jump lengths. While for bistable potential, the diffusion process is more complicated and it is not possible to reproduce the jump diffusion model results. The corresponding results of bistable periodic potential demonstrate that the diffusion process combines liquid-like behavior with lattice-like properties. Further, we have calculated the diffusion coefficient D in both potentials in wide range of theory parameters: the strength of friction and the amplitude of the potential barriers. At high friction, the diffusion coefficient obtained in bistable potential close exactly to the one obtained for metastable potential independently of the ratio D. While at low friction except for D = 0,1 the diffusion coefficient behaves differently in these two types of potential. The excellent coincidence and the greater deviation of D is well explained in terms of vibration frequency at the bottom and at the saddle point of the periodic potential. Also, in the framework of this work, we have made a direct comparison between the numerical diffusion coefficient D and the analytical approximation for two potentials shape in intermediate friction limit. The agreement between the numerical results and the analytical approximation is qualitatively good.
6. Summary and conclusion In this paper, we have discussed the diffusive motion of a particle subject to symmetric double-well potential (both in bi- and metastable potentials). The investigation is essentially based on the Brownian motion theory (Fokker-Planck model) for these processes, which can be applied equally well to diffusion in solids. The Fokker-Planck equation has been numerically solved by using the matrix-continued fraction method. In this way, the Green function of the Fokker-Planck equation and the dynamical structure factor can be obtained. The diffusion coefficient as well as diffusion mechanisms can be extracted from calculations of the dynamical structure factor Ss( q,x).
References [1] A. Bunde, Z. Phys., B Condens. Matter 44 (1981) 225 – 232. [2] S. Paik, S. DasSarma, Surf. Sci. Lett., L 53 (1989) 208; J.D. Doll, H.K. McDowell, J. Chem. Phys. 77 (1982) 479. [3] H. Risken, The Fokker-Planck Equation, Springer, Berlin, 1989. [4] G. Wahnsto¨m, Surf. Sci. 159 (1985) 311. [5] R. Ferrando, R. Spadacini, G.E. Tommei, Phys. Rev., E 51 (1995) 1. [6] R. Ferrando, R. Spadacini, G.E. Tommei, Surf. Sci. 311 (1994) 411 – 421. [7] A. Asaklil, M. Mazroui, Y. Boughaleb, Eur. Phys. J., B 10 (1999) 91. [8] M. Hairibi, A. Asaklil, M. Mazroui, Y. Boughaleb, M.J. Condens. Matter l1 (1999) 1.
A. Asaklil et al. / Solid State Ionics 159 (2003) 331–343 [9] H. Risken, K. Voigtla¨nder, J. Stat. Phys. 41 (1985) 825; K. Voigtla¨nder, H. Risken, Chem. Phys. Lett. 105 (1984) 506. [10] P. Ha¨nngi, P. Talkner, M. Brokovec, Rev. Mod. Phys. 62 (1990) 251. [11] V.I. Mel’nikov, Phys. Rep. 209 (1991) 1. [12] J.W.M. Frenken, B.J. Hinch, J.P. Toennies, Ch. Wo¨ll, Phys. Rev., B 41 (1990) 938. [13] G. Radons, T. Geisel, J. Kell, Solid State Ionics 13 (75) (1984) 80 (North-Holland, Amsterdam). [14] A.C. Levi, R. Spadacini, G.E. Tommei, Surf. Sci. 121 (1982) 357. [15] R. Kubo, Rep. Prog. Phys. 29 (1966) 255. [16] J.P. Hansen, I.R. McDonald, Theory of Simple Liquids, Academic Press, London, 1986. [17] R. Ferrando, R. Spadacini, G.E. Tommei, G. Caratti, Physica, A 195 (1993) 506. [18] G. Caratti, R. Ferrando, R. Spadacini, G.E. Tommei, Phys. Rev., E 55 (1997) 4. [19] W. Dieterich, P. Fulde, I. Peschel, Adv. Phys. 29 (3) (1980) 527 – 605.
343
[20] M. Mazroui, A. Asaklil, Y. Boughaleb, J. Phys., I Fr. 7 (1997) 675 – 690. [21] H.P. Weber, H. Schulz, J. Chem. Phys. 85 (1986) 475. [22] K. Funke, in: A.L. Laskar, S. Chandra (Eds.), Superionic Solids and Solid Electrolytes, Academic Press, New York, 1989. [23] F. Montalenti, R. Ferrando, Phys. Rev., B 59 (1999) 8. [24] M. Aniya, Solid State Ionics 70 (1994) 673. [25] M. Aniya, Physica, B 239 (1997) 144. [26] H.A. Kramers, Physica 7 (1940) 284. [27] M. Mazroui, A. Asaklil, Y. Boughaleb, Surf. Sci. 409 (1998) 528 – 540; M. Mazroui, Y. Boughaleb, Int. J. Mod. Phys., B 15 (2000) 2193. [28] M. Mazroui, M. Hairibi, Y. Boughaleb, Phys. Lett., A 253 (1999) 227. [29] H. Risken, H.D. Vollmer, Phys. Lett. 6 (1979) 69. [30] A.R. Bishop, W. Dieterich, I. Peschel, Z. Phys., B Condens. Matter 33 (1979) 187.