Soliton diffusion on chains of coupled nonlinear oscillators

Available online at www.sciencedirect.com

Physica A 334 (2004) 417 – 443

www.elsevier.com/locate/physa

Soliton di usion on chains of coupled nonlinear oscillators Edward Ar$evaloa;∗ , Yuri Gaidideib , Franz G. Mertensa a Physikalisches

Institut, Universitat Bayreuth, Bayreuth D-95440, Germany Institute for Theoretical Physics, Kiev 03143, Ukraine

b Bogolyuobov

Received 30 April 2003; received in revised form 29 September 2003

Abstract We study the di usivity of solitons in noisy-coupled nonlinear oscillators. We consider an anharmonic atomic chain and a network of coupled RLC circuits with nonlinear capacitance and linear inductances and resistances. The solitons propagate in the presence of either thermal noise and/or shot noise. We use a Langevin formalism in order to model the noisy systems, for that reason noise and damping are added to the discrete equations of motion. In the long-wave approximation both systems can be described by identical noisy and damped Korteweg–de Vries (KdV) equations. By using the results of the well-known adiabatic perturbation theory of the forced KdV equation we derive ordinary di erential equations (ODEs) for the relevant variables of the soliton, namely position and inverse of the width. We solve these ODEs by using standard perturbation theory to obtain analytical expressions of the variance and average of the soliton position. We perform Langevin dynamics simulations of the full discrete systems which con>rm our analytical results, namely superdi usivity of the solitons depending on the initial velocity. c 2004 Elsevier B.V. All rights reserved.  PACS: 05.10.Gg; 05.40.−a; 05.45.Yv; 05.50.+q Keywords: Lattice solitons; Noise; Langevin dynamics simulations; Chains of oscillators

1. Introduction Many processes in nature can be described in terms of >nite chains of coupled nonlinear oscillators. Examples include the monatomic chain of particles describing polypeptide chains in muscle proteins [1–3] or the DNA molecule [4]. ∗

Corresponding author. E-mail address: [email protected] (E. Ar$evalo).

c 2004 Elsevier B.V. All rights reserved. 0378-4371/$ - see front matter  doi:10.1016/j.physa.2003.10.084

418

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

Another example is a chain of coupled nonlinear electrical oscillators used as equivalent electrical circuit of a transmission line or long Josephson junctions [5]. When the amplitude of an excitation in the system (displacement of particles or voltages and currents) is high, the nonlinear e ects cannot be ignored. The exact balance between dispersion and nonlinearity leads to the concept of solitary wave. The solitary waves are pulse-like waves which can propagate with constant velocity and pro>le. In fact, these solitary waves propagating on chains, for simplicity called here lattice solitons, have been used to clarify many features of molecular chains [5–10]. For example, due to their robust character lattice solitons have been used to explain the energy transport in polypeptide chains in muscle proteins [1–3] or the energy transport in DNA [4]. Numerical simulations at realistic temperatures for transport in proteins have shown that lattice solitons can propagate over long distances in a chain with the Lennard–Jones potential [2]. In fact, those solitons in the presence of Gaussian white noise present a superdi usive behavior and tend to vanish when the system is close to its stationary state [11]. On the other hand, laboratory experiments on chains of electrical oscillators have shown that those chains can bear lattice solitons whose shape can be predicted very well in the long wave approximation [5]. In fact, for a polynomial interaction potential between particles, namely harmonic term plus cubic or=and quartic anharmonicity, there are exact soliton solutions for the partial di erential equations following from the discrete systems. However, for more realistic interaction potentials like Lennard–Jones or Morse there are no exact soliton solutions. In a real system, noise is a fundamental stochastic variable which is always present. For instance, at the atomic scale, thermal noise of thermodynamical origin causes stochastic Luctuations of the positions of the particles, whereas in the case of an electrical circuit it causes stochastic Luctuations in the electron density within a conductor. In fact, Nyquist’s theorem [12] states that the mean square value of a Luctuating voltage in a circuit is proportional to the electrical resistance of a circuit and the temperature [13]. Thermal noise may be important in oscillators where electrical resistances appear as basic circuit elements. This type of noise in a circuit is independent of the material of the resistors and its spectrum [12] is Lat, i.e., a constant function of frequency (white noise). Another random noise signal in electrical circuits is the shot noise, current Luctuations due to the discreteness of the electron’s charge. This type of noise arises from the presence of active elements, namely vacuum tubes or solid-state devices, in electrical circuits. In this case Schottky’s theorem states that the mean square value of a Luctuating current in a circuit is proportional to the average electrical direct current (DC) through the active element times the electron’s charge. Schottky’s result corresponds to the uncorrelated arrival of particles (electrons or holes) with a Poisson distribution function of the time intervals between arrival tunes [14–16]. Notice that for a large average of carrier population the Poisson distribution changes into the Gaussian distribution. In the case of >nite temperatures the distribution function of the carrier population at thermal equilibrium must >t the Gibbs (or Maxwell–Boltzmann) form [15]. Shot noise in electrical circuits is independent of frequency from a few kHz to several hundred MHz. Thermal and shot noise are uncorrelated, so the total noise in a

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

R

R

R

L

L

C

na

R

Vn

Vn-1

Vn+1 L

C

419

C

L

(n+1) a

Fig. 1. Equivalent circuit of an electrical network with linear inductance L and non-linear capacitance C(Vn ) in each unit section. The series of resistances R represents the losses of the system and can be active and=or passive.

system is the square root of the sum of the squares of all the incoherent noise sources. Notice that experimentally the shot noise can be adjusted to be much greater than the thermal noise [13]. The aim of the present work is to study in the framework of the Langevin formalism the di usion of lattice solitons on anharmonic chains of coupled oscillators. We consider the presence of intrinsic random perturbations, namely thermal and=or shot noise, depending on the system. These sources of noise are modeled as Gaussian white noise. We study numerically and analytically two systems, namely a monatomic chain with nearest-neighbor interaction and the electrical network shown Fig. 1. Notice that the con>guration shown in Fig. 1 is a simple equivalent electrical circuit of a transmission line. The soliton di usion on monatomic chains in contact with a thermal bath was already studied [11] in the case when the thermal Luctuations are very large with respect to the soliton amplitude. In this regime of large Luctuations (large temperatures) and also low soliton energy (soliton velocities very close to the sound velocity) the shape of the soliton is strongly distorted, so the noise cannot be considered as perturbation, i.e., perturbation theory fails. Here, we study the soliton di usion in the regime where the noise sources can be considered as perturbations, namely higher energy solitons. In the present work we consider that the interaction potential between particles in the monatomic chain is harmonic plus a cubic anharmonicity. The coupling to the thermal bath is provided by adding an additive thermal noise to the discrete equations of motion describing the chain. The energy input of the thermal noise is balanced by a damping term providing energy dissipation. We consider a dissipation which is due to irreversible processes arising from the >nite velocity of the internal motions taking place within the system. This type of dissipation can be modeled by the so-called inner or hydrodynamical damping [17,18,11], which is extensively used in elasticity theory. The corresponding noise term, which ful>lls the Luctuation dissipation theorem, takes the form of a discrete gradient of Gaussian white noise delta-correlated in space and time [11]. On the other hand, the electrical network, which is realized with identical components, has a large number N of identical unit sections. Here we restrict ourselves to the case where the inductance is independent of the current but the capacitance depends

420

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

on the voltage. The resistances in the circuit can be active or can represent internal losses of the inductances. Since the Nyquist and Schottky theorems state the relation between the resistance in a circuit and the random Luctuations of current or voltage, we introduce a random term in Kircho ’s equations of the circuit representing these Luctuations. This noise term ful>lls the Luctuation dissipation theorem and has the form of a discrete gradient of Gaussian white noise delta-correlated in space and time. Notice that in many cases the results following from Boltzmann–Langevin theories of shot noise in mesoscopic conductors turn out to be identical to exact the quantum results averaged over an ensemble. This makes Boltzmann–Langevin theories of shot noise credible even in those situations where quantum results are not available [16]. The Langevin approach has been used for a long time to model the dynamics of solitons coupled to a noisy environment [19]. For instance, these types of models have been used in the problem of nucleation of topological kinks [20–22], or in the problem of di usion of single solitons [11,23,24]. In the next section we present the noisy and damped KdV equation as an approximation of the discrete equations of motion of the anharmonic chain and the electrical network of Fig. 1. Afterwards we introduce the well-known results of the adiabatic perturbation theory of the perturbed KdV equation. From this initial point we derive stochastic ordinary di erential equations (ODEs) of the relevant soliton variables. By using perturbation theory we >nd and approximate the analytical solution of the ODEs. We compare these analytical results with the results following from Langevin dynamics simulations of the full discrete systems. Our conclusions are summarized in the last section. 2. The long-wave approximation In the long-wave approximation the discrete equations of motion of both systems considered here lead to noisy and damped Bq equations. These equations possess bi-directional wave solutions [5]. Since in the present work we consider a simple solitary wave traveling in one direction, we can reduce the noisy and damped Bq equation in order to get an unidirectional wave. In this respect one can use a reductive perturbation technique [25,26,5], so >nally one gets the most simpli>ed nonlinear partial di erential equation modeling these anharmonic chains, namely a noisy and damped KdV equation. 2.1. The monatomic chain Here we follow Ref. [11]. In the monatomic chain of particles with mass M we consider that the interaction between nearest-neighbors is harmonic plus a cubic anharmonicity. The Hamiltonian of this system reads  
 P2 A 1 n 2 3 H= : (1) +G (Yn+1 − Yn ) + (Yn+1 − Yn ) 2M 2 3 n

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

421

where Yn denotes the longitudinal displacement of the nth particle from its equilibrium position, and dYn Pn = M (2) dt is the momentum. Here G and A are the potential parameters. The discrete equations of motion of the system in relative displacements read d 2 rn 2 2 − 2rn2 + rn−1 ) M 2 = G(rn+1 − 2rn + rn−1 ) + GA(rn+1 dt   drn drn−1 drn+1 −2 + + M dt dt dt √ + D(n+1 (t) − 2n (t) + n−1 (t)) ; (3) where rn (t) = Yn+1 (t) − Yn (t). In Eq. (3) we have already added both stochastic and damping forces characterized by damping () and di usion (D) constants, respectively. Here, D = 2MkB T ;

(4)

where kB is the Boltzmann constant, and T is the temperature of the thermal bath. n (t) is the delta-correlated white noise, n (t)m (t  ) = nm (t − t  ) ;

(5)

n (t) = 0 :

(6)

The noise and damping terms ful>ll the Luctuation–dissipation theorem. If the relative displacement rn varies slowly from one bond of the chain to the next one, the continuum approximation can be applied [7,5], i.e., we expand rn±1 (t) and rn±2 (t) in a Taylor series around r(x; t) with x = na, where the equilibrium atomic spacing a is regarded as an expansion parameter. Substituting this expansion in Eq. (3) yields a noisy and damped Bq equation, √ (7) 92t r = c2 92x r + p9x (r9x r) + h94x r + a2 92x 9t r + a2 D92x (x; t) ; whose leading term is of order a4 . Here n+1 (t) − 2n (t) + n−1 (t) → 92x (x; t) ; a5=2 with the properties (x; t) = 0;

(x; t)(x ; t  ) = (x − x )(t − t  ) :

(8) (9)

Besides, D=

2kB T ; 

h=

a3 G ; 12

c2 =

Ga ; 

 = M=a :

p=

2a2 AG ;  (10)

422

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

By using a reductive perturbation technique [25,26,5] the noisy and damped Bq equation (7) can be reduced to a noisy and damped KdV equation (A.1), namely √ 9 u + 6u9s u + 93s u = 1 92s u − D1 9s (s; ) ; (11) where r(x; t) → u(s; ) and (s; ) is white noise delta-correlated in space (s) and time ( ). The damping (1 ) and di usion (D1 ) constants are de>ned in Eqs. (A.8). Notice that Eq. (11) is an approximation of Eq. (7) in a reference frame moving with the sound velocity c. 2.2. The electrical network The elementary network in Fig. 1 has linear inductors L and nonlinear capacitors C. For example, if the nonlinear capacitance is that of a reverse-biased pn diode junction the capacitance–voltage relation can be approximated by a polynomial expansion [5]. For a small enough voltage, one can keep the >rst two terms of the polynomial expansion and write C(Vn ) = C0 − 2C1 Vn ;

(12)

where Vn (t) is the voltage across the nth capacitor and C0 and C1 are constants. From Kircho ’s law one obtains   d 1 dVn C(Vn ) = (Vn−1 − 2Vn + Vn+1 ) dt dt L   1 dVn−1 dVn dVn+1 : (13) + −2 + R dt dt dt By inserting Eq. (12) into Eq. (13) we get C0

d 2 Vn 1 1 = (Vn−1 − 2Vn + Vn+1 ) + 2 dt L R + C1



dVn−1 dVn dVn+1 −2 + dt dt dt

d 2 Vn2 √ + D(n+1 (t) − 2n (t) + n−1 (t)) ; dt 2



(14)

where we have already added the random functions n which are delta-correlated in time and space. In the case of thermal noise the di usion constant takes the form [12] D=

2kB T : RL2

(15)

In the case of shot noise, i.e., the resistances (R) are from active elements, the Schottky formalism states that the power spectral density of current noise is proportional to a DC current traveling through the active elements. So in the case of the circuit in Fig. 1 the shot noise may be neglected since most of the time the electrical current traveling through R is zero. However, an alternative circuit can be considered here. In fact, one can introduce a noisy DC voltage source Vndc in each unit section as shown in Fig. 2.

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

R −

dc Vn-2

L

R −

+ Vn-1 C

dc Vn-1

L

R −

+ Vn

L

C

Vndc

423

R −

+ Vn+1 C

dc Vn+1

+

L

Fig. 2. Circuit with a voltage source, Vndc , in each unit section.

This noisy voltage may be de>ned as where √ Vndc = V dc + D(Wn (t) − Wn+1 (t)) ;

(16)

where V dc = Vndc 

and

Wn (t) = 0

n = 0; 1; : : : :

(17)

In Eq. (16) the Wn (t) are independent Wiener processes [12], i.e., dWn (t) = n (t) ; dt

(18)

where n (t) is Gaussian white noise. The equation of motion of each unit section of the network, Fig. (2), is identical to Eq. (14). The de>nition of the di usion constant in the case of shot noise is D=

eI dc ; L2

(19)

where I dc = V dc =R is a constant current traveling through R and e is the electron’s charge. The constant current I dc can be a function of temperature T , i.e., I dc = I dc (T ) [15,16]. It is interesting to notice that, for instance, in the case of a linear conductor R the shot noise magnitude converges to the thermal noise magnitude as the charge quantum tends to zero [15], namely lim eV dc (T ) = 2kB T :

e→0

(20)

De>nition (19) is intuitive rather than a result of formal derivation. However, a similar consideration, using the Langevin formalism, has been used successfully to describe shot noise in mesoscopic conductors [16]. In Eq. (14) the noise and dissipation terms satisfy the Luctuation–dissipation theorem. Notice that the nonlinear terms of Eqs. (14) and (3) are di erent from each other. By using the continuum approximation Eq. (14) yields a noisy and damped Bq equation, √ 92t v = c2 92x v + h94x v + p9t (v9t v) + a2 92x 9t v + a2 D92x (x; t) ; (21)

424

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

which is di erent from Eq. (7). Here, Vn (t) → v(x; t) ;

(22)

 2akB T     RL2 C 2 0

D D= 2 = C0  aeI dc    2 2 L C0

thermal noise; (23) shot noise

and =

1 ; RC0

c2 =

a2 ; LC0

p=

2C1 ; C0

h=

a4 : 12LC0

(24)

By using the reductive perturbation technique Eq. (21) can be reduced to a noisy and damped KdV equation identical to Eq. (11) (A.2), where v(x; t) → u(s; ) and (x; t) → (s; ). Thus both the monatomic chain and the electrical circuit can be approached by the same perturbed KdV equation. 3. The perturbed KdV equation This section is devoted to solve the noisy and damped KdV equation by using perturbation theory. We write Eq. (11) in the form 9r u + 6u9s u + 93s u = (F ; where F = 1 92s u −

√

D1 9s (s; ) :

(25)

(26)

In Eq. (25) we have introduced a small parameter ( for convenience. In the case ( = 0 we get the standard KdV equation, 9r u0 + 6u0 9s u0 + 93s u0 = 0 ;

(27)

whose one-soliton solution reads u0 (s; ) = 2*20 sech2 (*0 (s − 4*20 )) :

(28)

Notice that *0 =

1

3c(v0 − c) p

(29)

is the inverse of the soliton width written in terms of the soliton velocity v0 in the rest frame. The sound velocity c and the constant p are de>ned in Eqs. (10) and (24).

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

425

So our problem is reduced to study the evolution of a single KdV soliton a ected by a force F which is a function of space and time. In order to solve approximately Eq. (25) we express the solution in a perturbation series of the form u(s; ) = u0 (s; ) + (u1 (s; ) + · · · ;

(30)

u0 (s; ) = 2*2 ( )sech2 (z);

(31)

where z = *( )(s − S( )) :

Here (1. By using either inverse scattering transform [27], perturbation theory [28] or a multiple-scale perturbation expansion [29] one can derive the following equations for the variables *( ) and S( ):  ∞ d*( ) ( d z sech2 (z)F ; (32) = d 4* −∞  ∞ dS( ) ( 2 d z(z sech2 (z) + tanh(z) + tanh2 (z))F : (33) = 4* + 3 d 4* −∞ Eqs. (32) and (33) take into account the force-induced phonons (continuum eigenfunctions). From now on we set ( to the unity and we suppose that F is a perturbative force. By inserting Eq. (26) into Eq. (32) we get √  ∞ d*( ) D1 81 *3 ds 9s (sech2 (z))(s; ) ; (34) =− + d 15 4 −∞ √  ∞ D1 81 * dS( ) ds 9s = 4*2 + + d 15 4*2 −∞ × (z sech2 (z) + tanh(z) + tanh2 (z))(s; ) :

(35)

Eqs. (34) and (35) contain multiplicative noise which is interpreted in the Stratonovich sense. This means that we assume (s; ) to be a real noise with >nite correlation time, which is then allowed to become in>nitesimally small after calculating measurable quantities [30,31,12]. Notice that white noise means taking the limit of zero correlation time. The Fokker Planck equation (FPE) associated with Eqs. (34) and (35) in the Stratonovich sense reads   81 81 * D1 D1 9 = 9S  − 9S  9* (*3 ) − 4*2 + 9*  − 15 15 10 20*2   D1 D1 2  D1 + + (42 + -2 )92S  : (36) 9 (*) − 9* 9S 30 * 15 * 360*3 At this point, we are interested in getting a set of Langevin equations with uncorrelated noises and associated FPE (36). For that reason we de>ne a new set of Langevin equations, ˜ dY ˜Y ˜ ; ] + B[ ˜ ; ]˜( ) ; ˆY = A[ (37) d

426

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

which we have interpreted in the Stratonovich sense. The elements of both the vector ˜ and the matrix Bˆ have to be determined. Here, {Y1 ; Y2 } = {S; *} are the elements of A ˜ , and ˜ is a vector of uncorrelated noises i which does not depend on the vector Y the position variable s, so i ( )j (  ) = ( −  )ij ;

i; j = 1; 2 :

(38)

The corresponding FPE for this set (37) [32] is
1
9Yi (Ai ) + 9Yi (Bij 9Yk (Bkj )); 9 =− 2 i

i; j; k = 1; 2 :

(39)

ijk

If we compare Eq. (39) with Eq. (36), we can derive >ve algebraic equations for the terms Ai and Bij , whose solutions read

−81 *3 =15 + D1 =12 ˜ (40) A[Y; ] = 4*2 + 81 *=15 and

 



D1 *  15  ˆ  B[Y; ]=  D1  − 15*3

0



   : 2  D1 (30 + - )  180*3

(41)

Then the Langevin equations (37) read  d*( ) D1 D1 * 8*3 1 = − + 1 ( ) ; d 12 15 15 8*1 dS( ) = 4*2 + − d 15



(42)

D1 1 ( ) + 15*3



D1 (30 + -2 ) 2 ( ) : 180*3

(43)

Eqs. (42) and (43) ate statistically equivalent to Eqs. (34) and (35), because both share the same FPE (36). By integrating once with respect to time, Eqs. (42) and (43) take the form      D1 D1 *(  ) 8*3 (  )1  d +S *( ) = *(0) + (44) − dW1 (  ) ; 12 15 15 0 0   S( ) = S(0) +

0

  +S

0

8*(  )1 4* ( ) + 15 2





D1 (30 + -2 ) dW2 (  ) ; 180*3 (  )



d −S

  0

D1 dW1 (  ) 15*3 (  ) (45)

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

where the Stratonovich integrals read    1
S Bil dWi (  ) = Bil dWi (  ) + Bkj 9Yk Bij d 2 0 0 0



427

:

(46)

jk

Here Bij are the matrix elements of Eq. (41) and dWj ( )=j ( )d are Wiener processes where j ( )j
(  ) = ( −  )jj
:

(47)

Notice that (46) gives the connection between the Ito and Stratonovich representations [12]. Inserting Eq. (46) in Eqs. (44) and (45) yields      D1 D1 *(  ) 8*3 (  )1  d + (48) *( ) = *(0) + − dW1 (  ) 10 15 15 0 0  D1 8*(  )1 d  + 4* ( ) + S( ) = S(0) + 15 20*2 (  ) 0     D1 D1 (30 + -2 )  − ( ) + dW dW2 (  ) : 1 15*3 (  ) 180*3 (  ) 0 0  

2



(49)

In Appendix B we solve Eqs. (48) and (49) by using perturbation analysis. In fact, we consider the thermal terms as perturbations, so Eqs. (48) and (49) take the form   3   8* ( )1 d  *( ) = *(0) − 15 0 

  D1  D1 *(  )  +( (50) d + dW1 ( ) 15 0 10 0   S( ) = S(0) +

0





+(

0

  +

0

4*2 (  ) +

8*(  )1 15  

D1 d − 20*2 (  ) 2

0

 d



D1 dW1 (  ): 15*3 (  )



D1 (30 + - ) dW2 (  ) : 180*3 (  )

(51)

Now, we seek an asymptotic solution in the form of a small-noise expansion S( ) = s0 ( ) + (s1 ( ) + · · · ; *( ) = *0 ( ) + (*1 ( ) + · · · :

(52)

428

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

Here, ( is a factor introduced for convenience in the analytical calculations. Notice that the case ( = 0 reduces Eqs. (50) and (51) to the deterministic case. In order to interpret our perturbation theory we must set ( = 1 and assume that the terms on the r.h.s. of Eqs. (50) and (51) are suPciently small. So we must restrict ourselves to a regime of low temperatures of the thermal bath (D1 small). From the perturbation analysis (see Appendix B) we obtain the following analytical solutions for the average position S( ) and the variance of the position Var(S( )): S( ) = s0 ( ) + s1 ( ) √ 15 2−1 + = log(2) *0 (0) 41     √ 3D1 20(3 − 52 + 225=2 )*0 (0) + 8 2 − 152 + 723 1 − 6402*40 (0)21

(53)

and Var(S( )) = Var(s1 ( ))  3 −75 225 152 2  = D1 − + 32*30 31 22422 *30 31 28*30 31 345 − − 224*40 21

15 3

562 2 *40 21 3

1

152 4 22922 + − 14*40 21 224*40 21 5

(218 + 7-2 )2 2 9*0 + 1 458 + 7-2 24 + + √ 5 2 − + 5 5 5 3360*0 1 3360*0 1 42*0 1 24 2*0 1  7 2 4 (45*0 + 1 ) 315*20 + 35*0 1 − 221  ; + + 42*50 21 1122*50 31

(54)

where 2=1+16*20 (0) 1 =15. Notice that in the absence of damping, 1 =0, Var(S( ))=0 and S( ) = 4*20 (0) , where 4*20 (0) is the initial velocity of the KdV soliton. Since formula (54) is not amenable to physical interpretation, we expand this in a Taylor series in about = 0, namely  Var(S( )) D1 +

(42 + -2 ) ((34 + -2 )1 − 120*0 (0)) + 225*0 (0) 180*30 (0)

2

8*0 (0)((38 + -2 )21 + 2401 *0 (0) + 1800*20 (0)) 10 125

3

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

16(*30 (0)1 (27 000*20 (0) + 3270*0 (0)1 + (112 + -2 )21 )) 151 875  + O( 5 ) :

+

429 4

(55)

The term linear in derives from the direct e ect of the noise on the position of the soliton. The higher order terms in follow from an indirect e ect; the noise changes the soliton width, *−1 ( ), stochastically. This variable determines the velocity of the soliton, 4*2 ( ), and therefore the soliton position is also a ected. This indirect e ect of the noise is nonlinear in because the solitons here are nontopological [5], i.e., ∞ the quantity −∞ u(s; ) ds is not conserved. In fact, this quantity depends on *( ) which is nonlinearly a ected by the damping, 1 . This nonlinear e ect of the damping is reLected in the higher order terms in of expansion (55). For low initial soliton velocities (*0 (0) → 0) the linear term in is dominant at early times ( ¡ 3=1 ), so the variance tends to be linear here. At later times ( ¿ 3=1 ) the higher order terms in are dominant, so a superdi usive behavior takes place. Notice that for high-initial soliton velocities, (*−1 0 (0) → 0) the higher order terms also become dominant even at early times ( ¡ 3=1 ). In Appendix B we compare the variance given by Eq. (54) with the numerical solution of Eqs. (42) and (43). Notice that Eq. (29) gives the relation between * and the soliton velocity v0 in the rest frame.

4. Simulations Our Langevin dynamics simulations were performed for chains with 1500 oscillators. The time integration was carried out by using the Heun method [33], which has been successfully used in the numerical solution of partial di erential equations and di erence–di erential equations, coupled to either an additive or a multiplicative noise term [34,24,23,11]. In order to start the simulations at t = 0 we have used the one-soliton solution of the respective homogeneous Bq equation [7,5]. The average values have been calculated over 200 realizations up to a >nal time 5000. All values of the constants of Eqs. (56) and (62) are set at unity except the damping constant which is set at =1=RC0 =0:003. Notice that for lower values of damping the relaxation of the system energy would take more time in our simulations to reach a regime close to its stationary value. On the other hand, higher values of damping can strongly distort the soliton shape, namely the tail induced by the damping cannot be neglected when the value of the damping is high. In Appendix C we show the energy relaxation process of the system up to times when the system energy is close to its stationary value. In our simulations the values of the di usion constant, D, and initial soliton velocity, v0 (0), are parameters (see >gure captions). In order to check the accuracy of our codes we have computed some conserved quantities in both systems which are mentioned below.

430

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

4.1. The monatomic chain For this system we have performed Langevin dynamics simulation of Eq. (3). For that we have written Eq. (3) as follows; dPn 2 2 = G(rn+1 − 2rn + rn−1 ) + GA(rn+1 − 2rn2 + rn−1 ) dt   √ dPn−1 dPn+1 dPn + + D(n+1 (t) − 2n (t) + n−1 (t)) ; + −2 dt dt dt d n Pn (56) = M dt (compare these equations with Eqs. (62) below). Notice that the lattice solitons in relative displacements have a pulse shape, i.e., the amplitude of the soliton vanishes at in>nity. This characteristic allows us to use periodic boundary conditions which are necessary for long simulation times, because we want to avoid reLections at the boundaries. The periodic boundary conditions read d l r0 d l rN −1 d l rN d l r1 = ; = l ; l = 0; 1; l l l dt dt dt dt 0 (t) = N −1 (t);

N (t) = 1 (t) ;

(57)

where N is the number of particles of our chain and N − 1 is the number of bonds. A suitable method to detect the position of a pulse lattice soliton, rn (t), is to search for its maximum [18]. However, in the presence of stochastic perturbations this method is not useful since the pulse shape is masked by the noise. So from the data of our simulations we have taken snapshots of the system at di erent times, and from them we have generated the kink shape Yn (t) of the lattice soliton by using algorithm [11] n−1
ri (t); n = 2; 3; : : : ; N : (58) Yn (t) = Y1 (t) + i=1

The kink shape is less distorted by the noise than the pulse shape ri (t). In (58) Y1 (t) is a boundary condition that we have demanded to be N −1 1
ri (t) ; (59) Y1 (t) = − 2 i=1

since the kink is symmetrical about the origin (odd function) [7,5]. Notice that N −1
YN (t) − Y1 (t) = ri (t) (60) i=1

is a conserved quantity in this system, i.e., N −1 dYN (t) dY1 (t)
dVi (t) − = =0 : dt dt dt

(61)

i=1

We have checked Eq. (61) with a precision higher than 10−14 over the whole time range of our Langevin dynamics simulations.

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

431

Afterwards we have searched, using linear interpolation, the position where this kink intersects the x-axis. Finally, we have de>ned this intersection point as the soliton position. 4.2. Electrical network In this case we rewrite Eq. (14) so C0

d6n 1 = (Vn−1 − 2Vn + Vn+1 ) dt L   6n−1 6n 6n+1 1 −2 + + R 1 − 2(C1 =C0 )Vn−1 1 − 2(C1 =C0 )Vn 1 − 2(C1 =C0 )Vn+1 √ + D(n+1 (t) − 2n (t) + n−1 (t)) ; dVn 6n = 1 − 2(C1 =C0 )Vn dt

(62)

(compare these equations with Eqs. (56) above) under periodic boundary conditions. We note that quantity (60) is not conserved here. In fact, we have note in our simulations that N −1


(Vi (t) − Vi (0)) ¿ 0

(63)

i=1

for t ¿ 0. Therefore algorithm (58) for constructing the kink shape is not suitable here. On the other hand, the quantity  N −1 
C1 2 Vi (t) − QN (t) − Q1 (t) = V (t) (64) C0 i i=1

is conserved, i.e., N −1

dQN (t) dQ1 (t)
dVi (t) − = dt dt dt i=1



 C1 1−2 Vi (t) = 0 : C0

(65)

We have checked Eq. (65) in the deterministic case, namely D = 0, with a precision higher than 10−6 over the whole time range of our Langevin dynamics simulations. Here  n−1 
C1 2 Qn (t) = Q1 (t) + Vi (t) − Vi (t) (66) C0 i=1

is kink function. In Eq. (66) we have demanded  N −1  1
C1 2 Vi (t) − Vi (t) ; Q1 (t) = − 2 C0 i=1

(67)

432

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

70 60 C z (t)

50 40

B

30 20

A

10 1000

2000

3000

4000

5000

t Fig. 3. Averaged soliton position vs. time in the sound velocity frame z = x − t (c = 1), with  = 0:003 and T = 10−5 . Dotted and dash-dotted lines correspond to the Langevin dynamics simulations of both discrete systems, Eqs. (56) and (62), respectively. Solid lines correspond to the analytical solution, Eq. (53). A, B and C correspond to initial velocities v0 (0) = 1:007; 1:01 and 1:02, respectively. Notice that in some cases it is diPcult to distinguish between di erent lines.

so the kink function is odd at t = 0. Here, again we have taken snapshots of the system at di erent times, and from them we have generated the kink Qn (t). Afterwards we have de>ned the position of the lattice soliton for t ¿ 0 as the point where the kink intersects the x-axis. 4.3. Results In Fig. 3 we show a comparison of the soliton position obtained in both systems and the analytical prediction given by Eq. (53). Notice that Z(t) = x(t) − ct =

1 S(:t) ; 9

(68)

where S is de>ned in Eq. (53), the sound velocity is de>ned in Eqs. (10) and (24), and the constant 9 is de>ned in Eq. (A.4). The analytical prediction of the soliton position for low velocities, namely v0 =1:007; 1:01, >ts very well with the observed position of the soliton in the Langevin dynamics simulations. For higher velocities (v0 = 1:02) there is discrepancy between the soliton position observed in the electrical network and the analytical prediction. This discrepancy is not observed in the case of the monatomic chain. Notice that the nonlinear terms of the discrete equations of motion of both systems are di erent. Moreover, these nonlinearities are approached di erently in the reductive perturbation technique in order to get KdV-like equation (see Appendix A). In fact, we have made a change of variables, namely Eq. (A.3), where we have introduced a small coePcient in front of the time and position variables. The coePcient in front of the position variable, ;1=2 , turns out to be larger than the coePcient in front of the time variable, ;3=2 , so our KdV approximation is less sensible to time variations than to position variations of the >eld. Notice that the nonlinear term of the equations

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

0.06 Var (x)

Var (x)

0.05 0.04 0.03 0.02 0.01 0

0

1000

2000

(a)

3000

4000

5000

0.3 0.2

Var (x)

Var (x)

0.25 0.15 0.1 0.05 0

0

1000

2000

3000

4000

5000

t

0.4

Var (x)

Var (x)

0.5 0.3 0.2 0.1 0

0

1000

2000

3000 t

1000

2000

4000

5000 (f)

3000

4000

5000

3000

4000

5000

t 3.5 3 2.5 2 1.5 1 0.5 0

0

1000

2000

(d)

0.6

(e)

0

(b)

t

(c)

0.7 0.6 0.5 0.4 0.3 0.2 0.1 0

433

t 7 6 5 4 3 2 1 0

0

1000

2000

3000

4000

5000

t

Fig. 4. Variances of the soliton position: (a) E = 10−6 ; v0 (0) = 1:007, (b) E = 10−5 ; v0 (0) = 1:007, (c) E = 5 × 10−6 ; v0 (0) = 1:01, (d) E = 5 × 10−5 ; v0 (0) = 1:01, (e) E = 5 × 10−6 , v0 (0) = 1:02, (f) E = 5 × 10−5 ; v0 (0) = 1:02. (dashed lines) Eq. (54), (solid and dotted lines) Langevin dynamics simulations of Eqs. (56) and (62), respectively.

of motion of the monatomic chain is a ected by discrete spatial derivatives while in the equations of the electrical network it is a ected by a time derivative. In order to analyze our results we de>ne the Luctuation energy   thermal noise;  kB T (69) E = eV dc (T )   shot noise: 2 In the thermal noise case this term can be interpreted as temperature (E = T ) since in our simulations kB = 1 (see the values of E in the captions of Fig. 4).

434

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

In Fig. 4 we show the variance of the soliton position vs. time for di erent values of E and initial soliton velocity, v0 (0). We observe in all cases short-time linear behavior followed by a superdi usive behavior in both systems (continuous and dotted lines). The magnitude of the variance scales with the value of E. The behavior of the variances in both systems tends to be similar for low values of initial velocity (low-energy solitons: Figs. 4(a) and (b)) as well as for the intermediate velocities (Figs. 4(c) and (d)). This is because analytically both systems are similar for low amplitudes of the solitons (low-energy solitons). Indeed, Eq. (14) under the condition Vn (t)1 becomes C0

d 2 qn 1 C1 2 2 = (qn−1 − 2qn + qn+1 ) + − 2qn2 + qn+1 ) (q dt 2 L C0 L n−1   dqn dqn+1 1 dqn−1 −2 + + R dt dt dt √ + D(n+1 (t) − 2n (t) + n−1 (t)) ;

(70)

where qn = Vn −

C1 2 V ; C0 n

Vn qn +

C1 2 q ; C0 n

1 dVn2 (t)

0: R dt

(71)

Notice that Eq. (70) is identical to the equation of motion of the monatomic chain, Eq. (3). For higher velocities (high-energy solitons: Figs. 4(e) and (f)) the variances in the electrical network (dotted lines) are higher than in the monatomic chain (continuous line). This di erence, as we have mentioned above, is mainly due to the fact that the nonlinearities in both systems are di erent. This di erence becomes more evident for higher energy solitons. Another factor which makes a di erence between the results of both systems are the rounding errors of the simulations. In fact, the precision of the simulations in the electrical network is lower than in the monatomic chain as we have mentioned above (see comments below Eqs. (61) and (65)). Notice that in Eq. (62) divisions are involved but not in Eqs. (56). With respect to the analytical results (dashed lines), they predict qualitatively well the behavior of the numerical results (continuous and dotted lines). We observe for low-energy solitons that the variance is linear in time for times t ¿ 3= = 1000. At later times (t ¿ 1000) a superdi usive behavior is observed. Notice that for higher energy solitons (Figs. 4(e) and (f)) the behavior remains similar except that the superdi usive behavior becomes important even for times t ¡ 1000. This whole behavior agrees with our analysis of Eq. (55). The discrepancies between theory and simulations in the case of low-energy solitons (Figs. 4(a) and (b)) are mainly due to numerical artefacts generated by the method of detection of the lattice soliton in the chain. We have de>ned the position of the soliton as the point where the center of its kink shape is localized, as we mentioned in Sections 4.1 and 4.2. In the case of low-energy solitons, the width is broader, i.e., there are several site points close to the center of the kink. In this case with the presence of noise it happens often that there are more than one point which satis>es our de>nition of the position of the kink for a given time. To avoid this we interpolate the 10 near-neighbor site points to the center of the kink.

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

435

From this interpolation we make a >rst estimation where the center of the kink lies. Afterwards, we choose from the points which satisfy the de>nition of soliton position, the nearest one to the point estimated from the interpolation. This chosen point is called the soliton position. This criterion introduces an error in the detection of the soliton position which can be of the order of the lattice constant (a = 1) in the case of low-energy solitons. Unfortunately there is no way to overcome this numerical artefact, since we are working in a discrete system.

On the other hand, notice that for low-energy solitons Var(x) ¡ 1 (Figs. 4(a) and (b)). This means that we observe variations of the soliton position which are even smaller than the length of the lattice constant. So, in this regime the discrepancy of the numerical variance and the theoretical prediction can be relatively large, as we observe in Figs. 4(a) and (b). The analytical approximation is better for higher energy solitons (see Figs. 4(c) and (d)) where the random Luctuations of the soliton shape are comparatively smaller with respect to the soliton amplitude. Moreover, the soliton are more narrow, i.e., there are less site points close to the center of the soliton, so the soliton position gets less ambiguous, i.e., we do not meet more than one point satisfying the de>nition of the soliton position. That is the reason why our analytical approximation predicts better the di usion of high-energy solitons in the monatomic chain. Notice that we have mentioned above the reason why the KdV equation (11) approaches better the monatomic chain than the chain of electrical oscillators. It is also clear that our analytical result Eq. (54) itself contains a truncation error. However, the relative error of this truncation is lower than 10%, as we show in Appendix B. Finally, in order to compare the results of the monatomic chain presented here with the results from the collective variable approach presented in Ref. [11] we de>ne the reduced temperature TS = kB T=H (0), where H (0) is the initial soliton energy in the monatomic chain. The values of the reduced temperature for the cases presented in the panels of Fig. 4 are: TS = 3:39 × 10−4 (a), 3:39 × 10−3 (b), 9:83 × 10−4 (c), 9:83 × 10−3 (d), 3:34 × 10−4 (e) and 3:34 × 10−3 (f). These values are one order of magnitude lower than those presented in Ref. [11]. This is the reason why perturbation theory predicts well the results here. Notice that the theory here works well for higher soliton velocities, while the collective variable approach presented in Ref. [11] is valid only for velocities very close to the sound velocity (soliton velocities 6 1:007). 5. Summary and conclusions We have studied the di usion of lattice solitons on two systems of nonlinear coupled oscillators, namely the monatomic chain of atoms and a network of coupled electrical oscillators shown in Fig. 1. We have considered the cases where the solitons propagate in the presence of thermal and/or shot noise (in the electrical network). In the monatomic chain the interaction potential between the atoms is harmonic plus a cubic anharmonicity. The chain is coupled to a thermal bath at a given temperature. In the electrical network we have considered a chain of coupled RLC circuits where the

436

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

resistance R and the inductance L are linear. The capacitance C depends nonlinearly on the voltage. We have used a Langevin formalism, i.e., we have added noise and damping terms to the discrete equations of motion in both systems. The noise and the damping terms satisfy the Luctuation–dissipation theorem. The systems exhibit di erent nonlinear terms. In the long-wave approximation both systems can be approximated by noisy and damped Boussinesq (Bq) equations. By using a reductive perturbation technique these equations can be approximated by the same noisy and damped KdV equation. We have used the results of the well-known adiabatic perturbation theory of the forced KdV equation. We have derived ordinary di erential equations (ODEs) for the relevant variables of the soliton, namely position and inverse of the width. We have solved these ODEs by using standard perturbation theory to obtain analytical expressions for the variance and average of the soliton position. We have also performed Langevin dynamics simulations of the full discrete systems. We have observed the di usion of the lattice solitons when the system energy is far from and close to its stationary state. We have calculated numerically and analytically the variance of the position of the lattice soliton in the presence of noise and damping. We have observed that the variances scale with the value of the Luctuation energy (temperature in the case of the monatomic chain). There is a short-time linear behavior of the di usion followed by a superdi usive behavior. Our analytical results predict qualitatively well this behavior for di erent values of the soliton velocity and Luctuation energy. We have found that the behavior of the soliton di usion depends on the soliton velocity (soliton energy). For low-energy solitons a short-time linear behavior of the variance is observed followed by a superdi usive behavior. For higher energy solitons the superdi usive behavior becomes more dominant even for short time scales. The superdi usive behavior is always dominant for long time scales. The short-time linear behavior of the di usion is due to the direct e ect of the noise in the soliton position. The superdi usive behavior is due to the fact that the noise a ects the width of the soliton, and therefore its velocity and position. This e ect is nonlinear because the solitons here are nontopological. This means that the solitons here lose robustness against the stochastic Luctuations with the time due to the e ect of the damping on the velocity. We have found a discrepancy between our theory and the results from simulations for low-energy solitons. This discrepancy is due to artefacts associated with the discreteness of the systems. In fact, we have realized that the resolution of the detection of lattice solitons is limited by the lattice constant. So the numerical error in the detection of the soliton position can be of the order of the lattice constant in the case of low-energy solitons. Unfortunately, the standard deviation of the soliton position, due to the noise e ects, can also be of the order of the lattice constant for low-energy solitons, therefore the discrepancies between the theory and the simulations here are important. For higher energy solitons the theory describes well the behavior of the soliton di usion in both systems. In fact, our analytical results >t very well the soliton di usion in the monatomic chain. However, there is a discrepancy between the soliton di usion of the monatomic chain and the electrical network. This is due to the nonlinear terms which are di erent in both systems. On the other hand, the KdV approximation used in our analytics describes better the monatomic chain than the electrical network.

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

437

We do not observe in our numerical results any dependence on the size of the system. In general, our above results show that both systems, the monatomic chain and the electrical network, behave similarly for low-energy solitons and di er for high-energy solitons. The KdV approximation can predict well the soliton di usion in both systems, however, this is a better approximation for the monatomic chain. Acknowledgements Yuri Gaididei would like to acknowledge the hospitality of the University of Bayreuth where this work was performed. Appendix A. From the Bq equation to the KdV equation A.1. The Bq equation from the monatomic chain The noisy and damped Bq equation (7) follows from the continuum approximation of the monatomic chain. In order to reduce this equation to a noisy and damped KdV equation we use a reductive perturbation technique [25,26,5]. So we rewrite Eq. (7) in a perturbational form as follows: √ 92t r − c2 92x r − p9x (r9x r) − h94x r = ;1=2 (a2 92x 9t r) + a2 D92x (x; t) ; (A.1) where we have introduced a small parameter ; for convenience. Afterwards, we perform the following change of variables: s = ;1=2 9(x − ct);

= ;3=2 :t;

u = >r ;

(A.2)

so that 9x = ;1=2 99s

and

9t = ;3=2 :9 − ;1=2 9c9s ;

(A.3)

where 9 = 1;

:=

h ; 2c

>=

p : 6h

(A.4)

We have also expressed u and  in a perturbation series u = ;u1 + ;2 u2 + · · · ;

(A.5)

 = ;1 + ;2 2 + · · · :

(A.6)

Here the parameter ; indicates the magnitude of the rate of change, the co-ePcients ; and ;3 in (A.2) are chosen in order to balance the nonlinear and dispersive terms. Notice that the dispersive term and the time derivative are of the same order in ;. The small noise expansion (A.6) is de>ned such that the lowest order terms of noise and

438

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

damping are of the same order. Substituting Eqs. (A.3), (A.5) and (A.6) into (A.1), and keeping the lowest order terms, namely O(;3 ), we get √ (A.7) 9s 9 u1 + 392s (u12 ) + 94s u1 = 1 93s u1 − D1 92s 1 (s; ) ; where both u1 and 9s 1 vanish at s → −∞. By integrating once Eq. (A.7) with respect to s, we >nally obtain Eq. (11) where we have set u = u1 and  = 1 . Here 1 =

a2 c h

D1 =

and

a4 p 2 D : 72h3 c

(A.8)

A.2. The Bq equation from the electrical circuit The noisy and damped Bq equation (21) follows from the continuum approximation of the chain of electrical oscillators. In order to reduce this equation to a noisy and damped KdV equation we follow the same procedure used in the Section A.1. So we rewrite Eq. (21) in a perturbational form √ (A.9) 92t v = c2 92x v + h94x v + p9t (v9t v) + a2 D92x (x; t) + ;1=2 (a2 92x 9t v) : Afterwards we make the same change of variables (A.3) where the coePcients 9; : and > are exactly the same as in Eq. (A.4) if one replaces > by → c2 >. By expressing u and  in a perturbation series (A.5), and keeping the lowest order terms, namely O(;3 ), we >nd again Eq. (11), where 1 and D1 are exactly the same as in Eq. (A.8) if one replaces D1 by c5 D1 . Appendix B. Perturbation analysis In this Appendix we develop a perturbation approach to the equations   3   8* ( )1 *( ) = *(0) − d  15 0 

  D1  D1 *(  )  +( d + dW1 ( ) 15 0 10 0   S( ) = S(0) +

0





+(  +

0

0



4*2 (  ) +

8*(  )1 15  

D1 d − 20*2 (  ) 2

0

 d

(B.1)



D1 dW1 (  ) 15*3 (  )



D1 (30 + - ) dW2 (  ) : 180*3 (  )

(B.2)

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

439

We seek an asymptotic solution of the form S( ) = s0 ( ) + (s1 ( ) + · · · ; *( ) = *0 ( ) + (*1 ( ) + · · · :

(B.3)

Inserting Eqs. (B.3) into Eqs. (B.1) and (B.2) and collecting powers of ( we get (0 :   3   8*0 ( )1 *0 ( ) = *0 (0) − d  15 0    8*0 (  )1 2  d : 4*0 ( ) + (B.4) s0 ( ) = s0 (0) + 15 0 (1 :   *1 ( ) = *1 (0) +

0



D1 81 *20 (  ) − *1 (  ) 10 5



 

d +

0

D1 *0 (  ) dW1 (  ) 15 (B.5)

 

 D1 8*1 (  )1   + s1 ( ) = s1 (0) + + 8*0 ( )*1 ( ) d  15 20*20 (  ) 0     D1 D1 (30 + -2 )  − dW1 ( ) + dW2 (  ) : 3  3( ) 15* ( ) 180* 0 0

(B.6)

Solving Eqs. (B.4) we obtain  1 161 *0 ( ) = + 15 *20 (0)  s0 ( ) =

15 1 161 1 + log + − 2 15 *0 (0) 41 *0 (0)

(B.7)  1+

161 *20 (0) 15

 ;

(B.8)

where s0 (0) = 0. Inserting Eq. (B.7) in (B.5) and solving it, with the initial condition *1 (0) = 0, we >nd √ D1 (−225 15 + (15 + 16 *20 (0)1 )5=2 ) *1 ( ) = 400*20 (0)1 (15 + 16 *20 (0)1 )3=2

 D1 *0 (0) + 1=4 (15 + 16  *20 (0)1 )5=4 dW1 (  ) : (B.9) 15 (15 + 16 *20 (0)1 )3=2 0

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

2

0.2

1.5

0.15

Var (x)

Var (x)

440

1

0.1 0.05

0.5

0

0

0

1000

2000

4000

5000

t

0

2000

0.4

3

0.3

2 1

3000

4000

5000

3000

4000

5000

t

4

0

1000

(b)

Var (x)

Var (x)

(a)

3000

0.2 0.1

0

1000

2000

(c)

3000

4000

t

0

5000 (d)

0

1000

2000 t

Fig. 5. (a) E = 5 × 10−5 ; v0 (0) = 1:01, (b) E = 5 × 10−6 ; v0 (0) = 1:01, (c) E = 5 × 10−5 ; v0 (0) = 1:02, (d) E = 5 × 10−6 ; v0 (0) = 1:02. The parameter E is de>ned in Eq. (69). (Solid lines), Eq. (54). (Dashed lines) numerical solution of Eqs. (42) and (43).

Then inserting Eqs. (B.7) and (B.9) in Eq. (B.5) and integrating once we get √ 3D1 (20(3 − 52 + 225=2 )*0 (0) + (8 2 − 152 + 723 )1 ) s1 ( ) = − 6402*40 (0)21     D1 D1 (30 + -2 )  − dW1 ( ) + dW2 (  ) 3  3( ) ( ) 15* 180* 0 0  √   8 15*0 (0) + 21 
+ 153=4 d  25=4 dW1 (  ) ; (B.10) 22522 0 0 where 2 = 1 + 16*20 (0) 1 =15. To this order we have that the variance of the position Var(S( )) = (2 Var(s1 ( )) and thus we >nally obtain Eq. (54). Notice that the variance in the rest frame reads Var(x(t)) =

1 Var(S(:t)) : 92

(B.11)

In Fig. 5 we show some examples of Var(x) given by Eq. (B.11) where Var(S( )) is given by Eq. (54). We compare this analytical result with results from a numerical solution of Eqs. (42) and (43) for which we have used the Heun method [33]. The numerical variances have been obtained by averaging over 1000 realizations. The

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

441

relative error |(Var(x))numeric − (Var(x(t)))analytic | . 0:1 (Var(x(t)))numeric

(B.12)

for the cases shown in Fig. 5. Here, (Var(x(t)))analytic follows from Eq. (B.11) together with Eq. (54), and (Var(x(t)))numeric follows from the numerical solution of Eqs. (42) and (43). Appendix C. Energy relaxation process In this Appendix we show the energy relaxation process of the electrical network. In the case of the monatomic chain in contact with a thermal bath this process has been discussed already in Ref. [11]. In fact, the results for the electrical network shown here are similar to those obtained in the monatomic chain. The energy function of the circuit in Fig. 1 in the absence of losses (1=R = 0) reads 
L C0 2 2C1 3 Vn − Vn ; (C.1) h= In2 + 2 2 3 n where In is the current through the nth inductance L and Vn is the voltage across the nth capacitor. For suPciently small values of the voltage Vn , namely        3C0   C0 2   2C1 3      Vn    (C.2)   3 Vn  i:e:  4C1  |Vn | ;  2 the energy function (C.1) can be approximated by its harmonic limit. Notice that C0 = C1 = 1 in our simulations. Condition (C.2) is always satis>ed in the present work. In the harmonic limit the Hamiltonian of the circuit can be written as follows:
P2 (@n − @n−1 )2 n ; (C.3) H= + 2L 2C0 n where @n = C0

n
j=1

Vj

and

Pn = L

d@n dt

(C.4)

are canonical variables of the system. Notice that In =I0 −(d@n =dt), where I0 is constant which we have set at zero (I0 = 0). The generalized equipartition theorem states that [35]    NkB T  N thermal noise; 
9H (t) = NeV dc @n (t) (C.5)  9@n (t)  shot noise; n=1 2 where T is the temperature of an external bath and V dc is de>ned in Eq. (17). For >nite values of C1 relation (C.5) is a very good approximation to evaluate the temperature

442

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

1

1Heqp2 / NE

0.8 0.6 0.4 0.2 1000

2000

3000

4000

5000

t fl

Fig. 6. Heqp =NE vs. time t with E = 5 × 10−5 ;  = 0:003, solid line: the electrical network, dashed line: the monatomic chain.

of the system when the voltages Vn satisfy condition (C.2). So, in order to examine the relaxation process one can de>ne the ensemble average   N
9H (t) : (C.6) @n (t) Heqp (t) = 9@n (t) n=1

For >nite times Heqp (t) possesses two contributions, one due to the random Luctuafl sol tions (shot or thermal noise), Heqp (t), and the other one due to the soliton, Heqp (t). So fl sol Heqp (t) = Heqp (t) + Heqp (t) :

(C.7)

Notice that for a very large system (N 1500) the contribution of the soliton energy, sol (t), can be neglected. However, in our system of 1500 unit sections (particles) Heqp this contribution is appreciable. In fact, the time evolution of Heqp (t) for di erent initial soliton energies presents di erent values due to the soliton contribution. Here we remark that in saturation (thermal equilibrium), i.e., t → ∞, these di erences fl vanish. In Eq. (C.7), Heqp (t) gives us information about the time evolution of the energy (temperature) of the system in terms of the saturation energy (temperature of sol the thermal bath). The soliton contribution Heqp (t) can be evaluated numerically using Eq. (C.6) when the soliton propagates in the presence of the damping but without noise. fl sol Then one can perform the numerical subtraction Heqp (t)−Heqp (t) to obtain Heqp (t), which is shown in a normalized form in Fig. 6. Note the di erence between both processes of relaxation following from the electrical network and from the monatomic chain. This quantitative di erence is due to the e ect of the nonlinearities of both systems. In fact, for lower values of E both relaxation processes are identical. We remark that we observe the same result in systems with the double number of sites, namely 3000. The normalized relaxation process depends only on the damping constant which has the same value  = 1=RC0 = 0:003 in all our simulations. Notice that for times t . 3= = 1000 there is a fast relaxation process, but for larger times, t & 3=, the system energy approaches its stationary value very slowly. The saturation (thermal

E. Ar-evalo et al. / Physica A 334 (2004) 417 – 443

443

fl equilibrium) of the system corresponds to the case Heqp (t)=NE = 1. The parameter E is de>ned in Eq. (69).

References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35]

S. Yomosa, Phys. Rev. A 32 (1985) 1752. P. Perez, N. Theodorakopoulos, Phys. Lett. A 117 (1986) 405. P. Perez, N. Theodorakopoulos, Phys. Lett. A 124 (1987) 267. V. Muto, P.S. Lomdahl, P.L. Christiansen, Phys. Rev. A 42 (1990) 7452. M. Remoissenet, Waves called Solitons, Springer, Berlin, Heidelberg, 1996. M. Toda, Theory of Nonlinear Lattices, Springer, Berlin, 1981. Sp. Pnevmatikos (Ed.), Singularities and Dynamical Systems, Mathematical Studies, Vol. 103, North-Holland, Amsterdam, 1985. A.S. Davydov, Solitons in Molecular Systems, Reidel, Dordrecht, 1985. M.A. Collins, Chem. Phys. Lett. 77 (1981) 342. D. Hochstrasser, F.G. Mertens, H. BTuttner, Phys. Rev. A 40 (1989) 2602. E. Ar$evalo, F.G. Mertens, Yu. Gaididei, A.R. Bishop, Phys. Rev. E 67 (2003) 016 610. C.W. Gardiner, Stochastic Methods, Springer, Berlin, 1990. R.L. Stratonovich, Topics in the Theory of Random Noise, Vol. II, Gordon and Breach, New York, 1967. P. HTanggi, Z. Physik B 36 (1980) 271. J.L. Wyatt, G.J. Coram, IEEE Trans. Electron Devices 46 (1999) 184. Ya. M. Blanter, M. BTuttiker, Phys. Rep. 336 (2000) 1. L.D. Landau, E.M. Lifshitz, Theory of Elasticity, Pergamon Press, Oxford, 1986. E. Ar$evalo, Yu. Gaididei, F.G. Mertens, Eur. Phys. J. B 27 (2002) 63. F. Marchesoni, Phys. Lett. A 115 (1986) 29. G. Cattuto, F. Marchesoni, Europhys. Lett. 62 (2003) 363. S. Habib, G. Lythe, Phys. Rev. Lett. 84 (2000) 1070. P. HTanggi, F. Marchesoni, P. Sodano, Phys. Rev. Lett. 60 (1988) 2563. M. Meister, F.G. Mertens, A. S$anchez, Eur. Phys. J. B 20 (2001) 405. N.R. Quintero, A. S$anchez, F.G. Mertens, Eur. Phys. J. B 16 (2000) 361. T. Taniuti, C. Wei, Phys. Soc. Jpn. 21 (1968) 209. A.M. Samsonov, Strain Solitons in Solids and How to Construct them, Chapman & Hall, CRC, St. Petersburg, Russia, 2001. R.L. Herman, J. Phys. A 23 (1990) 1063. E. Mann, J. Math. Phys. 38 (1997) 3772. R. Grimshaw, H. Mitsudera, Stud. Appl. Math. 90 (1993) 75. E. Wong, M. Zakai, Ann. Math. Stat. 36 (1965) 1560. E. Wong, M. Zakai, Z. Wahrscheinlichkeitstheorie Verw. Geb. 12 (1969) 87. V. Konotop, L. Vazquez, Nonlinear Random waves, World Scienti>c, Singapore, 1994. P.E. Kloeden, E. Platen, Numerical solution of Stochastic Di erential Equations, Springer, New York, 1997. T. Kamppeter, F.G. Mertens, E. Moro, A. Sanchez, A.R. Bishop, Phys. Rev. B 59 (1999) 11 349. K. Huang, Statistical Mechanics, 2nd Edition, Wiley, New York, 1987.