Mauro Prencipe Dip. Scienze della Terra – Università di Torino Impact of Crystallography on Modern Science Accademia delle Scienze Torino, 25 giugno 2014

Accademia delle Scienze - Torino 25/06/2014 Outline of the presentation  What we know about the structure of the Earth  How do we know…  Problems and possible solutions  Contributions from theoretical methods  Some results

Accademia delle Scienze - Torino 25/06/2014 What do we know about the structure of the Earth?

Accademia delle Scienze - Torino 25/06/2014 What do we know about the structure of the Earth?  Average chemical composition  Mineralogical composition  Pressure and Temperature  Physical properties  Dynamics Olivine (Mg2 SiO4 ) monoclinic Wadsleite (Mg2 SiO4 ) orthorhombic Ringwoodite (Mg2 SiO4 ) cubic (spinel structure) Mg-perovskite (MgSiO3 ) + periclase (MgO) Increasing pressure

Accademia delle Scienze - Torino 25/06/2014 How do we are reasonably sure that the model is true? Force of gravity Total mass of the Earth (5.972·1024 kg) Dimension (radius: 6371 km) Chemical composition (from some meteorites reflecting the average composition of the early solar nebula) Average density: Mass/Volume (5.54 g/cm3)

Accademia delle Scienze - Torino 25/06/2014 How do we are reasonably sure that the model is true? Mineralogy and mineral’s stability (thermodynamics, experimental petrology) Dynamic phenomena Magnetic field Less direct sources of information

Accademia delle Scienze - Torino 25/06/2014 One of the most important sources of information… Earthquakes

Accademia delle Scienze - Torino 25/06/2014 Seismic Tomography = Vp : velocity of longitudinal waves B : bulk modulus ρ : density B and ρ do depend by pressure and temperature, and hence by the depth Vp is a function of depth; more in general, it is a function of position

Accademia delle Scienze - Torino 25/06/2014 Seismic Stations What we see…

Accademia delle Scienze - Torino 25/06/2014 Questions: Number of layers? Refraction indices? Depth of each discontinuity? …Therefore, we need a prior model… Detector line Discontinuity #1a, 1b Discontinuity #2a ray ray a ray b n1a n2a n1b Different models can lead to the same arrival point at the detector line… Warning! Q P

Seismic Tomography Experimental measurements; Quantum Mechanical Calculations Thermo-Elastic Properties Models of the Inner Earth Direct sources of Information Accademia delle Scienze - Torino 25/06/2014 How to build and check the model From the macroscopic to the microscopic scale Very accurate and precise at high P and low T, or at low P and high T. What about at HP and HT?

) , ( ) , ( R x E R x H    Hamiltonian Operator Wave Function, depending upon the electronic spatial and spin coordinates (x) and nuclear (R) Energy of the system ) , ( s r x  r: space coordinates; s: spin coordinates of all the electrons The key is in the time independent Schrödinger’s equation: Accademia delle Scienze - Torino 25/06/2014 Theoretical methods

                 k i i j i h k h hk k h ij ik k n i i k N k k r e Z Z r e r e Z m M H , , , 2 2 2 1 2 2 2 1 2 2 1 2 1 2 1 2   Kinetic energy of the nuclei Kinetic energy of the electrons Electrons-nuclei potential energy Interelectronic potential energy Internuclear potential energy The Hamiltonian for a multi-electronic (n electrons), multi-nuclear (N nuclei) system, reads: Accademia delle Scienze - Torino 25/06/2014 Theoretical methods

    h k h hk k h el r e Z Z R E R E , 2 2 1 ) ( ) ( ) (R E NN The total energy, at a fixed nuclear configuration R, is given by the sum of the electronic energy and of the internuclear potential ENN (R). The E(R) function defines the Born-Oppenheimer surface, whose minima (R0 ) identify stable or metastable phases (at least at T=0K). Born-Oppenheimer approximation The nuclear motion is generally solved within the harmonic approximation, in a small neighbourhood of R0 Accademia delle Scienze - Torino 25/06/2014

ON 1 2 N. Points Precision Geometry Optimization Phonons Step size Geometry Editor Input Cell Space Group Atoms Rem. Atoms Subst. Atoms Insert Coins Cost: too much Output HF B88 PW91 LDA None LYP VWN %HF Exchange Integrals 2 3 Shrink 4 6 Convergence Integration grid Fock Mixing Basis set quality RUN OFF ON Pseudopotentials Exchange Correlation Omphacite P2/n L H 0.001 0.1 10-5 10-10 10-6 10-10 LGRID XXLGRID 0% 80% 0% 100% VDZ V6Z* Off

Accademia delle Scienze - Torino 25/06/2014 The Quasi-harmonic approximation It is assumed that the frequency of vibration of each normal mode does not depend upon the temperature of the crystal, but does depend upon the volume. For each normal mode, a Grüneisen parameter () can be defined. V V i i i        The kinetic energy of the nuclei is obtained within the limit of the harmonic approximation: at any cell volume V, the nuclear displacements from the equilibrium positions are small; in these conditions, the nuclear motion can be described by a set of independent harmonic oscillators, each one being characterized by a frequency ν and an energy :    h T n T E         2 1 ) , ( ) ( Vibrational quantum number: number of phonons of frequency ν, at the temperature T, following the Bose- Einstein statistics.

Accademia delle Scienze - Torino 25/06/2014 Reliability of the calculated vibrational frequencies Low albite (NaAlSi3 O8 ) Triclinic, space group P-1 Calculated Raman spectrum (WC1LYP Hamiltonian) v.s. the experimental one An important ingredient for the calculation of thermal effects on various properties is the vibrational spectrum of the crystal Raman spectrum

Thermal and zero point pressures in the quasi-harmonic approximation                   j kT V vib vib vib T vib vib vib vib vib j Z kT T V F V F T V P T V TS T V U T V F / ) ( e , lnZ ) , ( ) , ( ) , ( ) , ( ) , (  Given the classical relation between F and P, and the statistical definition of F in terms of the partition function Z: Derivation of Fvib (V,T), with respect to V                j j j j j j j T st T n V h V h V E T V P      ) , ( 2 ) , ( 0 Total Pressure Static Pressure (Pst ) Zero Point Pressure (Ppz ) Thermal Pressure (Pth ) Number of phonons of the jth oscillator, at the temperature T Accademia delle Scienze - Torino 25/06/2014

Negative static pressure field Beryl (Be6 Al4 Si12 O36 ) WC1LYP Hamiltonian Pst (GPa) Est (hartree) V (Å3) Static pressure: no thermal and zero point effects V V E V P st st d ) ( d ) (   Static pressures (Pst ) at a given cell volumes (V) can easily be calculated as first derivatives, with respect to V, of the static energy curve Est (V) Accademia delle Scienze - Torino 25/06/2014

Thermal and zero point pressures: the quasi-harmonic approximation A whole vibrational spectrum is evaluated for a number of values of the unit cell volume in a given range, and the Grüneisen’s parameters are then estimated. Beryl, WC1LYP Hamiltonian. Distribution of Grüneisen’s parameters Accademia delle Scienze - Torino 25/06/2014

T (K) 0 100 200 300 400 500 P st , P th (GPa) -1.3 -1.2 0.0 Ppz Zero Point pressure (right P axis) Thermal pressure (left P axis) Static pressure (left P axis) 1.25 1.30 GPa Negative thermal expansion Positive thermal expansion Total pressure and thermal expansion Case of Beryl (WC1LYP) Thermal expansion against a zero external pressure. To keep the total pressure (P) of the phonon gas P=Pst +Pzp +Pth at the constant value of 0GPa, being Pzp almost constant in the investigated T range, Pst = -Pth At T lower than about 280K, Pst must increase with T; this is realized through a contraction of the unit cell volume: negative thermal expansion. The reverse is true at higher temperatures. Accademia delle Scienze - Torino 25/06/2014

Thermal expansion Case of beryl Thermal expansion coefficient (), as a function of temperature, at three different pressures. Accademia delle Scienze - Torino 25/06/2014

Compressibility Once the total pressure (P) has been obtained as a function of the cell volume (V) and of the temperature (T), the (in)compressibility can be calculated as a function of both P and T Beryl (WC1LYP calculation) CALC EXP V0 692.1 676.8 B0 180 179 B’ 4.0 3.7 Slope at different pressures (GPa/K) Accademia delle Scienze - Torino 25/06/2014 Volume (V) in Å3; bulk modulus (B) in GPa BT P 0 3 6 9 (dBT /dT)P -0.020 -0.022 -0.028 -0.048

Accademia delle Scienze - Torino 25/06/2014 Mineral name Chemical formula Bcalc Bexp Beryl Be6 Al4 Si12 O36 180 179 Diamond C 445 445 Diopside CaMgSi2 O6 106 106 Jadeite NaAlSi2 O6 131 134 Omphacite (ord.) (Na0.5 Ca0.5 )(Mg0.5 Al0.5 )Si2 O6 112 116 Low Albite NaAlSi3 O8 54 54 Mg-perovskite MgSiO3 241 244 Calcite CaCO3 76 77 Aragonite CaCO3 65 65 Periclase MgO 160 160-164 Some results (not only beryl…) Experimental (Bexp ) and calculated bulk moduli (Bcalc ) of several minerals, at T=300K (P=0GPa). All values are in GPa Hamiltonian: hybrid HF/DFT (WC1LYP)

Accademia delle Scienze - Torino 25/06/2014 Diamond Bulk modulus at T=300K Exp: 445(1) Gpa Calculated: 445 GPa

Accademia delle Scienze - Torino 25/06/2014 Conclusions Ab Initio hybrid HF/DFT methods can provide very reliable estimation of thermo-elastic parameters of crystalline solids, which are useful in:  the construction of a model for the inner Earth (or other bodies);  the understanding of mechanisms of compression and thermal expansion. Thank You for Your Attention!