Molecular Dynamics simulation of glass structures J.-M. Delaye 1 - PowerPoint PPT Presentation

Molecular Dynamics simulation of glass structures J.-M. Delaye 1 with the contributions of T. Charpentier 2 , L.-H. Kieu 1 , F. Pacaud 1 , M. Salanne 3 1 Service d’Etudes de Vitrification et procédés hautes Températures (SEVT), CEA Marcoule, France 2 Nanosciences et Innovation pour les Matériaux la Biomédecine et l'Énergie (NIMBE), CEA Saclay, France 3 Physicochimie des Electrolytes et Nanosystèmes interfaciaux (PHENIX), Université Pierre et Marie Curie, France Joint ICTP – IAEA Workshop November 10, 2017 CEA | 10 AVRIL 2012 | PAGE 1 6-10 November 2017, Trieste, Italy

OUTLINE Objective of this lecture: Describe the state of the art about simulation of silicate glass structure by classical molecular dynamics � Some fundamentals about classical molecular dynamics (MD) and interatomic potentials (15’) � Alumino silicate glasses: three examples (10’) � Boro silicate glasses: two examples (10’) � Conclusions � Perspectives: some words about very recent approaches (Reaxff, machine learning) (5’) PAGE 2

QUICK INTRODUCTION TO THE CLASSICAL MOLECULAR DYNAMICS � Classical molecular dynamics is able to represent the dynamics of several thousands of atoms (from 1000 to 10 6 or more) during several picoseconds (10ps - 1µs) � Computers are more and more powerful. � Larger systems are simulated g i v i n g a c c e s s t o n e w mechanisms: mechanical properties, longer relaxations … � 10 6 atoms < > cubic box of 25 nm of side (SiO 2 with some water molecules) November 10, 2017 PAGE 3

QUICK INTRODUCTION TO THE CLASSICAL MOLECULAR DYNAMICS � Representation of the atomic interactions → INTERATOMIC POTENTIALS � The first models were quite simple At short distance, interpenetration of the electronic clouds: repulsion At large distance, coulombic and dipolar attraction (dispersion term) 12 6 ⎛ ⎞ q q r C ⎛− ⎞ ⎛ ⎞ ⎛ ⎞ σ σ ⎜ ⎟ ( ) i j ij ij ( r ) B exp ⎜ ⎟ r 4 ⎜ ⎟ ⎜ ⎟ φ = + − φ = ε − ij ij ij ⎜ ⎟ ⎜ ⎟ 6 ⎜ ⎟ ⎜ ⎟ r r r r ρ ⎝ ⎠ ⎝ ⎠ ⎝ ⎠ ij ij ij ij ij ⎝ ⎠ Lennard-Jones potentials Buckingham potentials PAGE 4

QUICK INTRODUCTION TO THE CLASSICAL MOLECULAR DYNAMICS � Oxide glasses are also subjected to covalent interactions → orbital hybridation, developed local angular order O Tetrahedral order around the Si and Al ions (SiO 4 - AlO 4 ) [4] Si O Tetrahedral or triangular order around the B ions (BO 4 or [4] B [3] B BO 3 ) � This chemical property is taken into account using three body (angular) potentials Triplet <jik>: r ij , r ik , θ jik γ γ ( ) 2 r , r , exp( )(cos cos ) φ θ = λ + θ − θ 3 ij ik jik jik 0 r r r r − − ij c ik c A. Tilocca, N.H. de Leeuw, A.N. Cormack, Phys. Rev. B 73 (2006) 104209 Energy associated to the triplet PAGE 5 M.-T. Ha, S.H. Garofalini, J. Am. Ceram. Soc. 100 (2017) 563

QUICK INTRODUCTION TO THE CLASSICAL MOLECULAR DYNAMICS � The ions can have dipolar momenta → the polarisation is treated via shell models or polarisation terms Shell models The ions are represented as cores (massive) connected to charged shells (very small mass) representing the valence electrons The shift between the core and the shell creates a dipole. The core and its shell are connected by an harmonic spring potential In the silicate glasses, the shell model is used only to represent the O ions because the polarisability is larger for this species compared to the other ones. PAGE 6 B.G. Dick Jr., A.W. Overhauser, Physical Review 112 (1958) 90

QUICK INTRODUCTION TO THE CLASSICAL MOLECULAR DYNAMICS Polarizable potentials (PIM, AIM) � Charge ¡: ¡ V V V V V = + + + tot charge disp rep pol ! " * * + ! = $%&'() ' !# " ! !" > � Polarisa4on ¡ ¡: ¡ � Répulsion ¡: ¡ q r r q ⎡ ⎤ ⋅ µ µ ⋅ i ij j ij i ij j ij V g ( r ) g ( r ) ∑ = − ⎢ ⎥ !" # pol ij ij $ !" r 3 r 3 ) ( & ! " !" = ⎢ ⎥ i , j i ⎣ ⎦ ij ij > #&' !% " ! > ( r )( r ) 3 ⎡ ⎤ µ ⋅ µ ⋅ µ ⋅ µ i j ij i ij j ∑ + − ⎢ ⎥ � Dispersion ¡: ¡ 3 5 r r ⎢ ⎥ i , j i ⎣ ⎦ ij ij > !" !" * * ' $ 2 !" !" . - ,) + % - ,) + $ µ ! = ( + i % " ∑ α &!'( % !" $ !" + % $ ) ) % " 2 !# " ! & # !" !" > i i At each time step, the dipolar momenta are determined by minimizing the polarisation energy V pol . PAGE 7 F. Pacaud, J.-M. Delaye, T. Charpentier, L. Cormier, M. Salanne, J. Chem. Phys. 147 (2017) 161711

QUICK INTRODUCTION TO THE CLASSICAL MOLECULAR DYNAMICS � From the interatomic potentials to the forces: E ∂ ( ) ( ) E r r , r , ∑ ∑ tot = φ + φ θ F = − tot ij 3 ij ik jik i r ∂ i , j i , j , k i Total energy Force exerted on an atom � From the forces to the atomic displacements: discretization in time of the Newton’s equation ( r ) 2 d r ∂ φ ij i F m = ∑ − = i i 2 r dt ∂ j 2 i F ( ) t t b ( ) t δ ( ) 3 4 ( ) ( ) ( ) r t t r t v t t t t + δ = + δ + + δ + Ο δ m 2 6 ( ) F t ( ) ( ) ( ) ( ) 2 4 r t t 2 r t r t t t t + δ = − − δ + δ + Ο δ m 2 F ( ) t t b ( ) t δ ( ) 3 4 ( ) ( ) ( ) r t t r t v t t t t − δ = − δ + − δ + Ο δ m 2 6 November 10, 2017 Timestep = 1fs typically Positions are known at time t Calculation of the atomic displacements Calculation of the forces Positions are known at time t+ δ t PAGE 8

GLASS PREPARATION � A glass is prepared by equilibrating a liquid and by quenching it at ambient temperature Equilibration of the liquid Thermal quench Final relaxation to determine the equilibrium volume � The difficulty is to obtain a structure that corresponds to the minimum of the potential energy: one possibility is to prepare a glass in two steps PAGE 9

GLASS PREPARATION First step: using different initial densities, the potential energy - density curve is plotted (no modification of the volume during the glass preparation) Equilibration of the liquid Thermal quench Final relaxation Glass: 67,73%SiO 2 .18,04%B 2 O 3 .14,23%Na 2 O (mol%) The minimum of the potential energy and the corresponding density are determined PAGE 10 L. Deng, J. Du, Journal of Non-Crytalline Solids, 453 (2016) 177

GLASS PREPARATION Second step: using the equilibrium density, a complete glass preparation is performed using the complete thermal scheme Equilibration of the liquid Thermal quench Final relaxation to determine the equilibrium volume The initial density is taken equal to 1.05 the equilibrium density. The final configuration is an equilibrated glass ready for the analysis PAGE 11

POTENTIAL FIT � This step is very important: all the simulated glass properties (structural, dynamical, mechanical …) depend on the interatomic potentials � Two main methods to fit the interatomic potentials → fit on experimental data → fit on ab initio calculations (more and more popular) Reference configurations are prepared by ab initio methods The adjustable parameters are determined by a DFT: Optimisation of Classical MD: determination convergence loop. the wave functions - of the forces (and dipoles) determination of the with the trial interatomic forces and dipoles potentials Too large Shift between DFT and MD? N a 2 classique abinitio F F ∑ − 1 k 1 Small enough = Ξ = F N N * N a 2 abinitio a b F ∑ PAGE 12 Fitting procedure is finished k 1 =

ALUMINO SILICATE GLASSES: GUILLOT AND SATOR’S POTENTIALS � A work has been done to simulate complex systems that can be found in the Earth’s mantle (SiO 2 , TiO 2 , Al 2 O 3 , FeO, Fe 2 O 3 , MgO, CaO, Na 2 O, K 2 O) Buckingham type potential q O = -0.945 (from Matsui ’ s potential) The adjustable parameters have been fitted in order to reproduce the densities of 11 natural silicate melts - PAGE 13 B. Guillot, N. Sator, Geochim. and Cosmochim. Acta 71 (2007) 1249

ALUMINO SILICATE GLASSES: GUILLOT AND SATOR’S POTENTIALS � The interatomic potentials have been used to investigate other properties Molar volume vs. modifier concentration Distribution of local coordinations Na and Ca contributions to the electrical conductivity PAGE 14

ALUMINO SILICATE GLASSES SIMULATED BY AIM POTENTIALS (Y. ISHII ET AL.) � AIM (Aspherical Ion Model) potentials are used to simulate sodium alumino silicate glasses: ion shapes are introduced in the potentials With δσ i a deviation from the ionic radius and ν i the distorsion of the dipolar shape � Validation of the potentials on different structural characteristics: S(Q) Neutron X-Ray structural structural factors factors PAGE 15 Y. Ishii, M. Salanne, T. Charpentier, K. Shiraki, K. Kasahara, N. Ohtori, J. Phys. Chem. C 120 (2016) 24370