- Email: [email protected]
Correlation exponents for trajectories in the low-dimensional discrete selftrapping equation
Volume 110A, number 7,8
PHYSICS LETTERS
12 August 1985
CORRELATION EXPONENTS FOR TRAJECTORIES IN T H E L O W - D I M E N S I O N A L D I S C R E T E S E L F T R A P P I N G E Q U A T I O N ~ J.H. J E N S E N , P.L. C H R I S T I A N S E N , J.N. E L G I N t, J.D. G I B B O N a and O. S K O V G A A R D
Laboratoryof Applied MathematicalPhysics, The Technical Universityof Denmark, DK-2800 Lyngby, Denmark Received 30 May 1985; accepted for publication 4 June 1985
The non-linear dynamics of small polyatomic molecules has recently been modeled by the discrete selftrapping equation (DST) for which an alternative formulation is given. In the three-bond case of this hamiltonian system chaos occurs in the intermediate range of the non-linearity to dispersion ratio. Here calculations of the Grassberger-Procaccia correlation exponent exhibit noninteger values.
The recent progress in the field of bioenergetics was initiated in the early 70's by Davydov and his coworkers [ 1,2] who provided a theoretical explanation for the high efficiency in the transfer and storage of energy in proteins. They suggested that the amide-I vibrations of the peptide groups in an alpha-helix chain, through interaction with low frequency phonons along the chain, could form solitons. Davydov's model consists of a system of non.linear ordinary differential equations. In the continuous limit, where the number of peptide groups is very large, Davydov's equations become the non-linear Schr6dinger equation (NLS). This limiting case is studied in detail in ref. [3]. In 1984 a generalized model for the dynamics of small molecules of arbitrary geometry was proposed [4]. The model was denoted the discrete self-trapping equation (DST). In the linear case this system of ordinary differential equations was studied by Anderson [5]. It is noteworthy that DST describes the nonlinear dynamics of small polyatomic molecules such as water, ammonia, methane etc. [6]. The standard form of the DST equation is [4] JJi + 7 D ( I A I 2 ) A + eMA = 0 ,
(1)
Supported by the Danish Council for Scientific and Industrial Research and Phillips Foundation of 1958. i Permanent address: Department of Mathematics, Imperial College, London SW7 2BZ, UK. 0.375-9601/85/$ 03.30 © Elsevier Science Publishers B.V. (North-Holland Physics Publishing Division)
where A = col(A 1, A2 ..... An) is a complex n-component vector each component of which represents the probability amplitude of finding some conserved quantity on the nth subunit of the structure. In the nonlinear term, with coefficient 7, D denotes the diagonal matrix diag( IA 112, IA 212 ..... IAn 12). This term represents the tendency of A to self-trap through a nonlinear interaction with the adjacent structure. In the dispersion term, with coefficient e, M is a real symmetric matrix with zero diagonal elements. We shall assume all off-diagonal elements in M to be unity. The system (1) is hamiltonian with Z~ IAj 12 conserved. In terms of the density matrix, p, with elements ~g
Ojk = Aj A k , the DST equation (1) can be rewritten as
Ojk = i [O, H]jk, where H/k are the elements of the matrix H = - T D ( I A 12) - e M Using SU(n) notation and expanding p and bl in terms o f the corresponding generators we have reduced the DST equation for n = 2 to a pendulum equation in agreement with the fact that the system was found to be integrable with periodic solutions in ref. [4]. For n = 3 we get a complicated system of equations which will be reported elsewhere. In this case it was shown in ref. [4] that the solution becomes periodic in the limit 7/e -> 0 and quasiperiodic when 7/e ~ oo. For intermediate values of the non429
Volume 110A, number 7,8
PHYSICS LETTERS
linearity to dispersion ratio, "//el chaotic solutions were found computationaUy. In the present note we calculate the fractal dimensions for typical trajectories in this interval. We integrate eq. (1) from initial conditions (A 1 (0), A 2 (0), A3(0)) = (1, 0.001i, 0) using the IMSL routine DVERK [7]. The two parameters SCALE and TOL were chosen to 2 and 10 -8, respectively. The integration was carried out from t = 0 to t = 4500. The performance of the integration routine was checked by means of the conserved quantity 2~3=1 ~4/( t)12" The results agreed well with numerical results for the same initial conditions in ref. [4]. For a strange attractor Grassberger and Procaccia [8] computed the correlation integral given by N C(r)= lim N - 2 ~ O(r-tX i-Xjl) (2)
N--,*~
l (C) 0
'
7' =
000
~"
v
T ~ --5
~
3.75
0
a
ln( )
1.(c) 7
: 800 #°
i,/=l
where 0 is the Heaviside function and X i and Xj are points on the trajectories, and related the correlation integral to the exponent v through
5 -10 /
(3)
c ( r ) o: r v
They argued that u is close to (and never greater than) the fractal dimension of the trajectory. These authors [9] also proposed an effective algorithm for calculation of the correlation integral. Here we apply their method to the DST equation (1) with n = 3. Other hamiltonian systems have been treated in ref. [10]. It is clearly not fully understood how this method applies to hamiltonian systems, but here we present a second example. By means of a bit manipulating procedure the exponents of the real numbers IX i - X j l are extracted. In FORTRAN (IBM extension [11 ]) the function ISHFT performs a logical bit shift and is used for very fast collection of the number of distances having a fixed exponent. For the DST system the binning is narrowed by using IXi - Xjl 8 instead of IXi - Xil to get a suitable number of points in the diagrams showing In Cversus In r. Since the IBM 3081 computer operates in base 16 the distances are grouped in the intervals .... ] 16 -1/8, 16°], ] 160, 161/8], ] 161/8, 162/8 ] .... by the procedure. The correlation integral is now obtained as k
C(16k/8)=N -2 ~
430
12 August 1985
/
500
~ v
-15-
,-, 3 . 9 5
-5
0
in(r)
b
l.(C)
--
I//'/// / 6.00
I
7
=
vz % 3 . 8 7
5
C
0
Fig. 1. C versus r given by (2) for solutions to (1) with n = 3 and e = 1. (a): 3"= 0, 3; (b): 3' = 5, 8;and (c): 3"= 6. Correlation exponents v are indicated. where N i is the number of distances in the interval
] 16 (i-1)/8, 16i/8]. Ni,
(4)
Fig. 1 shows our computations for the dispersion
PHYSICS LETTERS
Volume 110A, number 7,8
4 3 2 1
0
~-~T_ ' ~ ' I ~ ' 4 8
[ 12
12 August 1985
six-dimensional phase space) and increases to a maxim u m value u - 3.95 for 3' ----5. At 3' = 6 two u-values occur (as shown in fig. lc). For higher values of 7 the correlation exponent approaches the value two corresponding to a trajectory on a two-dimensional surface in six-dimensional phase space. We conclude that the fractal dimension of the solution trajectory for the DST equation with n = 3 - calculated as the correlation exponent - exhibits a maxim u m for intermediate values of the ratio 7/e in accordance with the occurrence of chaos in this region as found in ref. [3].
7 Fig. 2. Correlation exponent v as a function of 7 for solution of (1) with n = 3 and e = 1. Error bonds indicated by bars.
P. Grassberger is acknowledged for drawing our attention to ref. [10].
References parameter e = 1 and the non-linearity parameter 7 = 0, 3, 5, 6, and 8. In fig. l a and l b the points approximately lie on a straight line - confirming eq. (3) the slope of which determines the correlation exponent v as indicated for each case. For 3' = 6 (fig. lc) the points approach lines with different slopes (u 1 ----2.01 and u 2 --- 3.87). This agrees with fig. 11~ in ref. [3] where the projection of the trajectory in sixdimensional solution space onto the Re(A 1), Im(A 1) plane shows a "two.band" structure. The inner and outer bands correspond to P2 and Vl, respectively. Fig. 2 shows the correlation exponent v as a function of the non-linearity parameter 3' for timed value of the dispersion parameter e. Error bonds, indicated by vertical bars, are obtained by means of the statistics routine BEMIRI from ref. [7]. The correlation exponent is unity for 3' = 0 (where DST is integrable and the solution trajectory becomes a closed curve in
[1 ] A.S. Davydov and N.I. Kislukha, Phys. Stat. Sol. 59b (1973) 465. [2] A.S. Davydov, J. Theor. Biol. 38 (1973) 559. [3] A.C. Scott, Phys. Scr. 25 (1982) 651. [4] J.C. Eilbeck, P.S. Lomdahl and A.C. Scott, to be published in Physica D. [5] P.W. Anderson, Phys. Rev. 109 (1958) 1492. [6] A.C. Scott, P.S. Lomdahl and J.C. Eilbeck, to be published. [7] International Mathematical and Statistical Library, Ed. 9 (Houston, 1982). [8] P. Grassberger and I. Procaccia, Phys. Rev. Lett. 50 (1983) 346. [9] P. Grassberger and I. Procaccia, Physica 9D (1983) 189. [10] M. Pettini and A. Vulpiani, Phys. Lett. 106A (1984) 207. [11] IBM. VS FORTRAN Application Programming: Language Reference GC26-3986-3. Release 3.0 (San Jose. 1983).
431
PHYSICS LETTERS
12 August 1985
CORRELATION EXPONENTS FOR TRAJECTORIES IN T H E L O W - D I M E N S I O N A L D I S C R E T E S E L F T R A P P I N G E Q U A T I O N ~ J.H. J E N S E N , P.L. C H R I S T I A N S E N , J.N. E L G I N t, J.D. G I B B O N a and O. S K O V G A A R D
Laboratoryof Applied MathematicalPhysics, The Technical Universityof Denmark, DK-2800 Lyngby, Denmark Received 30 May 1985; accepted for publication 4 June 1985
The non-linear dynamics of small polyatomic molecules has recently been modeled by the discrete selftrapping equation (DST) for which an alternative formulation is given. In the three-bond case of this hamiltonian system chaos occurs in the intermediate range of the non-linearity to dispersion ratio. Here calculations of the Grassberger-Procaccia correlation exponent exhibit noninteger values.
The recent progress in the field of bioenergetics was initiated in the early 70's by Davydov and his coworkers [ 1,2] who provided a theoretical explanation for the high efficiency in the transfer and storage of energy in proteins. They suggested that the amide-I vibrations of the peptide groups in an alpha-helix chain, through interaction with low frequency phonons along the chain, could form solitons. Davydov's model consists of a system of non.linear ordinary differential equations. In the continuous limit, where the number of peptide groups is very large, Davydov's equations become the non-linear Schr6dinger equation (NLS). This limiting case is studied in detail in ref. [3]. In 1984 a generalized model for the dynamics of small molecules of arbitrary geometry was proposed [4]. The model was denoted the discrete self-trapping equation (DST). In the linear case this system of ordinary differential equations was studied by Anderson [5]. It is noteworthy that DST describes the nonlinear dynamics of small polyatomic molecules such as water, ammonia, methane etc. [6]. The standard form of the DST equation is [4] JJi + 7 D ( I A I 2 ) A + eMA = 0 ,
(1)
Supported by the Danish Council for Scientific and Industrial Research and Phillips Foundation of 1958. i Permanent address: Department of Mathematics, Imperial College, London SW7 2BZ, UK. 0.375-9601/85/$ 03.30 © Elsevier Science Publishers B.V. (North-Holland Physics Publishing Division)
where A = col(A 1, A2 ..... An) is a complex n-component vector each component of which represents the probability amplitude of finding some conserved quantity on the nth subunit of the structure. In the nonlinear term, with coefficient 7, D denotes the diagonal matrix diag( IA 112, IA 212 ..... IAn 12). This term represents the tendency of A to self-trap through a nonlinear interaction with the adjacent structure. In the dispersion term, with coefficient e, M is a real symmetric matrix with zero diagonal elements. We shall assume all off-diagonal elements in M to be unity. The system (1) is hamiltonian with Z~ IAj 12 conserved. In terms of the density matrix, p, with elements ~g
Ojk = Aj A k , the DST equation (1) can be rewritten as
Ojk = i [O, H]jk, where H/k are the elements of the matrix H = - T D ( I A 12) - e M Using SU(n) notation and expanding p and bl in terms o f the corresponding generators we have reduced the DST equation for n = 2 to a pendulum equation in agreement with the fact that the system was found to be integrable with periodic solutions in ref. [4]. For n = 3 we get a complicated system of equations which will be reported elsewhere. In this case it was shown in ref. [4] that the solution becomes periodic in the limit 7/e -> 0 and quasiperiodic when 7/e ~ oo. For intermediate values of the non429
Volume 110A, number 7,8
PHYSICS LETTERS
linearity to dispersion ratio, "//el chaotic solutions were found computationaUy. In the present note we calculate the fractal dimensions for typical trajectories in this interval. We integrate eq. (1) from initial conditions (A 1 (0), A 2 (0), A3(0)) = (1, 0.001i, 0) using the IMSL routine DVERK [7]. The two parameters SCALE and TOL were chosen to 2 and 10 -8, respectively. The integration was carried out from t = 0 to t = 4500. The performance of the integration routine was checked by means of the conserved quantity 2~3=1 ~4/( t)12" The results agreed well with numerical results for the same initial conditions in ref. [4]. For a strange attractor Grassberger and Procaccia [8] computed the correlation integral given by N C(r)= lim N - 2 ~ O(r-tX i-Xjl) (2)
N--,*~
l (C) 0
'
7' =
000
~"
v
T ~ --5
~
3.75
0
a
ln( )
1.(c) 7
: 800 #°
i,/=l
where 0 is the Heaviside function and X i and Xj are points on the trajectories, and related the correlation integral to the exponent v through
5 -10 /
(3)
c ( r ) o: r v
They argued that u is close to (and never greater than) the fractal dimension of the trajectory. These authors [9] also proposed an effective algorithm for calculation of the correlation integral. Here we apply their method to the DST equation (1) with n = 3. Other hamiltonian systems have been treated in ref. [10]. It is clearly not fully understood how this method applies to hamiltonian systems, but here we present a second example. By means of a bit manipulating procedure the exponents of the real numbers IX i - X j l are extracted. In FORTRAN (IBM extension [11 ]) the function ISHFT performs a logical bit shift and is used for very fast collection of the number of distances having a fixed exponent. For the DST system the binning is narrowed by using IXi - Xjl 8 instead of IXi - Xil to get a suitable number of points in the diagrams showing In Cversus In r. Since the IBM 3081 computer operates in base 16 the distances are grouped in the intervals .... ] 16 -1/8, 16°], ] 160, 161/8], ] 161/8, 162/8 ] .... by the procedure. The correlation integral is now obtained as k
C(16k/8)=N -2 ~
430
12 August 1985
/
500
~ v
-15-
,-, 3 . 9 5
-5
0
in(r)
b
l.(C)
--
I//'/// / 6.00
I
7
=
vz % 3 . 8 7
5
C
0
Fig. 1. C versus r given by (2) for solutions to (1) with n = 3 and e = 1. (a): 3"= 0, 3; (b): 3' = 5, 8;and (c): 3"= 6. Correlation exponents v are indicated. where N i is the number of distances in the interval
] 16 (i-1)/8, 16i/8]. Ni,
(4)
Fig. 1 shows our computations for the dispersion
PHYSICS LETTERS
Volume 110A, number 7,8
4 3 2 1
0
~-~T_ ' ~ ' I ~ ' 4 8
[ 12
12 August 1985
six-dimensional phase space) and increases to a maxim u m value u - 3.95 for 3' ----5. At 3' = 6 two u-values occur (as shown in fig. lc). For higher values of 7 the correlation exponent approaches the value two corresponding to a trajectory on a two-dimensional surface in six-dimensional phase space. We conclude that the fractal dimension of the solution trajectory for the DST equation with n = 3 - calculated as the correlation exponent - exhibits a maxim u m for intermediate values of the ratio 7/e in accordance with the occurrence of chaos in this region as found in ref. [3].
7 Fig. 2. Correlation exponent v as a function of 7 for solution of (1) with n = 3 and e = 1. Error bonds indicated by bars.
P. Grassberger is acknowledged for drawing our attention to ref. [10].
References parameter e = 1 and the non-linearity parameter 7 = 0, 3, 5, 6, and 8. In fig. l a and l b the points approximately lie on a straight line - confirming eq. (3) the slope of which determines the correlation exponent v as indicated for each case. For 3' = 6 (fig. lc) the points approach lines with different slopes (u 1 ----2.01 and u 2 --- 3.87). This agrees with fig. 11~ in ref. [3] where the projection of the trajectory in sixdimensional solution space onto the Re(A 1), Im(A 1) plane shows a "two.band" structure. The inner and outer bands correspond to P2 and Vl, respectively. Fig. 2 shows the correlation exponent v as a function of the non-linearity parameter 3' for timed value of the dispersion parameter e. Error bonds, indicated by vertical bars, are obtained by means of the statistics routine BEMIRI from ref. [7]. The correlation exponent is unity for 3' = 0 (where DST is integrable and the solution trajectory becomes a closed curve in
[1 ] A.S. Davydov and N.I. Kislukha, Phys. Stat. Sol. 59b (1973) 465. [2] A.S. Davydov, J. Theor. Biol. 38 (1973) 559. [3] A.C. Scott, Phys. Scr. 25 (1982) 651. [4] J.C. Eilbeck, P.S. Lomdahl and A.C. Scott, to be published in Physica D. [5] P.W. Anderson, Phys. Rev. 109 (1958) 1492. [6] A.C. Scott, P.S. Lomdahl and J.C. Eilbeck, to be published. [7] International Mathematical and Statistical Library, Ed. 9 (Houston, 1982). [8] P. Grassberger and I. Procaccia, Phys. Rev. Lett. 50 (1983) 346. [9] P. Grassberger and I. Procaccia, Physica 9D (1983) 189. [10] M. Pettini and A. Vulpiani, Phys. Lett. 106A (1984) 207. [11] IBM. VS FORTRAN Application Programming: Language Reference GC26-3986-3. Release 3.0 (San Jose. 1983).
431