Statistical-mechanical theory of nonlinear density fluctuations near the glass transition

Physica A 395 (2014) 31–47

Contents lists available at ScienceDirect

Physica A journal homepage: www.elsevier.com/locate/physa

Statistical-mechanical theory of nonlinear density fluctuations near the glass transition Michio Tokuyama ∗ High Efficiency Rare Elements Extraction Technology Area, IMRAM, Tohoku University, Sendai 980-8577, Japan

highlights • • • •

A new time-convolutionless equation is derived for the density fluctuations. A new equation is found for the Debye–Waller factor. The critical temperature is expected to be much lower than that of ideal MCT. The non-Gaussian parameter is shown to be zero at an initial time.

article

info

Article history: Received 15 May 2013 Received in revised form 25 September 2013 Available online 24 October 2013 Keywords: Critical temperature Debye–Waller factor Glass transition Nonlinear density fluctuations Projection-operator methods Supercooled liquids

abstract The Tokuyama–Mori type projection-operator method is employed to study the dynamics of nonlinear density fluctuations near the glass transition. A linear non-Markov timeconvolutionless equation for the scattering function Fα (q, t ) is first derived from the Newton equation with the memory function ψα (q, t ), where α = c for the coherent–intermediate scattering function and s for the self–intermediate scattering function. In order to calculate ψα (q, t ), the Mori type projection-operator method is then used and a linear non-Markov time-convolution equation for ψα (q, t ) is derived with the memory function ϕα (q, t ). In order to calculate ϕα (q, t ), the same binary approximation as that used in the mode-coupling theory (MCT) is also employed and hence ϕα (q, t ) is shown to be identical with that obtained by MCT. Thus, the coupled equations are finally derived to calculate the scattering functions, which are quite different from the so-called ideal MCT equation. The most important difference between the present theory and MCT appears in the Debye–Waller factor fα (q). In MCT it is given by fα (q) = Γα (q)/(Γα (q) + 1), where Γα (q) is the long-time limit of the memory function ϕα (q, t ). On the other hand, in the present theory it is given by fα (q) = exp[−1/Γα (q)]. Thus, it is expected that the critical temperature Tc of the present theory would be much lower than that of MCT. The other differences are also discussed. © 2013 Elsevier B.V. All rights reserved.

1. Introduction The main purpose of the present paper is to formulate a general scheme for deriving basic equations for density fluctuations near the glass transition, and thus to clarify the role of nonlinear density fluctuations in supercooled liquids from first principles. A well-known example of this kind is the mode-coupling theory (MCT) equation [1,2]. A rigorous formulation has been established to predict a critical temperature Tc , below which the MCT solutions reduce to non-zero values even in the long-time limit [3]. In fact, the MCT equations have been solved numerically for various systems [4–12] to confirm this prediction. As long as the system is in equilibrium, this situation is exactly the same as that in critical phenomena,

∗

Tel.: +81 22 217 5327; fax: +81 22 217 5327. E-mail addresses: [email protected], [email protected].

0378-4371/$ – see front matter © 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.physa.2013.10.028

32

M. Tokuyama / Physica A 395 (2014) 31–47

except that the crossover from an equilibrium state to a nonequilibrium state occurs at the glass transition temperature Tg , which should be higher than Tc . Thus, MCT has been applied for wide varieties of glass-forming liquids to understand the mechanism of the glass transition [13–27]. Although MCT was the origin of all later studies on the glass transition, we still find a few unsatisfactory points. First, the critical temperature Tc of MCT is close to the melting temperature Tm , which is much higher than Tg . As a result, the peak height of the non-Gaussian parameter α2 (t ) [28] is always much smaller than 1 since the system is in a liquid state. This is a first problem (i). Secondly, the initial value of α2 (t ) is not zero but α2 (t = 0) = −2/3. This is a second problem (ii). Thirdly, the MCT solutions do not satisfy such a universality in dynamics that all the dynamics in different systems coincide with each other if the values of their self-diffusion coefficients are identical [29]. In fact, the MCT solutions have been shown not to agree with the simulation results in the cage region [30,11]. This is a third problem (iii). In the previous paper [31,32], therefore, we have proposed the alternative MCT (AMCT) equations to solve those problems. As is shown in Section 2, the problems (i) and (iii) are partially solved. In the present paper, therefore, we propose a theoretical approach completely different from the previous ones and show how one could avoid the above three problems safely within the MCT binary approximation [3]. Our approach is mainly based on the fact that the self–intermediate scattering function Fs (q, t ) is well-known to be described by the following cumulant expansion form [33]:



Fs (q, t ) = exp −q

2

t



Ψ (q, s)ds



0

 = exp −

q2 6

M2 (t ) +

q4



2

M2 (t )

2

6

 α2 (t ) + · · ·

(1)

with the non-Gaussian parameter [28]

α2 (t ) =

3M4 (t ) 5M2 (t )2

− 1,

(2)

where M2n (t ) = ⟨|Xj (t ) − Xj (0)|2n ⟩, Xj (t ) being a position vector of a particle j at time t. By taking the time derivative of Eq. (1), one can then find

∂ Fs (q, t ) = −q2 Ψ (q, t )Fs (q, t ). ∂t

(3)

This equation is different from the so-called MCT equation. In fact, the MCT equations for the scattering functions Fα (q, t ) are given by

 t 2 ∂2 q2 vth ∂ Fα (q, t ) = − Fα (q, t ) − ϕα (q, t − s) Fα (q, s)ds, (4) ∂t2 Sα (q) ∂ s 0 where ϕα (q, t ) is the memory function, vth the thermal velocity, and Sα (q) = Fα (q, 0). Here α = c stands for the coherent–intermediate scattering function and s for the self–intermediate scattering function, where Ss (q) = 1. Then, MCT [3] shows that an ergodic to nonergodic transition occurs at a critical temperature Tc , below which the long-time solution reduces to a non-zero value fα (q). For T ≤ Tc , one can thus find the so-called Debye–Waller factor fα (q) as fα (q) = lim

t →∞

Fα (q, t ) Sα (q)

=

Γα (q) Γα (q) + 1

(5)

with the long-time limit of the memory function

Γα (q) = lim

z →0

z ϕα [q, z ] 2 q2 vth

Sα (q),

(6)

where ϕα [q, z ] is a Laplace transform of ϕα (q, t ). One of the most important results in the present paper is to predict that the critical temperature Tc would be much lower than that of MCT. In fact, the present theory leads to the coupled equations

 t ∂ Fα (q, t ) = −q2 ψα (q, s)dsFα (q, t ), (7) ∂t 0  t ∂ ψα (q, t ) = − ϕα (q, s)ψα (q, t − s)ds, (8) ∂t 0 where the memory function ϕα (q, t ) coincides with that of MCT within the MCT binary approximation. Similarly to Eq. (5), by using Eqs. (7) and (8), one can then find for T ≤ Tc   1 fα (q) = exp − . (9) Γα (q)

M. Tokuyama / Physica A 395 (2014) 31–47

33

Here the value of Γα (q) in the present theory is always larger than that in MCT as a function of fα (q). As T decreases, Γα (q) increases. Hence this suggests that the critical temperature Tc of the present theory is lower than that of MCT. Within the MCT binary approximation, one can also transform Eq. (4) into

∂ Fα (q, t ) = −q2 ∂t

t



ψα (q, t − s)Fα (q, s)ds.

(10)

0

This equation should be compared with Eq. (7). Thus, we emphasize that the problems (i)–(iii) all depend on whether the memory term is convolution or convolutionless in time. We begin in Section 2 by reviewing the theoretical background. We discuss two types of projection-operator methods and compare the MCT equations with the AMCT equations. In Section 3, we first discuss how to derive the time-convolutionless generalized Langevin equation for the density fluctuations with the memory function ψα (q, t ). Then, we discuss how to derive the time-convolution equation for ψα (q, t ) with the memory function ϕα (q, t ). Thus, we compare the present results with the previous ones. We conclude in Section 4 with a summary. 2. Theoretical background As already pointed out in the previous paper [31], the first problem (i) strongly depends on which types of projectionoperator methods one employs, the Mori type [M] [34] or the Tokuyama–Mori type [T] [35]. In the present section, we briefly summarize and discuss the basic equations and main results, which are compared with the present ones. 2.1. Two types of projection-operator methods We first summarize two types of projection-operator methods, [M] and [T]. Let A(t ) denote a column vector of dynamical variables whose time evolution is described by the Heisenberg equation d

A(t ) = iLA(t ), dt where L is the Liouville operator. We now introduce a projection operator ℘ onto A(0) subspace by Ref. [34]

℘ B(t ) = ⟨B(t )A(0)∗ ⟩ · ⟨A(0)A(0)∗ ⟩−1 · A(0),

(11)

(12)

where the brackets denote the equilibrium ensemble average, B(t ) an arbitrary dynamical variable, and the asterisk the complex conjugate. Then, there are two types of projection operator methods to derive a generalized Langevin equation from Eq. (11), Mori type [34] and Tokuyama–Mori type [35]. As discussed in Refs. [34,35], starting from Eq. (11) and transforming it by employing two types of projection operator methods, one obtains

 t    ϕ(s) · A(t − s)ds + f (t ), [M] iω · A(t ) − d A(t ) = 0 t  dt  iω · A(t ) − ψ(s)ds · A(t ) + g (t ), [T]

(13)

0

where iω ≡ ⟨(iLA)A∗ ⟩ · ⟨AA∗ ⟩−1 , and

ϕ(t ) = ⟨f (t )f (0)∗ ⟩ · ⟨A(0)A(0)∗ ⟩−1 ,

(14)

ψ(t ) = ⟨g (t )g (0) ⟩ · ⟨A(t )A(0) ⟩ .

(15)

∗

∗ −1

In the following, we use the abbreviation [M] for Mori type equation and [T] for Tokuyama–Mori type equation. Here the fluctuating forces f (t ) and g (t ) are given by f (t ) = etQ iL Q iLA(0), g (t ) = e

(16)

[1 − Q {1 − e e }] Q iLA(0), = A˙ (t ) − ⟨A˙ (t )A∗ ⟩ · ⟨A(t )A∗ ⟩−1 · A(t ), ˙ (t ) = iLA(t ), and Q = 1 − ℘ . The fluctuating forces f (t ) and g (t ) satisfy the orthogonality conditions where A tQ iL

−tiL tQ iL

⟨f (t )A(0)∗ ⟩ = ⟨g (t )A(0)∗ ⟩ = 0.

−1

(17) (18) (19)

We should mention here that both equations [M] and [T] are mathematically identical unless any approximations are made for the fluctuating forces. Once some approximations are made for them, however, the time evolution of the correlation function ⟨A(t )A(0)∗ ⟩ obtained by solving both equations of types [M] and [T] shows a quite different behavior from each other [35]. 2.2. Ideal MCT equations In this paper, we consider the three-dimensional equilibrium glass-forming system, which consists of N particles with mass m and diameter σ in the total volume V at temperature T . Near the glass transition, the slowly-varying variables are

34

M. Tokuyama / Physica A 395 (2014) 31–47

then given by the density fluctuations ρc (q, t ) and ρs (q, t ) 1

ρc (q, t ) =

N 1/2

 N 

 ρs (q, t ) − N δq,0 ,

(20)

j =1

ρs (q, t ) = eiq·Xj (t ) , (21) where Xj (t ) denotes the position vector of the jth particle at time t. Then, the coherent–intermediate scattering function Fc (q, t ) and the self–intermediate scattering function Fs (q, t ) are given by Fα (q, t ) = ⟨ρα (q, t )ρα (q, 0)∗ ⟩,

(22)

where α = c and s. Here Fc (q, 0) = S (q) and Fs (q, 0) = 1, where S (q) is a static structure factor. We first briefly review the ideal MCT equations and its main results. From Eq. (11), the density fluctuations ρα (q, t ) obey

∂ ρc (q, t ) = iLρc (q, t ) = iqJα (q, t ) ∂t with the current Jα (q, t ) given by N 

1

Jc (q, t ) =

N 1/2 j=1

Js (q, t ),

(23)

Js (q, t ) = qˆ · X˙ j (t )eiq·Xj (t ) ,

(24)

where qˆ = q/q and q = |q|. Following MCT, we take two variables ρα (q, t ) and Jα (q, t ) as A(t ) and then introduce the projection operator ℘α′′ by

℘α′′ B(t ) =

⟨B(t )ρα (q, 0)∗ ⟩ ⟨B(t )Jα (q, 0)∗ ⟩ ρα (q, 0) + Jα (q, 0), 2 Sα (q) vth

(25)

2 where ⟨Jα (q, 0)ρα (q, 0)∗ ⟩ = 0 and ⟨Jα (q, 0)Jα (q, 0)∗ ⟩ = vth . Here vth = (kB T /m)1/2 . Then, use of Eq. (13) of type [M] leads to 2 ∂ q2 vth Jα (q, t ) = − ρα (q, t ) − ∂t Sα (q)

t



ϕα (q, t − s)Jα (q, s)ds + fα (q, t )

(26)

0

with the memory function 2 , ϕα (q, t ) = ⟨fα (q, t )fα (q, 0)∗ ⟩/vth where the fluctuating force fα (q, t ) is given by

(27)

′′

fα (q, t ) = etQα iL Qα′′ iLJα (q, 0). Here

Qα′′

=1−

℘α′′ ,

(28)

and Sc (q) = S (q) and Ss (q) = 1. From Eq. (26), one can then obtain the so-called ideal MCT equations

2 q2 vth ∂ F ( q , t ) = − Fα (q, t ) − α ∂t2 Sα (q) 2

t



ϕα (q, t − s) 0

∂ Fα (q, s)ds. ∂s

It is convenient to introduce the Laplace transform of ϕα (q, t ) by ϕα [q, z ] = obtain

(29)

∞ 0

ϕα (q, t )e−zt dt. From Eq. (29), we then

 −1 v2 1 = z + q2 th . Sα (q) Sα (q) z + ϕα [q, z ]

Fα [q, z ]

(30)

By using Eq. (30), one can now discuss the asymptotic forms for a short and a long times. For a short time, one finds 2 Fα [q, z ] ≃ Sα (q)/[z + q2 vth /(zSα (q))],

(31)

Fα (q, t ) ≃ Sα (q) cos[qvth t /Sα (q)1/2 ].

(32)

For a long time, one obtains Fα [q, z ] Sα (q)

≃

2 z ϕα [q, z ]Sα (q)/(q2 vth ) 2 z ϕα [q, z ]Sα (q)/(q2 vth )+1

.

(33)

Here we should note that Eq. (32) is an origin of the second problem (ii). In fact, it leads to α2 (t = 0) = −2/3 since 2 2 4 4 M2 (t ) ≃ 3vth t and M4 (t ) ≃ 5vth t . 2.2.1. Nonlinear memory function We now show how to calculate the memory function ϕα (q, t ) within the MCT binary approximation [3]. In addition to the randomly-varying processes, the fluctuating force fα (q, t ) still contains nonlinear terms of density fluctuations since

M. Tokuyama / Physica A 395 (2014) 31–47

35

those are the slowly-varying variables. Hence it can be formally split up into two parts, the random force Rα (q, t ) and the nonlinear fluctuating force ξα (q, t ): fα (q, t ) = ξα (q, t ) + Rα (q, t ).

(34)

Here we assume that the random force Rα (q, t ) is uncorrelated with any functions of dynamical variables and obeys the Gaussian-Markov processes with zero mean, which satisfy 2 ⟨Rα (q, t )Rα (q, t ′ )∗ ⟩ = 2γα vth δ(t − t ′ ),

(35)

where γα is a positive constant. Then, the memory function ϕα (q, t ) can be written as [36–38]

ϕα (q, t ) = 2γα δ(t ) + 1ϕα (q, t )

(36)

with the nonlinear memory term 2 1ϕα (q, t ) = ⟨ξα (q, t )ξα (q, 0)∗ ⟩/vth .

(37)

We next calculate the nonlinear memory function 1ϕα (q, t ) formally. Following MCT [3], we introduce the projection operator ℘2α by

℘2c B =

 ⟨Bρc (q − k )∗ ρc (k )∗ ⟩ 2S (k)S (|q − k |)

k

℘2s B =

 ⟨Bρs (q − k )∗ ρc (k )∗ ⟩ S (k)

k

ρc (k )ρc (q − k ),

(38)

ρc (k )ρs (q − k ),

(39)

where ρα (k ) = ρα (k , 0). By using the operator identities ′′

′′

′′

etiL = etiL ℘2α + etiL Q2α , ′′

′′

etiL Q2α = etQ2α iL Q2α +

(40) t



e(t −s)iL ℘2α iL′′ esQ2α iL Q2α ds, ′′

′′

(41)

0

and Eq. (28), one finds ′′

ξα (q, t ) = etiL ℘2α [iL′′ Jα (q)], Rα (q, t ) = e

tQ2α iL′′

(42)

Q2α [iL′′ Jα (q)],

where Q2α = 1 − ℘2α and iL ⟨Rα (q, t )gα (q, 0)∗ ⟩ = 0.

′′

(43)

Qα′′ iL.

=

Here the second term of Eq. (41) becomes zero since ⟨Rα (q, t )ρα (q, 0) ⟩ = ∗

We now calculate the fluctuating force ξα (q, t ). By using the operator identity e

tiL′′

=e

tiL

t

 − 0

e(t −s)iL ℘α′′ iLesiL ds, ′′

(44)

one can first write the memory function 1ϕα (q, t ) as

1ϕα (q, t ) = 1φα (q, t ) +

t



1ϕα (q, t − s) 0

⟨Jα (q, s)ξα(0) (q, 0)∗ ⟩ ds, 2 vth

(45)

where

1φα (q, t ) =

⟨ξα(0) (q, t )ξα(0) (q, 0)∗ ⟩ , 2 vth

(46)

ξα(0) (q, t ) = etiL ℘2α [iL′′ Jα (q)].

(47)

Taking the Laplace transform of Eqs. (23) and (26), one can find

⟨Jα [q, z ]ξα(0) (q, 0)∗ ⟩ = 2 vth z+

1φα [q, z ] 2 q2 vth

zSα (q)

+ γα + 1φα [q, z ]

.

(48)

Use of Eqs. (45) and (48) then leads to z 1φα [q, z ]2

1ϕα [q, z ] = 1φα [q, z ] + z2 +

2 q2 vth

Sα (q)

+ z (γα − 1φα [q, z ])

.

(49)

Hence one can obtain the memory function 1ϕα [q, z ] if one can find 1φα [q, z ] in terms of Fα (q, t ). We discuss this next.

36

M. Tokuyama / Physica A 395 (2014) 31–47

Use of Eqs. (38) and (39) leads to 2 vth

℘2c [Qc′′ iLJc (q)] =



2iN 1/2

2  vth

℘2s [Qs′′ iLJs (q)] =

iN 1/2

vc (q, k )ρc (k )ρc (q − k ),

(50)

k

vs (q, k )ρc (k )ρs (q − k )

(51)

k

with the vertex amplitudes

vc (q, k ) = qˆ · kc (k) + qˆ · (q − k )c (|q − k |),

(52)

vs (q, k ) = qˆ · kc (k),

(53)

where c (k) = 1 − 1/S (k) and we have used the fact that ⟨ρc (q)ρc (k )∗ ρc (q − k )∗ ⟩ ≃ S (q)S (k)S (|q − k |)/N 1/2 . From Eq. (47), we then obtain

ξc(0) (q, t ) = ξs(0) (q, t ) =

2 vth



2iN 1/2

2  vth

iN 1/2

vc (q, k )ρc (k , t )ρc (q − k , t ),

(54)

k

vs (q, k )ρc (k , t )ρs (q − k , t ).

(55)

k

In order to calculate the nonlinear memory function 1φα (q, t ), we next employ the following factorization approximation used in MCT:

⟨ρc (k , t )ρα (q − k , t )ρc (k ′ , 0)∗ ρα (q − k ′ , 0)∗ ⟩ ≃ ⟨ρc (k , t )ρc (k ′ , 0)∗ ⟩ ⟨ρα (q − k , t )ρα (q − k ′ , 0)∗ ⟩ + ⟨ρc (k , t )ρα (q − k ′ , 0)∗ ⟩ ⟨ρα (q − k , t )ρα (k ′ , 0)∗ ⟩ = 21−nα δk ,k ′ Fc (q, t )Fα (|q − k ′ |, t ),

(56)

where nc = 0 and ns = 1, and we have used the thermodynamic limit N → ∞ and V → ∞ with ρ = N /V being fixed. Thus, using Eqs. (46), (54), (55), and (56), we finally obtain

v2 1φc (q, t ) = th 2ρ



2 vth ρ



1φs (q, t ) =

dk

(2π )3 dk

(2π )3

vc (q, k )2 Fc (k, t )Fc (|q − k |, t ),

(57)

vs (q, k )2 Fc (k, t )Fs (|q − k |, t ),

(58)

where we have used the fact that V −1 k → dk /(2π )3 . Those are identical to the nonlinear memory functions obtained by MCT. We should mention here that in MCT the second term of Eq. (49) is usually neglected. This approximation is ′′ equivalent to putting etiL ≃ etiL in ξα (q, t ). As mentioned in Ref. [39] as Q approximation, it may cause trouble in some cases. In the present system, this might be reasonable because the second term is of order O(1φα2 ). In fact, the density fluctuations are small as compared to their average values even near the glass transition [40]. This situation is quite different from that in the critical phenomena, where the density fluctuations become the same order as that of their average values near the critical point [41]. In the following, therefore, we just neglect the second term of Eq. (49) to discuss the MCT solutions within the MCT approximation.





2.2.2. Debye–Waller factor Next, we discuss the solutions of the MCT equations (29) at long times. The most important prediction of MCT is the ergodic to nonergodic transition at a critical temperature Tc , below which the long-time solutions reduce to non-zero values. In fact, following MCT [3], from Eq. (33), one can find the Debye–Waller factor fα (q) as Fα (q, t )

fα (q) = lim

Sα (q)

t →∞

= lim

zFα [q, z ]

z →0

Sα (q)

=

Γα (q) Γα (q) + 1

(59)

with the long-time limit of the memory function

Γα (q) = lim

z ϕα [q, z ]

z →0

=

2 q2 vth

1 2(2π )3



Sα (q),

dkVα(2) (q, k, |q − k |)fc (k)fα (|q − k |),

(60)

(61)

M. Tokuyama / Physica A 395 (2014) 31–47

37

where the vertex Vα(2) is given by Vα(2) (q, k, |q − k |) = 2nα Sα (q)Sc (k)Sα (|q − k |)vα (q, k )2 /(ρ q2 ).

(62)

One can also solve Eq. (59) for Γα (q) as

Γα (q) =

fα (q)

.

1 − fα (q)

(63)

Then, MCT predicts that there exist non-zero solutions fc (q) > 0 in the long time limit. This is possible if the vertices are sufficiently large enough. When the temperature decreases, the value of vertices increases. Then, the highest temperature at which the non-zero solution appears first is the so-called critical temperature Tc of MCT. 2.2.3. A single-particle dynamics Finally, we discuss the single-particle dynamics. Expanding Fs (|q − k |) in powers of q, from Eq. (58) we first obtain

1ϕs (q, t ) = 1ϕs(0) (t ) + q2 1ϕs(2) (t ) + · · · ,

(64)

where

1ϕs(0) (t ) =

2 vth 6π 2 ρ

∞



dkk4 c (k)2 Fc (k, t )Fs (k, t ),

(65)

0

2 vth 1ϕs (t ) = 20π 2 ρ

∞



(2)

∂2 dkk c (k) Fc (k, t ) + 2 3k ∂ k ∂k 4

0



2

2 ∂



Fs (k, t ).

(66)

By expanding Eq. (22) in powers of q as Fs (q, t ) = 1 − q2 M2 (t )/3! + q4 M4 (t )/5! + · · ·

(67)

and using Eq. (29), one also finds d2 dt

2 M2 (t ) = 6vth − 2

d2 dt

M4 (t ) = 20v 2

t



d

ϕs(0) (t − s)

ds

0

2 th M2

(t ) −

t



M2 (s)ds,

(0)

ϕs (t − s) 0

d ds

(68)

M4 (s)ds + 20

t



1ϕs(2) (t − s)

0

d ds

M2 (s)ds.

(69)

The Laplace transform of Eqs. (68) and (69) then leads to M2 [z ] =

M4 [z ] =

2 6vth

(0)

z 2 (z + γs + 1ϕs [z ]) 10 3

 zM2 [z ]

2

,

(2)

z 1ϕs [z ] 2 vth

(70)

 +1 .

(71)

From Eq. (68), one can then find the long-time self-diffusion coefficient D as D=

2 vth

γs + 1ϕs(0) [0]

.

(72)

From Eqs. (68) and (69), one can also find α2 (t = 0) = −2/3. This is the second problem (ii). 2.3. Alternative MCT equations Next, we briefly review the AMCT equations and its main results. As shown in the previous papers [31,32], we here take Eq. (13) of type [T]. In fact, using Eqs. (13) and (25) and employing the same binary approximation as that employed in MCT, one can derive the following alternative time-convolutionless MCT (AMCT) equations: 2 ∂2 q2 vth Fα (q, t ) = − Fα (q, t ) − 2 ∂t Sα (q)

 0

t

ϕα (q, s)ds

∂ Fα (q, t ). ∂t

(73)

Here we note that Eq. (73) has the same form as that of Eq. (29), except that the memory term is convolutionless in time. Within the MCT binary approximation, the memory function ϕα (q, t ) is also shown to have the same form as that of Eq. (29) [31,32]. However, we should mention here that their numerical values are different from each other because the starting Eqs. (29) and (73) are different. We now discuss the asymptotic forms of Eq. (73) for a short time and also for a long time. Since the second term of Eq. (73) is neglected for a short time, Eq. (73) then leads to the short time solutions given by Eq. (32). Hence the second problem (ii) still exists here, even though the memory function is convolutionless in time. On the other hand, since the second derivative

38

M. Tokuyama / Physica A 395 (2014) 31–47

term is neglected for a long time, one finds

 Fα (q, t ) ≃ Sα (q) exp −

2 q2 vth

t



Sα (q)

s

0



1 0

ϕα (q, τ )dτ

ds .

(74)

We next find the non-zero solutions fα (q) from Eq. (74). Similarly to Eq. (59), one can write Eq. (74) as

   ezt fα (q) = lim exp − dz K [q, z ] ,

(75)

z

t →∞

where K [q, z ] = z

∞



2 q2 vth

t



dte−zt

Sα (q)

0

0

s 0

ϕα (q, τ )dτ

ds.

(76)

In the limit of z → 0, K [q, z ] can be written as lim K [q, z ] = lim

z →0

z →0

= lim

z →0

2 q2 vth

∞



Sα (q)

ds  ∞

0

0

v

2 q2 th

z ϕα [q, z ]Sα (q)

=

e−zs dτ ϕα (q, τ )θ (s − τ ) 1

Γα (q)

.

(77)

Use of Eqs. (75) and (77) then leads to



fα (q) = exp −



1

Γα (q)

.

(78)

One can also solve Eq. (78) for Γα (q) as

Γα (q) = −

1 ln[fα (q)]

.

(79)

Eq. (78) is a new non-zero long time solution which should be compared with the MCT non-zero solution given by Eq. (59). This will be discussed in the next subsection. Finally, we discuss the single-particle dynamics. Use of Eqs. (64), (67), and (73) then leads to d2 dt

M2 (t ) = 6v − 2 2 th

t

 0

d2 dt

ϕs(0) (s)ds

2 M4 (t ) = 20vth M2 (t ) − 2

t

 0

d dt

M2 (t ),

ϕs(0) (s)ds

d dt

(80) M4 (t ) + 20

t

 0

ϕs(2) (s)ds

d dt

M2 (t ).

(81)

From Eq. (80), one can easily show that D has the same form as that of Eq. (72), although their temperature dependence is different from each other. From Eqs. (80) and (81), one can also find α2 (t = 0) = −2/3. This is the second problem (ii). 2.4. Critical temperature As shown in Eq. (78), the time-convolutionless equation (73) leads to a new non-zero solution, which is different from the MCT solution given by Eq. (59). Hence we now discuss how the difference between those solutions changes the critical temperature Tc . As is shown in Fig. 1, the value of Γα (q) in AMCT is always larger than that in MCT as a function of fα (q). Since Γα increases as temperature decreases, this suggests that Tc of AMCT is always lower than that of MCT. Even though the difference between two results is of order 0.5, this small difference might solve the serious problem that Tc of MCT is always higher than that of the simulation results. In order to prove this, we have first performed the extensive molecular-dynamics simulations on the one-component Lennard-Jones system with the potential given by U (r ) = 4ϵ[(σ /r )12 − (σ /r )6 ] under the NVT ensemble [42], where ϵ is energy unit and σ length unit. Length, time, and temperature are scaled by σ , σ /v0 , and ϵ/kB , respectively, where v0 = (ϵ/m)1/2 . Here the number density ρ is set as ρ = 0.9σ −3 . The reason why we choose this onecomponent system is mainly because this system becomes crystalline at the melting temperature Tm . Hence one can show that the MCT solutions predict a fictive supercooled state since there exists no supercooled state below Tm in the system. In Fig. 2, the long-time self-diffusion coefficient D(T ) versus scaled inverse temperature 1/T is shown for the simulation results, where the simulation results predict Tm ≃ 0.5. Then, Eqs. (29) and (73) have been solved numerically by using the static structure factor S (q) obtained by the molecular-dynamics simulations [42,43]. In Fig. 2, their long-time self-diffusion coefficients are also shown for comparison. In order to analyze each data from a unified point of view, it is convenient to introduce the master curve for fragile liquids recently proposed by the present author [44,45] D(T ) = d0

(1 − x)10/3 x

exp[62x13/3 (1 − x)10/3 ]

(82)

M. Tokuyama / Physica A 395 (2014) 31–47

39

Fig. 1. (Color online) A plot of Γα (q) versus fα (q). The dotted line indicates the numerical results obtained by Eq. (63) and the solid line by Eq. (79). The dashed line indicates a difference between those two lines.

Fig. 2. (Color online) A log plot of D(T )/σ v0 versus scaled inverse temperature 1/T . The symbols () indicate the simulation results on the one-component Lennard-Jones system at the number density ρσ 3 = 0.9 from Ref. [42], () the numerical solutions of Eq. (29) from Ref. [43], and (⊙) the numerical solutions of Eq. (73) from Ref. [43]. The dotted line indicates the master curve given by Eq. (82) at d0 = 0.0118σ v0 and 1/Tf = 5.806, the dashed line at d0 = 0.0118σ v0 and 1/Tf = 2.464, and the solid line at d0 = 0.0039σ v0 and 1/Tf = 5.806. The symbols (•) indicate the numerical results of Eq. (73) times 3.

with x = 1 − Tf /T , where Tf is a fictive singular temperature to be determined and d0 a positive constant. By analyzing many data in various glass-forming liquids [44–46], as long as the system is in equilibrium, all the diffusion data have been shown to be described by Eq. (82) well up to a deviation point Tn , below which the system becomes out of equilibrium and the data start to deviate from the master curve. Thus, Tf and d0 are obtained by fitting the date with Eq. (82) and are listed in Table 1. As is shown in Table 1, the singular temperature Tf ≃ 0.172 for the AMCT solutions coincides with that for the simulation

40

M. Tokuyama / Physica A 395 (2014) 31–47 Table 1 Fitting values of Tf and d0 for different cases. System

Tf

d0 /σ v0

Simulations AMCT MCT

0.172 0.172 0.406

0.0118 0.0039 0.0118

Fig. 3. (Color online) A log–log plot of M2 (t )/σ 2 versus scaled time T 1/2 t. The symbols (◦) indicate the simulation results on the one-component LennardJones system at T = 0.667 from Ref. [42] and () the β -relaxation time given by log10 (tβ ) = 0.825. The dotted line indicates the solution of MCT at T = 2.0 and the solid line the solution of AMCT at T = 1.67 from Ref. [43].

results, while that for the MCT solutions is Tf ≃ 0.406. Thus, the singular temperature of MCT is much higher than that of the simulation results and the AMCT solutions. The non-Gaussian parameter α2 (t ) of MCT is always much smaller than 1. This is reasonable because the system is in a liquid state for T > Tm (≃ 0.5). Thus, D(T ) of MCT turns out to predict a fictive supercooled state near Tm , even though the present system is in a equilibrium liquid state for T > Tm and becomes crystalline for T ≤ Tm . On the other hand, since Tf of AMCT is the same as that of the simulation results, α2 (t ) of AMCT is expected to become larger than 1 if the original system is chosen from one of glass-forming liquids, which shows a supercooled state. Next, we discuss the value of d0 . As is shown in Table 1, it is given by 85σ v0 for the MCT solutions. This value coincides with that for the simulation results. This is reasonable because the mechanism of MCT equations should be the same as that of simulations [44,45]. On the other hand, d0 for the AMCT solutions is nearly three times larger than others. The main reason why this difference appears is not known yet, although the singular temperature of AMCT coincides with that of the simulations. Since Tf of AMCT is identical to that of the simulations, however, all the self-diffusion data for the AMCT solutions should coincide with the simulation results if they are multiplied by 3. In fact, this is also shown in Fig. 2, where both results completely coincide with each other. Hence this fact suggests that the AMCT solutions do not have the third problem (iii). In fact, this is seen in Fig. 3, where the scaled self-diffusion coefficient D(T )/T 1/2 has the same value for three different cases; the simulation at T = 0.667, AMCT at T = 1.67, and MCT at T = 2.0. Thus, the solution of AMCT for M2 (t ) at T = 1.67 agrees with the simulation results at T = 0.667, while the solution of MCT for M2 (t ) at T = 2.0 does not agree with others in the cage region up to tβ . Here we note that the deviation from the simulation results in a cage region always appears above the simulation results in a liquid state, while it appears below them in a supercooled state [11,30]. Finally, we should mention that AMCT still have two unsatisfactory points. One is the second problem (ii), where α2 (t = 0) = −2/3 ̸= 0. The other is a problem that the value of d0 is not identical to that of the simulations. Hence Eqs. (73) are still not adequate equations to study the dynamics of supercooled liquids correctly, even though the problems (i) and (iii) are solved partially. In order to overcome those problems, therefore, we next propose time-convolutionless MCT (TMCT) equations, which are quite different from Eqs. (73). Thus, the problems (i)–(iii) are shown to be safely solved. However, one still needs to solve the TMCT equations numerically for glass-forming liquids to check whether the value of d0 coincides with that of simulations.

M. Tokuyama / Physica A 395 (2014) 31–47

41

3. TMCT equation In the present section, we discuss how to derive a new type of time-convolutionless equations for the density fluctuations

ρα (q, t ).

3.1. Basic equation for density fluctuations As shown in the previous sections, the initial value of the non-Gaussian parameter α2 (t ) is given by α2 (t = 0) = −2/3 in both MCT and AMCT cases. This is due to the fact that in both cases the starting equations are the second-order differential equations in time, although the memory term is convolution in time for MCT while it is convolutionless in time for AMCT. In both cases, the density fluctuation ρα (q, t ) and the current Jα (q, t ) were taken as the relevant variable A(t ) and the projection operator ℘α′′ given by Eq. (25) was used to derive the second-order differential equations for Fα (q, t ). Near the glass transition, however, the slowly-varying variables are only the density fluctuations. Hence we now take the density fluctuation ρα (q, t ) as the relevant variable A(t ) and introduce the projection operator ℘α from Eq. (12) by

℘α B(t ) =

⟨B(t )ρα (q, 0)∗ ⟩ ρα (q, 0). Sα (q)

(83)

Then, using Eq. (13) of type [T] exactly leads to

∂ ρα (q, t ) = −q2 ∂t

t



ψα (q, s)dsρα (q, t ) + iqgα (q, t )

(84)

0

with the memory function

ψα (q, t ) = ⟨gα (q, t )gα (q, 0)∗ ⟩/Fα (q, t ),

(85)

where iωα = 0. The fluctuating force gα (q, t ) is given by gα (q, t ) = etQα iL [1 − Qα {1 − e−tiL etQα iL }]−1 Jα (q, 0),

(86)

where Qα = 1 − ℘α , and gα (q, t ) satisfies the orthogonality condition give by Eq. (19). By using ⟨gα (q, t )ρα (q, 0) ⟩ = 0, from Eqs. (22) and (84), we then obtain ∗

∂ Fα (q, t ) = −q2 ∂t

t



ψα (q, s)dsFα (q, t ).

(87)

0

This equation is easily solved to give the formal solution



Fα (q, t ) = Sα (q) exp −q

2

t





(t − s)ψα (q, s)ds .

(88)

0

For a short time, Eq. (88) then reduces to Fα (q, t )

  2 q2 vth ≃ exp − t2 , Sα (q) 2Sα (q)

(89)

2 where ψα (q, t = 0) = vth /Sα (q). This leads to α2 (t = 0) = 0. Hence Eq. (87) does not have the second problem (ii).

3.2. Equation for memory function In order to calculate the memory function ψα (q, t ), one has to derive an equation for the fluctuating force gα (q, t ) given by Eq. (86). We note here that gα (q, t ) consists of randomly-varying variables and slowly-varying variables, such as nonlinear terms of density fluctuations. Since it is not an easy task to treat it, one has to find an approximate way. In fact, by taking a time derivative of Eq. (86), one can obtain d

gα (q, t ) = Qα iLgα (q, t ) − iqψα (q, t )Qα ρα (q, t ). (90) dt Next we show that the second term of Eq. (90) is negligible as long as the time scale of Qα ρα (q, t ) is much larger than that of ψα (q, t ). We first introduce the following two different quantities by Dα (q, t ) =

t



ψα (q, s)ds = − 0

Ωα (q, t ) =

t



Dα (q, s)ds = − 0

1 ∂ q2 ∂ t 1 q2

ln[Fα (q, t )],

 ln

Fα (q, t ) Sα (q)



.

(91)

(92)

42

M. Tokuyama / Physica A 395 (2014) 31–47

Fig. 4. (Color online) A plot of ψα (q, t ) and Fα (q, t )/Sα (q) versus scaled time t /t0 for A particle in the Stillinger–Weber binary mixtures A80 B20 [48] at scaled temperature T = 0.714. The solid lines indicate the simulation results from Ref. [47] for ψα (q, t ) and the dashed lines for Fα (q, t )/Sα . The symbols (+) indicate the simulation results for ψα (q, t )(Fα (q, t )/Sα (q)).

Fig. 5. (Color online) A log–log plot of Dα (q, t ) and Ωα (q, t ) versus scaled time t /t0 for A particle in the Stillinger–Weber binary mixtures A80 B20 [48] at scaled temperature T = 0.714. The solid lines indicate the simulation results from Ref. [47] for Dα (q, t ) and the dashed lines for Ωα (q, t ).

Then, the characteristic time scales of ψα (q, s), Dα (q, s), and Ωα (q, s) are different from each other. In order to see such differences numerically, it is convenient to use the simulation results [47] on the Stillinger–Weber [48] binary mixtures A80 B20 . In Fig. 4 the simulation results for ψα (q, t ) and Fα (q, t )/Sα (q) are shown versus scaled time t /t0 at scaled temperature T = 0.714, where t0 = σ /v0 . In Fig. 5 the simulation results for Dα (q, t ) and Ωα (q, t ) are also shown versus scaled time t /t0 . Then, Ωα (q, t ) (or Fα (q, t )) shows three crossover times. One is a mean-free time τf from a ballistic motion to a caging process, the second is a β -relaxation time τβ from a caging process to a α -relaxation process, and the last is a α -relaxation time τα from a α -relaxation process to a long-time diffusion process, where τα ≫ τβ ≫ τf . Then,

M. Tokuyama / Physica A 395 (2014) 31–47

43

ψα (q, t ) has a minimum around τf . Dα (q, t ) has a maximum for t < τf and decreases rapidly until τβ . It decreases slowly until τα and becomes constant for t ≫ τα . Hence τf is the characteristic time for ψα (q, t ), τβ for Dα (q, t ), and τα for Ωα (q, t ). In order to check whether the denominator Fα (q, t ) in Eq. (85) can be replaced by Fα (q, 0)(=Sα (q)) or not since τα ≫ τf , in Fig. 4 the simulation results for ψα (q, t )(Fα (q, t )/Sα (q)) are also shown versus scaled time t /t0 . They show that ψα (q, t )(Fα (q, t )/Sα (q)) coincides with ψα (q, t ). Thus, Fα (q, t ) can be safely replaced by Fα (q, 0) in Eq. (85). Hence this suggests that the time t dependence of the slower variables are negligible if they are multiplied by ψα (q, t ). This situation is generally true for lower temperatures near the glass transition since τα ≫ τf . By solving Eq. (84) formally, one can find  t gα (q, s)/Fα (q, s)ds. (93) Qα ρα (q, t ) = iqFα (q, t ) 0

Then, this leads to

⟨Qα ρα (q, t )gα (q, 0)∗ ⟩ = iqFα (q, t )Dα (q, t ).

(94)

Hence the time scale of Qα ρα (q, t ) turns out to be much larger than that of ψα (q, t ). As discussed above, therefore, the second term of Eq. (90) can be safely replaced by iqψα (q, t )Qα ρα (q, 0) = 0. Thus, Eq. (90) reduces to d dt

gα (q, t ) ≃ Qα iLgα (q, t ).

(95)

This is easily solved to give the formal solution gα (q, t ) = etQα iL Jα (q, 0),

(96)

which satisfies the orthogonality condition given by ⟨gα (q, t )ρα (q, 0) ⟩ = 0. Thus, one can approximately write the memory function ψα (q, t ) given by Eq. (85) as ∗

ψα (q, t ) ≃

⟨gα (q, t )gα (q, 0)∗ ⟩ , Sα (q)

(97)

where gα (q, t ) is given by Eq. (96) and the denominator of Eq. (85) is safely replaced by Sα (q). We now discuss how to derive an equation for the memory function ψα (q, t ) given by Eq. (97). By taking gα (q, t ) as A(t ), we introduce the projection operator ℘α′ by

℘α′ B(t ) =

⟨B(t )gα (q, 0)∗ ⟩ gα (q, 0), 2 vth

(98)

2 where ⟨gα (q, 0)gα (q, 0)∗ ⟩ = ⟨Jα (q, 0)Jα (q, 0)∗ ⟩ = vth . Since gα (q, t ) obeys Eq. (95), one can also use Eq. (13) to obtain the generalized Langevin equations for gα (q, t ), where iL is now all replaced by Qα iL. As shown in Fig. 4, however, ψα (q, t ) has a negative minimum around the mean-free time τf . This situation is the same as that for any other temperatures. Hence one cannot use type [T] here because as shown in Eq. (15) the memory function contains the correlation function ⟨gα (q, t )gα (q, 0)∗ ⟩ in the denominator, which becomes zero at a time shorter than τf . In order to derive a linear equation for gα (q, t ), therefore, one may use type [M] here. From Eq. (13), one then obtains

∂ gα (q, t ) = − ∂t

t



ϕα′ (q, s)gα (q, t − s)ds + fα′ (q, t )

0

(99)

with the memory function 2 ϕα′ (q, t ) = ⟨fα′ (q, t )fα′ (q, 0)∗ ⟩/vth ,

where the fluctuating force

fα′ (q, t )

(100)

is given by

′

fα′ (q, t ) = etiL iL′ gα (q, 0).

(101)

℘α′ . ∗

fα′ (q, t )

= 1− The fluctuating force is uncorrelated only Here iL = and iL gα (q, 0) = Qα iLJα (q, 0), where with the linear functions of ρα (q, 0) and gα (q, 0) as ⟨fα′ (q, t )ρα (q, 0) ⟩ = ⟨fα′ (q, t )gα (q, 0)∗ ⟩ = 0. From Eq. (99), we thus obtain ′

′

Qα′ Qα iL

∂ ψα (q, t ) = − ∂t

 0

Qα′

t

ϕα′ (q, s)ψα (q, t − s)ds.

(102)

Eqs. (87) and (102) are coupled equations to obtain the solutions for the scattering functions Fα (q, t ). In order to solve Eq. (102) formally, it is convenient to take the Laplace transform of Eq. (102). From (102), we then obtain

ψα [q, z ] =

2 vth 1 . Sα (q) z + ϕα′ [q, z ]

(103)

44

M. Tokuyama / Physica A 395 (2014) 31–47

From Eqs. (91) and (92), we also find

Ωα [q, z ] = Dα [q, z ]/z = ψα [q, z ]/z 2 .

(104)

3.3. Nonlinear memory function In addition to the randomly-varying processes, the fluctuating force fα′ (q, t ) still contains nonlinear terms of density fluctuations since those are the slowly-varying variables. This situation is exactly the same as that of fα (q, t ). Hence it can be formally split up into two parts, the random force Rα (q, t ) and the nonlinear fluctuating force ξα′ (q, t ): fα′ (q, t ) = ξα′ (q, t ) + Rα (q, t ),

(105)

where Rα (q, t ) satisfies Eq. (35). Then, the memory function ϕα′ (q, t ) can be written as

ϕα′ (q, t ) = 2γα δ(t ) + 1ϕα′ (q, t )

(106)

with the nonlinear memory term 2 1ϕα′ (q, t ) = ⟨ξα′ (q, t )ξα′ (q, 0)∗ ⟩/vth .

(107)

1ϕα′ (q, t )

Next, we show that the nonlinear memory function exactly reduces to 1ϕα (q, t ) within the MCT binary approximation. In fact, by using ℘2α , one can find the operator identities ′

′

′

etiL = etiL ℘2α + etiL Q2α , ′

′

etiL Q2α = etQ2α iL Q2α +

(108) t



e(t −s)iL ℘2α iL′ esQ2α iL Q2α ds. ′

′

(109)

0

From Eqs. (101) and (105), one then obtains ′

ξα′ (q, t ) = etiL ℘2α [Qα iLJα (q)], Rα (q, t ) = e

tQ2α iL′

(110)

Q2α [Qα iLJα (q)],

(111)

where the second term of Eq. (109) becomes zero since ⟨Rα (q, t )ρα (q, 0) ⟩ = ⟨Rα (q, t )gα (q, 0) ⟩ = 0. We now calculate the fluctuating force ξα′ (q, t ). Similarly to Eqs. (50) and (51), use of Eqs. (38) and (39) leads to ∗

℘2c [Qc iLJc (q)] = ℘2s [Qs iLJs (q)] =

2 vth



2iN 1/2

vc (q, k )ρc (k )ρc (q − k ),

(112)

k

2  vth

iN 1/2

∗

vs (q, k )ρc (k )ρs (q − k ).

(113)

k

From Eq. (110), we then obtain

ξα′ (q, t ) = ξα(0) (q, t ),

(114) tiL′

where we have used the fact that for arbitrary nonlinear function h(ρα ) of ρα (k , 0), one can find e h(ρα (k , 0)) = etiL h(ρα (k , 0)) = h(ρα (k , t )). Thus, the nonlinear memory function 1ϕα′ (q, t ) exactly coincides with 1φα (q, t ). In the following, therefore, we just write ϕα′ (q, t ) as ϕα (q, t ) within the MCT binary approximation. By using Eq. (103), one can also write the MCT solutions given by Eq. (30) as Fα [q, z ] Sα (q)

=

1 z + q2 ψα [q, z ]

.

(115)

This leads to

∂ Fα (q, t ) = −q2 ∂t

t



ψα (q, t − s)Fα (q, s)ds,

(116)

0

where ψα (q, t ) obeys Eq. (102) within the MCT approximation. Eq. (116) should be compared to Eq. (87). Thus, the main difference between TMCT and MCT is that the memory term in Eq. (87) is convolutionless in time, while that in Eq. (116) is convolution in time. Hence the short-time solution of MCT is given by Eq. (32), while that of TMCT is given by Eq. (89), 2 where Ω 2 = q2 vth /Sα (q). The former leads to α2 (t = 0) = −2/3, while the latter to α2 (t = 0) = 0.

M. Tokuyama / Physica A 395 (2014) 31–47

45

3.4. Debye–Waller factor We next show that Eq. (87) also has non-zero solutions fα (q) at long times. By using Eq. (104), one can write Eq. (88) as



fα (q) = lim exp −q2



dzezt

t →∞

ψα [q, z ]



z2

.

(117)

In the long time limit, use of Eqs. (60) and (103) then leads to lim q2

ψα [q, z ] z

z →0

2 q2 vth

=

z ϕ[q, z ]Sα (q)

=

1

Γα (q)

.

(118)

From Eq. (117), we thus obtain

 fα (q) = exp −



1

Γα (q)

.

(119)

One can also solve Eq. (119) for Γα (q) as

Γα (q) = −

1 ln[fα (q)]

.

(120)

Eqs. (119) and (120) have the same forms as those of Eqs. (78) and (79), respectively. Similarly to the AMCT results, therefore, the TMCT results are also expected to make Tc lower than that of MCT. Thus, the present formulation is expected to overcome the problems (i)–(iii) safely. 3.5. Single-particle dynamics We finally discuss the single-particle dynamics. Use of Eqs. (64) and (103) then leads to

ψs [q, z ] = ψs(0) [z ] − q2 ψs(2) [z ] + · · ·

(121)

with

ψs(0) [z ] =

2 vth

(0)

z + γs + 1ϕs [z ]

,

(122)

2 ψs(2) [z ] = 1ϕs(2) [z ]ψs(0) [z ]2 /vth .

(123)

By using Eqs. (1), (87), (97), and (121), we obtain, up to O(q2 ), d dt

M2 (t ) = 6D(t ) = 6

t



ψs(0) (s)ds,

(124)

0

where D(t ) denotes the time-dependent self-diffusion coefficient. The formal solution of Eq. (124) leads to M2 (t ) = 6

t



(t − τ )ψs(0) (τ )dτ ,

(125)

0

M2 [z ] = 6ψs(0) [z ]/z 2 .

(126)

The long-time self-diffusion coefficient D is then given by ∞

 D=

ψs(0) (τ )dτ =

0

2 vth

γs + 1ϕs(0) [0]

.

(127)

This form coincides with Eq. (72), although their temperature dependence is different from each other. Similarly, one can also find, up to O(q4 ),

α2 (t ) =

72 M 2 (t )2

t



(t − τ )ψs(2) (τ )dτ .

(128)

0

As is easily shown, Eq. (128) leads to α2 (t = 0) = 0. Thus, the present formulation does not have the second problem (ii). We next compare the single-particle dynamics of TMCT with those of MCT. By expanding Eq. (115) in powers of q and using Eqs. (122) and (123), one can rewrite Eqs. (70) and (71) as M2 [z ] = 6ψs(0) [z ]/z 2 ,

(129)

  M4 [z ] = 120 ψs(2) [z ]/z 2 + M2 [z ]ψs(0) [z ]/(6z ) ,

(130)

46

M. Tokuyama / Physica A 395 (2014) 31–47

respectively. Here we note that although Eqs. (126) and (127) have the same forms as Eqs. (129) and (72), respectively, their numerical values are different from each other since the starting equations for Fs (q, t ) are different from each other. The inverse Laplace transform of Eq. (130) is given by M4 (t ) = 120

t



  1 (t − τ ) ψs(2) (τ ) + M2 (τ )ψs(0) (τ ) dτ . 6

0

(131)

From Eq. (131), one then obtains for MCT

α2 (t ) =

72

t



M2 (t )2



(2)

1

(0)



(t − τ ) ψs (τ ) + M2 (τ )ψs (τ ) dτ − 1. 6

0

(132)

This should be compared with Eq. (128) for TMCT. As is easily shown, Eq. (132) leads to α2 (t = 0) = −2/3. Finally, we should mention that although the memory functions ψα (q, t ) and ϕα (q, t ) have the same forms in MCT and TMCT, their numerical results are completely different from each other since Fα (q, t ) obeys different equations given by Eqs. (87) and (116). 4. Summary The main purpose of the present paper was to overcome three difficulties (i)–(iii), which appear in the ideal MCT. This was done by using the time-convolutionless projection-operator method and by taking the density fluctuations as relevant variables. Thus, we have first derived the time-convolutionless equations for Fα (q, t )

∂ Fα (q, t ) = −q2 ∂t

t



ψα (q, s)dsFα (q, t ),

(133)

0

and then derived the time-convolution equations for ψα (q, t )

∂ ψα (q, t ) = − ∂t

t



ϕα (q, s)ψα (q, t − s)ds.

(134)

0

Then, we have shown that the memory function ϕα (q, t ) coincides with that obtained by MCT within the MCT binary approximation. These coupled equations are valid to describe not only the dynamics of liquids but also that of supercooled liquids as long as the systems are in equilibrium. On the other hand, one needs to derive new equations to study the dynamics of glasses because the systems are out of equilibrium there. By using Eq. (134), we have also rewritten the MCT equations as

∂ Fα (q, t ) = −q2 ∂t

t



ψα (q, t − s)Fα (q, s)ds.

(135)

0

Here Eqs. (133) differ from Eqs. (135), where the memory term in the former is convolutionless in time, while that in the latter convolution in time. In this paper, we have shown that this difference leads to quite different results. Both equations predict that the transport coefficients obey singular functions as long as the system is in equilibrium. This situation is exactly the same as that in critical phenomena, except that the crossover from an equilibrium state to a nonequilibrium state occurs at the deviation temperature Tn which is higher than the critical temperature Tc . Thus, the first result appears in the critical temperature. In fact, the Debye–Waller factor fα (q) is given from Eq. (133) by fα (q) = exp[−1/Γα (q)],

(136)

while Eq. (135) leads to fα (q) = Γα (q)/[1 + Γα (q)].

(137)

Hence Eqs. (136) and (137) suggest that the critical temperature Tc of TMCT is lower than that of MCT. As discussed in Section 2.4, this has already been seen in the alternative MCT equations, in which the memory term is convolutionless in time. The second result is that Eq. (133) leads to α2 (t = 0) = 0, while Eq. (135) leads to α2 (t = 0) = −2/3. Since the long-time dynamics of TMCT is qualitatively the same as that of AMCT, the singular temperature must be identical to that of the simulations as shown in Section 2.4. Thus, the problems (i)–(iii) have been shown to be naturally solved. Finally, we note that one needs to solve the TMCT equations numerically for glass-forming liquids to check whether d0 has the same value as that of simulations or not. This will be discussed elsewhere. Acknowledgments The author wishes to thank H. Mori and I. Oppenheim for valuable discussions and also for their hospitality and continual encouragement. He also thanks S. Enda, Y. Kimura, and T. Narumi for preparing the simulation data by using the SGI Altix3700Bx2 in Advanced Fluid Information Research Center, Institute of Fluid Science, Tohoku University. This work was partially supported by High Efficiency Rare Elements Extraction Technology Area, Institute of Multidisciplinary Research for Advanced Materials (IMRAM), Tohoku University and also by Institute of Fluid Science, Tohoku University.

M. Tokuyama / Physica A 395 (2014) 31–47

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] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48]

U. Bengtzelius, W. Götze, A. Sjölander, J. Phys. C 17 (1984) 5915. E. Leutheusser, Phys. Rev. A 29 (1984) 2765. W. Götze, in: J.P. Hansen, D. Levesque, J. Zinn-Justen (Eds.), Liquids, Freezing and Glass Transition, North-Holland, Amsterdam, 1991. M. Fuchs, W. Götze, I. Hofacker, A. Latz, J. Phys.: Condens. Matter 3 (1991) 5047. M. Nauroth, W. Kob, Phys. Rev. E 55 (1997) 657. M.H. Nauroth, Ph.D. Thesis, Institute für Theoretische Physik der Technischen Universität München, 1999. W. Götze, Th. Voigtmann, Phys. Rev. E 67 (2003) 021502. G. Foffi, W. Götze, F. Sciortino, P. Tartaglia, T. Voigtmann, Phys. Rev. E 69 (2004) 011505. E. Flenner, G. Szamel, Phys. Rev. E 72 (2005) 031508. Th. Voigtmann, J. Horbach, Europhys. Lett. 74 (2006) 459. T. Narumi, M. Tokuyama, Phys. Rev. E 84 (2011) 022501. M. Domschke, M. Marsilius, T. Blochowicz, Th. Voigtmann, Phys. Rev. E 84 (2011) 031506. J.L. Barrat, W. Götze, A. Latz, J. Phys.: Condens. Matter 1 (1989) 7163. W. van Megen, S.M. Underwood, Phys. Rev. Lett. 70 (1993) 2766. P. Lunkenheimer, A. Pimenov, A. Loidl, Phys. Rev. Lett. 78 (1997) 2995. T. Franosch, M. Fuchs, W. Götze, M.R. Mayr, A.P. Singh, Phys. Rev. E 55 (1997) 7153. M. Fuchs, W. Götze, M.R. Mayr, Phys. Rev. E 58 (1998) 3384. W. van Megen, T.C. Mortensen, S.R. Williams, J. Müller, Phys. Rev. E 58 (1998) 6073. C. Beck, W. Härtl, R. Hempelmann, J. Chem. Phys. 111 (1999) 8209. L. Fabbian, A. Latz, R. Schilling, F. Sciortino, P. Tartaglia, C. Theis, Phys. Rev. E 60 (1999) 5768. C. Theis, F. Sciortino, A. Latz, R. Schilling, P. Tartaglia, Phys. Rev. E 62 (2000) 1856. E. Bartsch, T. Eckert, C. Pies, H. Sillescu, J. Non-Cryst. Solids 307–310 (2002) 802. T. Eckert, E. Bartsch, Faraday Discuss. 123 (2003) 51. S.P. Das, Rev. Modern Phys. 76 (2004) 785. K. Binder, W. Kob, Glassy Materials and Disordered Solids, World Scientific, Singapore, 2005. X.J. Han, H. Teichler, Phys. Rev. E 75 (2007) 061501. W. Götze, Complex Dynamics of Glass Forming Liquids: A Mode Coupling Theory, Oxford Science, Oxford, 2009. A. Rahman, K.S. Singwi, A. Sjölander, Phys. Rev. 126 (1962) 986. M. Tokuyama, Phys. Rev. E 80 (2009) 031503. M. Tokuyama, Y. Kimura, Physica A 387 (2008) 4749. M. Tokuyama, Physica A 387 (2008) 1926. M. Tokuyama, Physica A 389 (2010) 951. R. Kubo, J. Math. Phys. 4 (1963) 174. H. Mori, Progr. Theoret. Phys. 33 (1965) 423. M. Tokuyama, H. Mori, Progr. Theoret. Phys. 55 (1976) 411. M.S. Green, J. Chem. Phys. 20 (1952) 1281. 22 (1954) 398. K. Kawasaki, Ann. Phys. 61 (1970) 1. H. Mori, H. Fujisaka, Progr. Theoret. Phys. 49 (1973) 764. T. Koide, T. Kodama, Phys. Rev. E 78 (2008) 051107. M. Tokuyama, Phys. Rev. E 62 (2000) R5915. H. Mori, M. Tokuyama, T. Morita, Progr. Theoret. Phys. Suppl. 64 (1978) 50. T. Narumi, M. Tokuyama, Rep. Inst. Fluid Sci. 19 (2007) 73. T. Narumi, Ph.D’s Thesis, Tohoku University, 2008. M. Tokuyama, Phys. Rev. E 82 (2010) 041501. M. Tokuyama, J. Chem. Phys. B 115 (2011) 14030. M. Tokuyama, AIP Conf. Proc. 1518 (2013) 47. M. Tokuyama, S. Enda, Physica A 392 (2013) 2999. T.A. Weber, F.H. Stillinger, Phys. Rev. B 31 (1985) 1954.

47