Icosahedral order in liquid and glassy phases of cyclohexane
IIcosahedral order in liquid and glassy phases of cyclohexane
Tomoko Mizuguchi and Soichi Tatsumi and Susumu Fujiwara
Faculty of Materials Science and Engineering, Kyoto Institute of Technology, Matsugasaki,Sakyo-ku, Kyoto 606-8585, Japan
ARTICLE HISTORY
Compiled January 15, 2020
ABSTRACT
We performed all-atom molecular dynamics simulations for bulk cyclohexane andanalysed the short- and medium-range structures in supercooled and glassy statesby using the Voronoi tessellation technique. From the analyses of both the poten-tial energy of the system and the radial distribution function of molecules, cyclo-hexane was found to be vitrified as the temperature decreased. Furthermore, theicosahedral-like structures are dominant at all temperatures and grow in a super-cooled liquid, whereas the face-centred cubic structures do not grow when the tem-perature decreases. It was also ascertained that the icosahedral-like structure is moredominant than the full-icosahedral one. The network of the distorted icosahedronspreads throughout the system at low temperatures. Our simulation demonstratesthe stability of the icosahedral structure even in a non-spherical molecule such ascyclohexane.
KEYWORDS cyclohexane; supercooled and glassy state; Voronoi tessellation; icosahedral localorder
1. Introduction
In the glass transition process, the motion of particles slows down drastically as thetemperature is decreased. This is accompanied by the divergence of relaxation time orviscosity in a small range of temperatures. Although many experimental, theoretical,and computational techniques have been employed in an attempt to elucidate theorigin of glass transition [1–5], none of them has succeeded in explaining all abnormalbehaviours related to the glass transition.To avoid crystallization, many experiments on glass transition have been conductedon molecules with complex structures, such as polymers [6,7]. In recent years, advancesin experimental techniques have made it possible to vitrify even simple molecules suchas propene, benzene, and cyclohexane [8–11]. Interestingly, cyclohexane, which exhibitsthe glass transition by confinement in nanopores with a diameter of 1.9 to 2.9 nm,demonstrates a heat capacity anomaly caused by the first-order phase transition ata higher temperature than the glass transition point [10]. This means that certainstructural changes should occur in supercooled liquid cyclohexane, but the changesare small and thus difficult to capture experimentally.
CONTACT T. Mizuguchi. Email: [email protected] a r X i v : . [ c ond - m a t . d i s - nn ] J a n dditionally, many computer simulations of supercooled liquids and glasses havebeen done with simple isotropic models. With recent developments in computer tech-nology, it has become possible to realistically model and simulate realistic moleculesat the atomic level. We use an all-atom model of cyclohexane as a realistic model of amolecular liquid for analysis in supercooled and glassy systems.Cyclohexane is a small molecule, with the chemical formula is C H (see Fig. 1(a)),and it mainly interacts via the van der Waals (vdW) force with like molecules. The meltis transformed at 280 K into a plastic phase, which has a cubic cell, as well as at 186 Kinto a crystalline phase, which has a monoclinic cell [12]. However, crystallization canbe prevented by confinement within a nano-sized pore, as described above. The focusof this work is the short- and medium-range structures in supercooled cyclohexane. Weexamine the change in the local structures using molecular dynamics (MD) simulationsand show that icosahedral structures grow in the supercooled state.In Section 2, we explain the methods used in the MD simulations as well as thestructural analysis. In Section 3, we analyse the short- and medium-range structure ofsupercooled and glassy cyclohexane using Voronoi tessellation. In Section 4 we closewith concluding remarks.
2. Methods
We carried out all-atom MD simulations of bulk cyclohexane using the CHARMMDrude force field developed by Vorobyov et al. [13] in the NAMD2.10 program pack-age [14]. Periodic boundary conditions were applied to the MD unit cell which contains1372 cyclohexane molecules. The vdW interaction is represented by the Lennard-Jonespotential which was truncated by applying a switching function [15] with a range of10–12 ˚A. The long-range electrostatic interaction was calculated using the smoothparticle-mesh Ewald method [16]. All simulations were performed in the isothermal-isobaric ensemble with a time step of 1 fs. The pressure was maintained at 1 atm by theLangevin piston Nos´e-Hoover method [17,18] with barostat oscillation and dampingtime constants of 100 and 50 fs, respectively. The temperature was controlled with theLangevin dynamics at a damping coefficient of 1.0 ps − . The SHAKE algorithm [19]was used to fix the bond lengths involving the hydrogen atoms.First, we prepared liquid cyclohexane at 360 K and quenched the system to 10 Kat a cooling rate of 10 K/s. Starting from the configurations dumped during thecooling process, we conducted 5-ns runs for structural relaxation and subsequent 1-nsproduction runs at each temperature.
The intermolecular structures were analysed in terms of Voronoi polyhedra. First, wecalculated the centre of mass for each molecule from the MD snapshots. Next, theVoronoi diagram was computed for the centre-of-mass coordinates using a Voro++library [20]. Finally, for each Voronoi polyhedron, each face with an area less than 1%of the total surface area was removed in order to avoid overcounting [21]. We used theVoronoi index (cid:104) n , n , n , n (cid:105) to identify the type of each polyhedron, where n i is thenumber of faces with i vertices. To characterise the medium-range order structures,we define the cluster as the nearest neighbour shell obtained by Voronoi tessellation.2 a) (b) Figure 1. (colour online) (a) A cyclohexane molecule represented by balls and sticks. Cyan and white spheresrepresent carbon and hydrogen atoms, respectively. (b) An example of the connection of the molecular clusters.The sphere represents a cyclohexane molecule. Cyan, blue, and red spheres represent nearest-neighbor shells,central molecules, and shared molecules, respectively. The bond between spheres is not a real bond and issimply shown for clarity. g m o l ( r ) r [Å] T = 320 K T = 60 K decreasing T Figure 2. (colour online) Radial distribution functions between the centre of mass of cyclohexane moleculesat various temperatures ranging from 60 K to 320 K with 20 K step. Each line is shifted 0.1 units with respectto the next temperature for readability.
The network is defined as the connection of the clusters: two clusters are consideredto be connected if they share at least one molecule [22]. An example of the connectionof the clusters is shown in Fig. 1(b).
3. Results and Discussion
First, we calculated the radial distribution functions for the centre of mass of themolecules g mol ( r ) at each temperature. Figure 2 shows g mol ( r ) at various temperatures.When the temperature was lowered, the peak positions of g mol ( r ) shifted to smaller r values, reflecting the increase in the density of the system. The second peak of g mol ( r )has a small sub-peak around r = 12 ˚A at lower temperatures. This corresponds to thewell-known observation that the second peak of g mol ( r ) splits in a glassy state [23–25]. The split second peak of g mol ( r ) is considered to be related to the icosahedrallocal order, which can be seen in disordered systems composed of isotropic particlessuch as Lennard-Jones systems [23], bulk metallic glass [24], and random hard-spherepacking [26]. In our systems also, icosahedral local order was observed, as describedlater.Potential energy per molecule is shown in Fig. 3 as a function of temperature. Theslope changes around 110 K but there is no jump and thus it does not crystallise and4 T [K] E [ k ca l/ m o l ] Figure 3.
Temperature dependence of potential energy per molecule. Solid lines are fitting functions. instead forms a glass at low temperatures. We can evaluate T g from the figure as thepoint where the slope changes [24]. T g is determined by dividing the curve into twopieces and fitting a straight line to each; the intersection of the two lines is definedas T g . From Fig 3, T g was estimated to be approximately 112 K. This value is inagreement with the experimental value obtained from heat capacity measurements ofconfinement cyclohexane, in which T g was 80–100 K [10]. Note that the higher T g in simulations than experiments can be attributed to the faster cooling rate in thesimulations than that in experiments.To examine the local intermolecular structures, we performed Voronoi tessellationfor the centre-of-mass coordinates. Figure 4 shows the fractions of the ten most pop-ulous types of polyhedra with the Voronoi index (cid:104) n , n , n , n (cid:105) at T = 300 K, 200K, and 80 K. The (cid:104) , , , (cid:105) type of polyhedron is the most populous at all threetemperatures, although its fraction at T = 300 K is not much different from that ofother polyhedra. It is an icosahedral-like cluster and is frequently observed in metal-lic glasses [22,27–33]. The full icosahedron (cid:104) , , , (cid:105) was also observed; however, thefraction was smaller than that of (cid:104) , , , (cid:105) . This implies that an icosahedral-like clus-ter is preferred as a local structure rather than a full icosahedron in the supercooledand glassy state of cyclohexane.The plastic phase of cyclohexane has a face-centred cubic (fcc) lattice, whose Voronoipolyhedron is (cid:104) , , , (cid:105) . However, this polyhedron, as well as other crystal polyhedrasuch as (cid:104) , , , (cid:105) in a body-centred cubic lattice are not seen at all temperatures. Con-versely, we do see (cid:104) , , , (cid:105) polyhedra at all temperatures, which is considered to bean fcc-like cluster [31,34,35]. However, the fraction of (cid:104) , , , (cid:105) polyhedra is small andalmost constant at 2–3 % at all temperatures, and importantly it does not grow when5emperature decreases. The temperature dependence of the fractions of (cid:104) , , , (cid:105) , (cid:104) , , , (cid:105) and (cid:104) , , , (cid:105) polyhedra is shown in Figure 5(a). As seen in the figure,the icosahedral structures (cid:104) , , , (cid:105) and (cid:104) , , , (cid:105) grow in a supercooled liquid. Weshow the simplified diagrams of (cid:104) , , , (cid:105) and (cid:104) , , , (cid:105) clusters in Fig. 5(b). Theindex (cid:104) , , , (cid:105) represents a full icosahedron but the (cid:104) , , , (cid:105) cluster in our sys-tem is not necessarily a perfect icosahedron. Hirata et al. , reported that icosahedralclusters observed in Zr Pt metallic glass are not ideal icosahedra but are distorteddue to the geometric frustration of icosahedral order [28]. Additionally, in our system,the (cid:104) , , , (cid:105) cluster is distorted, and the (cid:104) , , , (cid:105) cluster has an icosahedral-likesymmetry but has one extra molecule (indicated by an arrow in Fig. 5(b)).To investigate the medium-range order, we calculated the average network size andthe number of networks of (cid:104) , , , (cid:105) clusters as shown in Figure 6. The tempera-ture dependence of the average network size is similar to that of the fraction of the (cid:104) , , , (cid:105) polyhedron (Fig. 5(a)). The number of networks becomes approximately 1below 200 K, meaning that the network spreads throughout the system. The spatialdistribution of (cid:104) , , , (cid:105) clusters is shown in Figure 7 as a projection on a two-dimensional plane. A cyclohexane molecule is represented by a sphere that is locatedat the centre-of-mass coordinate of the molecule. Red spheres represent the centremolecules in (cid:104) , , , (cid:105) clusters, and orange spheres are molecules that are composedof (cid:104) , , , (cid:105) clusters excluding the centre molecule. At 300 K, several networks existseparately. At 250 K, the size of each network becomes larger but some of them arestill separate. At 200 K, the networks grow and merge into nearly a single network. At100 K, the density of network increases and the (cid:104) , , , (cid:105) cluster spreads throughoutthe system.
4. Summary and Conclusions
We performed all-atom MD simulations of the liquid and glassy phases of cyclohex-ane as a model of a molecular liquid. The glassy state was achieved by quenchingfrom the melt and confirmed by the energy change and radial distribution function.The local structure was analysed using Voronoi tessellation and icosahedral struc-tures were found in the supercooled and glassy state. The dominant structure is anicosahedral-like (cid:104) , , , (cid:105) structure rather than a full-icosahedral (cid:104) , , , (cid:105) struc-ture. The (cid:104) , , , (cid:105) polyhedral network grows in three dimensions as the temperaturedecreases and spreads throughout the system. The icosahedral structure is consideredto be important in supercooled liquids and glasses because it is highly close-packed [36]and causes geometrical frustration owing to inconsistencies in the translational sym-metry. In fact, it has been observed in many metallic glasses. In this study, we foundthe icosahedral order in supercooled cyclohexane even though the shape of cyclohexaneis not necessarily spherical. This observation reinforces the idea that the icosahedralorder plays a role in glass formation. The relation between the experimentally ob-served transition in confined cyclohexane and its local structures are currently underinvestigation. Furthermore, an all-atom model can incorporate even intramolecularstructural changes, but it has a disadvantage for investigating slow dynamics nearthe glass transition point. Long-time simulations using a coarse-grained model are forfuture research. 6 V o r ono i I nd e x < n , n , n , n > fraction [%]300 K(a) fraction [%] <0,1,10,2><0,0,12,0><0,1,8,4><0,2,8,4><0,0,12,2><0,0,10,2><0,1,10,3><0,3,6,4><0,0,10,3><0,1,8,3,1>
200 K(b) fraction [%]80 K(c) <0,1,10,2><0,0,12,0><0,0,12,2><0,1,10,3><0,2,8,4><0,1,8,4><0,0,10,3><0,0,10,2><0,0,10,4><0,3,6,4>
Figure 4.
Fractions of the ten most populous Voronoi polyhedra at (a) 300 K, (b) 200 K, and (c) 80 K. <0,1,10,2><0,0,12,0><0,3,6,4> T [K] fr ac ti on [ % ] (a) (b) <0, 0, 12, 0> <0, 1, 10, 2> Figure 5. (colour online) (a) Temperature dependence of the fractions of (cid:104) , , , (cid:105) , (cid:104) , , , (cid:105) and (cid:104) , , , (cid:105) polyhedra. (b) Structures of (cid:104) , , , (cid:105) and (cid:104) , , , (cid:105) polyhedra. In the upper panels, only carbon atoms aredisplayed and the centre molecule is coloured red. In the lower panels, a cyclohexane molecule is representedby balls and the bond between balls is shown only for clarity (not real bonds). T [K] nu m b e r o f n e t w o r k s a v e r a g e n e t w o r k s i ze Figure 6. (colour online) Average network size (purple square, left scale) and the number of networks (greencircle, right scale) of (cid:104) , , , (cid:105) clusters as a function of temperature. a) T = 300 K (b) T = 250 K(c) T = 200 K (d) T = 100 K Figure 7. (colour online) Spatial distributions of (cid:104) , , , (cid:105) clusters at 300, 250, 200, and 100 K are shown asa projection on a two-dimensional plane. The centre-of-mass coordinate of a cyclohexane molecule is representedby a sphere. Red spheres represent the centre molecules in (cid:104) , , , (cid:105) clusters, and orange spheres are moleculesthat are composed of (cid:104) , , , (cid:105) clusters excluding the centre molecules. Others are shown in light grey. cknowledgements This research used computational resources of the Supercomputer Center, the Insti-tute for Solid State Physics, the University of Tokyo, and of the MEXT Joint Usage/ Research Center ”Center for Mathematical Modeling and Applications”, Meiji Uni-versity, Meiji Institute for Advanced Study of Mathematical Sciences (MIMS).
Funding
This work was partially supported by the Japan Society for the Promotion of ScienceKAKENHI grant no. JP19K05209.
References [1] Angell CA, Ngai KL, McKenna GB, et al. Relaxation in glassforming liquids and amor-phous solids. J Appl Phys. 2000;88:3113–3157.[2] Berthier L, Biroli G. Theoretical perspective on the glass transition and amorphous ma-terials. Rev Mod Phys. 2011;83:587–645.[3] Kob W. Computer simulations of supercooled liquids and glasses. J Phys: Condens Matter.1999;11:R85.[4] Tarjus G, Kivelson SA, Nussinov Z, et al. The frustration-based approach of supercooledliquids and the glass transition: a review and critical assessment. J Phys: Condens Matter.2005;17:R1143–R1182.[5] Cavagna A. Supercooled liquids for pedestrians. Physics Reports. 2009;476(4-6):51–124.[6] Cangialosi D. Dynamics and thermodynamics of polymer glasses. Journal of Physics Con-densed Matter. 2014;26(15):153101.[7] Roth CB, editor. Polymer glasses. CRC Press; 2016.[8] Tatsumi S, Aso S, Yamamuro O. Thermodynamic study of simple molecular glasses:Universal features in their heat capacity and the size of the cooperatively rearrangingregions. Phys Rev Lett. 2012;109:045701.[9] Xia Y, Dosseh G, Morineau D, et al. Phase diagram and glass transition of confinedbenzene. J Phys Chem B. 2006;110:19735.[10] Tatsumi S, Uehara T, Oguni M. A new phase transition and involved structural changeof confined cyclohexane. Netsu Sokutei. 2015;42(4):142–147. [in Japanese].[11] Yamamuro O, Matsuo T, Yamamuro NO, et al. Neutron diffraction and thermal studiesof amorphous cs2 realised by low-temperature vapour deposition. Europhys Lett. 2003;63:368–373.[12] Kahn R, Fourme R, Andr´e D, et al. Crystal Structures of Cyclohexane I and II. ActaCryst. 1973;B29:131.[13] Vorobyov I, Anisimov VM, Greene S, et al. Additive and classical drude polarizable forcefields for linear and cyclic ethers. Journal of Chemical Theory and Computation. 2007;3(3):1120–1133.[14] Phillips JC, Braun R, Wang W, et al. Scalable molecular dynamics with namd. J ComputChem. 2005;26(16):1781–1802.[15] Steinbach PJ, Brooks BR. New spherical-cutoff methods for long-range forces in macro-molecular simulation. J Comput Chem. 1994;15(7):667–683.[16] Essmann U, Perera L, Berkowitz ML, et al. A smooth particle mesh ewald method. JChem Phys. 1995;103(19):8577–8593.[17] Martyna GJ, Tobias DJ, Klein ML. Constant pressure molecular dynamics algorithms. JChem Phys. 1994;101(5):4177–4189.
18] Feller SE, Zhang Y, Pastor RW, et al. Constant pressure molecular dynamics simulation:The langevin piston method. J Chem Phys. 1995;103(11):4613–4621.[19] Ryckaert JP, Ciccotti G, Berendsen HJ. Numerical integration of the cartesian equationsof motion of a system with constraints: molecular dynamics of n-alkanes. J Comput Phys.1977;23(3):327–341.[20] Rycroft CH. Voro++: A three-dimensional voronoi cell library in c++. Chaos. 2009;19:041111.[21] Brostow W, Chybicki M, Laskowski R, et al. Voronoi polyhedra and Delaunay sim-plexes in the structural analysis of molecular-dynamics-simulated materials. Physical Re-view B. 1998;57(21):13448–13458. Available from: http://link.aps.org/doi/10.1103/PhysRevB.57.13448 .[22] Ward L, Miracle D, Windl W, et al. Structural evolution and kinetics in Cu-Zr metallicliquids from molecular dynamics simulations. Physical Review B. 2013;88(13):134205.[23] Kob W, Andersen HC. Testing made-coupling theory for a supercooled binary lennard-jones mixture i: The van hove correlation function. Phys Rev E. 1995;51(5):4626–4641.[24] Bailey NP, Schiøtz J, Jacobsen KW. Simulation of Cu-Mg metallic glass: Thermodynamicsand structure. Physical Review B. 2004;69(14):1–11.[25] Chen DZ, An Q, Goddard WA, et al. Ordering and dimensional crossovers in metallicglasses and liquids. Physical Review B. 2017;95(2):1–8.[26] Clarke AS, J´onsson H. Structural changes accompanying densification of random hard-sphere packings. Physical Review E. 1993;47(6):3975–3984.[27] Mauro NA, Wessels V, Bendert JC, et al. Short- and medium-range order in Zr80Pt20liquids. Physical Review B - Condensed Matter and Materials Physics. 2011;83(18):1–8.[28] Hirata A, Kang LJ, Fujita T, et al. Geometric frustration of icosahedron in metallicglasses. Science. 2013;341(6144):376–379. Available from: .[29] Zhang J, Chen C, Pei Q, et al. Ab initio molecular dynamics study of the lo-cal atomic structures in monatomic metallic liquid and glass. Materials & Design.2015;77:1–5. Available from: .[30] Yang MH, Li JH, Liu BX. Comparatively studying the local atomic structures of metallicglasses upon cyclic-loading by computer simulations. RSC Advances. 2017;7(30):18358–18365.[31] Shimono M, Onodera H. Dynamics and Geometry of Icosahedral Order in Liquid andGlassy Phases of Metallic Glasses. Metals. 2015;5(4):1163–1187. Available from: .[32] Fukunaga T, Itoh K, Otomo T, et al. Voronoi Analysis of the Structure of Ni-Zr-AlTernary Metallic Glass. Materials Transactions. 2007;48(7):1698–1702. Available from: http://joi.jlc.jst.go.jp/JST.JSTAGE/matertrans/MJ200750?from=CrossRef .[33] Ganesh P, Widom M. Ab initio simulations of geometrical frustration in supercooledliquid Fe and Fe-based metallic glass. Phys Rev E. 2008;77:014205.[34] Cape JN, Finney JL, Woodcock LV, et al. An analysis of crystallization by homogeneousnucleation in a 4000-atom soft-sphere model in a 4000-atom soft-sphere model. J ChemPhys. 1998;75(5):2366–2373.[35] Jiang SQ, Wu ZW, Li MZ. Effect of local structures on crystallization in deeply under-cooled metallic glass-forming liquids. Journal of Chemical Physics. 2016;144(15).[36] Frank FC, Kasper JS. Complex alloy structures regarded as sphere packings. i. definitionsand basic principles. Acta Cryst. 1958;11:184..[33] Ganesh P, Widom M. Ab initio simulations of geometrical frustration in supercooledliquid Fe and Fe-based metallic glass. Phys Rev E. 2008;77:014205.[34] Cape JN, Finney JL, Woodcock LV, et al. An analysis of crystallization by homogeneousnucleation in a 4000-atom soft-sphere model in a 4000-atom soft-sphere model. J ChemPhys. 1998;75(5):2366–2373.[35] Jiang SQ, Wu ZW, Li MZ. Effect of local structures on crystallization in deeply under-cooled metallic glass-forming liquids. Journal of Chemical Physics. 2016;144(15).[36] Frank FC, Kasper JS. Complex alloy structures regarded as sphere packings. i. definitionsand basic principles. Acta Cryst. 1958;11:184.