- Email: [email protected]
The dynamics of supercooled silica: acoustic modes and boson peak
Journal of Non-Crystalline Solids 235±237 (1998) 320±324
The dynamics of supercooled silica: acoustic modes and boson peak J urgen Horbach, Walter Kob *, Kurt Binder Institut f ur Physik, Johannes-Gutenberg-Universit at, Staudinger Weg 7, D-55099 Mainz, Germany
Abstract Using molecular dynamics computer simulations we investigate the dynamics of supercooled silica in the frequency ÿ1 . We ®nd that for small wave-vectors the dispersion relations range 0.5±20 THz and the wave-vector range 0.13±3.5 A are in very good agreement with the ones found in experiments. In particular we discuss the wave-vector dependence of the longitudinal and transverse acoustic modes and of the boson peak. Ó 1998 Elsevier Science B.V. All rights reserved. PACS: 61.43.Fs; 61.20.Lc; 02.70.Ns; 64.70.Pf
1. Introduction In the last few years a signi®cant eort was undertaken to understand the nature of the so-called boson peak, a prominent dynamical feature at around 1 THz which is observed in strong glassformers [1±4]. Various theoretical approaches have been proposed to explain this peak, such as localized vibrational modes or scattering of acoustic modes, but so far no clear picture has emerged yet. Recently also computer simulations have been used in order to gain insight into the mechanism that gives rise to this peak, but due to the high cooling rates with which the samples were prepared (on the order of 1012 K/s) and small system the results of these investigations sizes (20±40 A) were not able to give a ®nal answer either [2,4]. In the present work we present the results of a large scale computer simulation of supercooled sil-
* Corresponding author. Tel.: 49 6131 393641; fax: 49 6131 395441; e-mail: [email protected].
0022-3093/98/$19.00 Ó 1998 Elsevier Science B.V. All rights reserved. PII: S 0 0 2 2 - 3 0 9 3 ( 9 8 ) 0 0 5 9 3 - 6
ica. At the temperature investigated we are still able to fully equilibrate the sample, thus avoiding the problem of the high cooling rates, and by using a large system we can minimize the possibility of ®nite size eects [5]. Thus, by making a large computational eort, we are able to study the dynamics of this strong glassformer in a frequency and wave-vector range which is not accessible to real experiments and can therefore investigate the properties of the boson peak in greater detail than was possible so far. 2. Details of the simulations The silica model we use for our simulation is the one proposed by van Beest et al. [6] and has been shown to give a good description of the static properties of silica glass [7] as well as the dynamical properties of the supercooled melt, such as the activation energies of the diusion constants [8]. In this model the potential between ions i and j is given by
J. Horbach et al. / Journal of Non-Crystalline Solids 235±237 (1998) 320±324
/ rij
qi qj e 2 Cij Aij eÿBij rij ÿ 6 : rij rij
1
The values of the parameters qi , Aij , Bij and Cij can be found in the original publication [6]. The nonCoulombic part of the potential was truncated The simulations were done and shifted at 5.5 A. at constant volume and the density of the system was ®xed to 2.37 g/cm3 . The system size was 3, 8016 ions, giving a size of the box of (48.37 A) and the equations of motion were integrated over 4 106 times steps of 1.6 fs, thus over a time span of 6.4 ns. This time is suciently long to fully equilibrate the system at 2900 K, the temperature considered in this study [8]. More details on the simulations can be found in Ref. [9]. 3. Results In the present work we study the dynamics of the system by means of JL q; m and JT q; m, the longitudinal and transverse current±current correlation functions for wave-vector q at frequency m [10]. These are de®ned as the longitudinal and transverse part of the current±current correlation function, i.e.
Ja q; m N X
ÿ1
321
Z1 dt exp i2pmt ÿ1
huk t ul 0 exp iq rk t ÿ rl 0i;
2
kl
where uk t is equal to q r_ k t=q for a L and equal to q r_ k t=q for a T . JL q; m is directly connected to the dynamical structure factor S q; m via the expression JL q; m
m2 S q; m; q2
3
whereas no similar relation exists for JT q; m. Note that S q; m can be measured in scattering experiments. In the following we will focus on the silicon±silicon correlation only, but we have found that the oxygen±oxygen correlation function behave very similarly. In Fig. 1 we show the frequency dependence of ÿ1 , the JL q; m for wave-vectors between 0.13 A smallest wave-vector compatible with our box, ÿ1 , a wave-vector which is at the location and 1.6 A of the ®rst sharp diraction peak in S q. From the ®gure we recognize that this correlation function has a peak at a frequency, mL q, which increases with increasing q and which corresponds to the
Fig. 1. Frequency dependence of the longitudinal current±current correlation for dierent wave-vectors q.
322
J. Horbach et al. / Journal of Non-Crystalline Solids 235±237 (1998) 320±324
Fig. 2. Frequency dependence of the dynamical structure factor for dierent wave-vectors q. The vertical line is the location of the boson peak at 1673 K as determined from the neutron scattering experiments [3].
longitudinal acoustic modes. For the larger waveÿ1 ) we observe that a plateau like vectors (>0.4 A region is formed for frequencies less than the corresponding acoustic mode peak. The latter feature is seen in the dynamical structure factor S q; m as the boson peak. In Fig. 2 S q; m is shown for sevÿ1 . Since S q; m is coneral values of q up to 1.6 A nected to JL q; m via the exact relation given by Eq. (3) it contains the same information as JL q; m, with the dierence, that the intensity at small frequencies is increased. Hence, although the acoustic modes appear for wave-vectors less ÿ1 as peaks, for larger values of q they than 0.22 A appear only as a plateau like feature, due to their overlap with the boson peak. Also included in Fig. 2 is the location of the boson peak as determined from neutron scattering experiments at 1673 K which is around 1.5 THz [3]. We will comment more on this work below. We also mention that a qualitatively similar behavior is observed for the transverse current±current correlation function [9]. The dispersion relations mL q and mT q are presented in Fig. 3 which are determined from JL q; m and JT q; m, respectively. We see that, as expected,
for small wave-vectors mL depends linearly on q. Also included in the ®gure is a line with slope cL 6370 m/s, the experimental value of the longitudinal sound velocity of silica at around 1600 K [3]. We see that the data points for mL for small q are very close to this line and thus we conclude that the sound velocity of this system is very close to the one of real silica, thus giving further support for the validity of the model potential. For the transverse acoustic modes the agreement between the experiment and the simulation data is a bit inferior, but still good. Also included in the ®gure is mBP , the location of the boson peak which is determined from S q; m. We ®nd that for large wave-vectors mBP is around 1.8 THz, a value that is a bit larger than the experimental value of 1.5 THz reported by Wischnewski et al. at 1673 K [3]. However, these authors also found that mBP increases with increasing temperature and a rough extrapolation of their data for mBP to T 2900 K shows that a value of 1.8 THz is quite reasonable, thus giving further support for the validity of our model. In the inset of Fig. 3 we show the dispersion curves at small values of q. We note that at small
J. Horbach et al. / Journal of Non-Crystalline Solids 235±237 (1998) 320±324
323
Fig. 3. Wave-vector dependence of mL (open circles), mT (®lled circles) and mBP (open triangles). The bold solid lines are the dispersion relations for the longitudinal and transverse acoustic modes (Ref. [3]). Inset: enlargement of the curves at small q.
wave-vectors the curve for mBP q seems to join smoothly the one for mL . Within the accuracy of our data it is not clear, whether the boson peak and the longitudinal acoustic mode become identical or whether the boson peak ceases to exist for wave-vectors smaller than approximately 0.2 ÿ1 . We also note that for wave-vectors larger A ÿ1 , a bit less than the location of the ®rst than 1.4 A sharp diraction peak, mL and mT do not increase anymore. The reason for this is likely the fact that at this q value the system has a quasi-Brillouin zone [2].
4. Conclusions To summarize we can say that our simulation allows to investigate the dynamics of supercooled silica in a wave-vector and frequency range which is not accessible to real experiments. We ®nd that the dispersion relations for the longitudinal and transverse acoustic modes agree very well with the experimental values. Furthermore we ®nd that the boson peak appears in the dynamical structure factor as well as in the longitudinal and transverse part of the current±current correlation function ÿ1 whereas for lowabove a q value of about 0.2 A
er values of q it cannot be seen in the latter quantities.
Acknowledgements We thank U. Buchenau, W. G otze, G. Ruocco and F. Sciortino for valuable discussions and the DFG, through SFB 262, and the BMBF, through grant 03N8008C, for ®nancial support.
References [1] U. Buchenau, M. Prager, N. N ucker, A.J. Dianoux, N. Ahmad, W.A. Phillips, Phys. Rev. B 34 (1986) 5665; P. Benassi, M. Krisch, C. Masciovecchio, V. Mazzacurati, G. Monaco, G. Ruocco, F. Sette, R. Verbeni, Phys. Rev. Lett. 77 (1996) 3835; M. Foret, E. Courtens, R. Vacher, J.-B. Suck, Phys. Rev. Lett. 77 (1996) 3831. [2] S.N. Taraskin, S.R. Elliott, Europhys. Lett. 39 (1997) 37; S.N. Taraskin, S.R. Elliott, preprint, 1997. [3] A. Wischnewski, U. Buchenau, A.J. Dianoux, W.A. Kamitakahara, J.L. Zarestky, preprint, 1997. [4] R. Dell'Anna, G. Ruocco, M. Sampoli, G. Viliani, Phys. Rev. Lett. 80 (1998) 1236. [5] J. Horbach, W. Kob, K. Binder, C.A. Angell, Phys. Rev. E 54 (1996) R5897.
324
J. Horbach et al. / Journal of Non-Crystalline Solids 235±237 (1998) 320±324
[6] B.W.H. van Beest, G.J. Kramer, R.A. van Santen, Phys. Rev. Lett. 64 (1990) 1955. [7] K. Vollmayr, W. Kob, K. Binder, Phys. Rev. B, 54 (1996) 15808; K. Vollmayr and W. Kob, Ber. Bunsenges. Phys. Chemie 100 (1996) 1399.
[8] J. Horbach, W. Kob, K. Binder, Philos. Mag. B 77 (1998) 297. [9] J. Horbach, W. Kob, unpublished. [10] J.P. Boon, S. Yip, Molecular Hydrodynamics, Dover, New York, 1980.
The dynamics of supercooled silica: acoustic modes and boson peak J urgen Horbach, Walter Kob *, Kurt Binder Institut f ur Physik, Johannes-Gutenberg-Universit at, Staudinger Weg 7, D-55099 Mainz, Germany
Abstract Using molecular dynamics computer simulations we investigate the dynamics of supercooled silica in the frequency ÿ1 . We ®nd that for small wave-vectors the dispersion relations range 0.5±20 THz and the wave-vector range 0.13±3.5 A are in very good agreement with the ones found in experiments. In particular we discuss the wave-vector dependence of the longitudinal and transverse acoustic modes and of the boson peak. Ó 1998 Elsevier Science B.V. All rights reserved. PACS: 61.43.Fs; 61.20.Lc; 02.70.Ns; 64.70.Pf
1. Introduction In the last few years a signi®cant eort was undertaken to understand the nature of the so-called boson peak, a prominent dynamical feature at around 1 THz which is observed in strong glassformers [1±4]. Various theoretical approaches have been proposed to explain this peak, such as localized vibrational modes or scattering of acoustic modes, but so far no clear picture has emerged yet. Recently also computer simulations have been used in order to gain insight into the mechanism that gives rise to this peak, but due to the high cooling rates with which the samples were prepared (on the order of 1012 K/s) and small system the results of these investigations sizes (20±40 A) were not able to give a ®nal answer either [2,4]. In the present work we present the results of a large scale computer simulation of supercooled sil-
* Corresponding author. Tel.: 49 6131 393641; fax: 49 6131 395441; e-mail: [email protected].
0022-3093/98/$19.00 Ó 1998 Elsevier Science B.V. All rights reserved. PII: S 0 0 2 2 - 3 0 9 3 ( 9 8 ) 0 0 5 9 3 - 6
ica. At the temperature investigated we are still able to fully equilibrate the sample, thus avoiding the problem of the high cooling rates, and by using a large system we can minimize the possibility of ®nite size eects [5]. Thus, by making a large computational eort, we are able to study the dynamics of this strong glassformer in a frequency and wave-vector range which is not accessible to real experiments and can therefore investigate the properties of the boson peak in greater detail than was possible so far. 2. Details of the simulations The silica model we use for our simulation is the one proposed by van Beest et al. [6] and has been shown to give a good description of the static properties of silica glass [7] as well as the dynamical properties of the supercooled melt, such as the activation energies of the diusion constants [8]. In this model the potential between ions i and j is given by
J. Horbach et al. / Journal of Non-Crystalline Solids 235±237 (1998) 320±324
/ rij
qi qj e 2 Cij Aij eÿBij rij ÿ 6 : rij rij
1
The values of the parameters qi , Aij , Bij and Cij can be found in the original publication [6]. The nonCoulombic part of the potential was truncated The simulations were done and shifted at 5.5 A. at constant volume and the density of the system was ®xed to 2.37 g/cm3 . The system size was 3, 8016 ions, giving a size of the box of (48.37 A) and the equations of motion were integrated over 4 106 times steps of 1.6 fs, thus over a time span of 6.4 ns. This time is suciently long to fully equilibrate the system at 2900 K, the temperature considered in this study [8]. More details on the simulations can be found in Ref. [9]. 3. Results In the present work we study the dynamics of the system by means of JL q; m and JT q; m, the longitudinal and transverse current±current correlation functions for wave-vector q at frequency m [10]. These are de®ned as the longitudinal and transverse part of the current±current correlation function, i.e.
Ja q; m N X
ÿ1
321
Z1 dt exp i2pmt ÿ1
huk t ul 0 exp iq rk t ÿ rl 0i;
2
kl
where uk t is equal to q r_ k t=q for a L and equal to q r_ k t=q for a T . JL q; m is directly connected to the dynamical structure factor S q; m via the expression JL q; m
m2 S q; m; q2
3
whereas no similar relation exists for JT q; m. Note that S q; m can be measured in scattering experiments. In the following we will focus on the silicon±silicon correlation only, but we have found that the oxygen±oxygen correlation function behave very similarly. In Fig. 1 we show the frequency dependence of ÿ1 , the JL q; m for wave-vectors between 0.13 A smallest wave-vector compatible with our box, ÿ1 , a wave-vector which is at the location and 1.6 A of the ®rst sharp diraction peak in S q. From the ®gure we recognize that this correlation function has a peak at a frequency, mL q, which increases with increasing q and which corresponds to the
Fig. 1. Frequency dependence of the longitudinal current±current correlation for dierent wave-vectors q.
322
J. Horbach et al. / Journal of Non-Crystalline Solids 235±237 (1998) 320±324
Fig. 2. Frequency dependence of the dynamical structure factor for dierent wave-vectors q. The vertical line is the location of the boson peak at 1673 K as determined from the neutron scattering experiments [3].
longitudinal acoustic modes. For the larger waveÿ1 ) we observe that a plateau like vectors (>0.4 A region is formed for frequencies less than the corresponding acoustic mode peak. The latter feature is seen in the dynamical structure factor S q; m as the boson peak. In Fig. 2 S q; m is shown for sevÿ1 . Since S q; m is coneral values of q up to 1.6 A nected to JL q; m via the exact relation given by Eq. (3) it contains the same information as JL q; m, with the dierence, that the intensity at small frequencies is increased. Hence, although the acoustic modes appear for wave-vectors less ÿ1 as peaks, for larger values of q they than 0.22 A appear only as a plateau like feature, due to their overlap with the boson peak. Also included in Fig. 2 is the location of the boson peak as determined from neutron scattering experiments at 1673 K which is around 1.5 THz [3]. We will comment more on this work below. We also mention that a qualitatively similar behavior is observed for the transverse current±current correlation function [9]. The dispersion relations mL q and mT q are presented in Fig. 3 which are determined from JL q; m and JT q; m, respectively. We see that, as expected,
for small wave-vectors mL depends linearly on q. Also included in the ®gure is a line with slope cL 6370 m/s, the experimental value of the longitudinal sound velocity of silica at around 1600 K [3]. We see that the data points for mL for small q are very close to this line and thus we conclude that the sound velocity of this system is very close to the one of real silica, thus giving further support for the validity of the model potential. For the transverse acoustic modes the agreement between the experiment and the simulation data is a bit inferior, but still good. Also included in the ®gure is mBP , the location of the boson peak which is determined from S q; m. We ®nd that for large wave-vectors mBP is around 1.8 THz, a value that is a bit larger than the experimental value of 1.5 THz reported by Wischnewski et al. at 1673 K [3]. However, these authors also found that mBP increases with increasing temperature and a rough extrapolation of their data for mBP to T 2900 K shows that a value of 1.8 THz is quite reasonable, thus giving further support for the validity of our model. In the inset of Fig. 3 we show the dispersion curves at small values of q. We note that at small
J. Horbach et al. / Journal of Non-Crystalline Solids 235±237 (1998) 320±324
323
Fig. 3. Wave-vector dependence of mL (open circles), mT (®lled circles) and mBP (open triangles). The bold solid lines are the dispersion relations for the longitudinal and transverse acoustic modes (Ref. [3]). Inset: enlargement of the curves at small q.
wave-vectors the curve for mBP q seems to join smoothly the one for mL . Within the accuracy of our data it is not clear, whether the boson peak and the longitudinal acoustic mode become identical or whether the boson peak ceases to exist for wave-vectors smaller than approximately 0.2 ÿ1 . We also note that for wave-vectors larger A ÿ1 , a bit less than the location of the ®rst than 1.4 A sharp diraction peak, mL and mT do not increase anymore. The reason for this is likely the fact that at this q value the system has a quasi-Brillouin zone [2].
4. Conclusions To summarize we can say that our simulation allows to investigate the dynamics of supercooled silica in a wave-vector and frequency range which is not accessible to real experiments. We ®nd that the dispersion relations for the longitudinal and transverse acoustic modes agree very well with the experimental values. Furthermore we ®nd that the boson peak appears in the dynamical structure factor as well as in the longitudinal and transverse part of the current±current correlation function ÿ1 whereas for lowabove a q value of about 0.2 A
er values of q it cannot be seen in the latter quantities.
Acknowledgements We thank U. Buchenau, W. G otze, G. Ruocco and F. Sciortino for valuable discussions and the DFG, through SFB 262, and the BMBF, through grant 03N8008C, for ®nancial support.
References [1] U. Buchenau, M. Prager, N. N ucker, A.J. Dianoux, N. Ahmad, W.A. Phillips, Phys. Rev. B 34 (1986) 5665; P. Benassi, M. Krisch, C. Masciovecchio, V. Mazzacurati, G. Monaco, G. Ruocco, F. Sette, R. Verbeni, Phys. Rev. Lett. 77 (1996) 3835; M. Foret, E. Courtens, R. Vacher, J.-B. Suck, Phys. Rev. Lett. 77 (1996) 3831. [2] S.N. Taraskin, S.R. Elliott, Europhys. Lett. 39 (1997) 37; S.N. Taraskin, S.R. Elliott, preprint, 1997. [3] A. Wischnewski, U. Buchenau, A.J. Dianoux, W.A. Kamitakahara, J.L. Zarestky, preprint, 1997. [4] R. Dell'Anna, G. Ruocco, M. Sampoli, G. Viliani, Phys. Rev. Lett. 80 (1998) 1236. [5] J. Horbach, W. Kob, K. Binder, C.A. Angell, Phys. Rev. E 54 (1996) R5897.
324
J. Horbach et al. / Journal of Non-Crystalline Solids 235±237 (1998) 320±324
[6] B.W.H. van Beest, G.J. Kramer, R.A. van Santen, Phys. Rev. Lett. 64 (1990) 1955. [7] K. Vollmayr, W. Kob, K. Binder, Phys. Rev. B, 54 (1996) 15808; K. Vollmayr and W. Kob, Ber. Bunsenges. Phys. Chemie 100 (1996) 1399.
[8] J. Horbach, W. Kob, K. Binder, Philos. Mag. B 77 (1998) 297. [9] J. Horbach, W. Kob, unpublished. [10] J.P. Boon, S. Yip, Molecular Hydrodynamics, Dover, New York, 1980.