Prediction of crystal properties by numerical simulation

Last updated

The prediction of crystal properties by numerical simulation has become commonplace in the last 20 years as computers have grown more powerful and theoretical techniques more sophisticated. High accuracy prediction of elastic, electronic, transport and phase properties is possible with modern methods.

Contents

Ab Initio Calculations

Ab initio or first principles calculations are any of a number of software packages making use of density functional theory to solve for the quantum mechanical state of a system. Perfect crystals are an ideal subject for such calculations because of their high periodicity. Since every simulation package will vary in the details of its algorithms and implementations, this page will focus on a methodological overview.

Basic theory

Density functional theory seeks to solve for an approximate form of the electronic density of a system. In general, atoms are split into ionic cores and valence electrons. The ionic cores (nuclei plus non-bonding electrons) are assumed to be stable and are treated as a single object. Each valence electron is treated separately. Thus, for example, a Lithium atom is treated as two bodies – Li+ and e- – while oxygen is treated as three bodies, namely O2+ and 2e.

The “true” ground state of a crystal system is generally unsolvable. However, the variational theorem assures us that any guess as to the electronic state function of a system will overestimate the ground state energy. Thus, by beginning with a suitably parametrized guess and minimizing the energy with respect to each of those parameters, an extremely accurate prediction may be made. The question as to what one's initial guess should be is a topic of active research. [1]

In the large majority of crystal systems, electronic relaxation times are orders of magnitude shorter than ionic relaxation times. Thus, an iterative scheme is adopted. First, the ions are considered fixed and the electronic state is relaxed by considering the ionic and electron-electron pair potentials. Next, the electronic states are considered fixed and the ions are allowed to move under the influence of the electronic and ion-ion pair potentials. When the decrease in energy between two iterative steps is sufficiently small, the structure of the crystal is considered solved.

Boundary conditions

A key choice that must be made is how many atoms to explicitly include in one's calculation. In Big-O notation, calculations general scale as O(N3) where N is the number of combined ions and valence electrons. [2] For structure calculations, it is generally desirable to choose the smallest number of ions that can represent the structure. For example, NaCl is a bcc cubic structure. At a first guess, one might construct a cell of two interlocked cubes – 8 Na and 8 Cl – as one's unit cell. This will give the correct answer but is computationally wasteful. By choosing appropriate coordinates, one might simulate it with just two atoms: 1 Na and 1 Cl.

Crystal structure calculations rely on periodic boundary conditions. That is, the assumption is that the cell you have chosen is in the midst of an infinite lattice of identical cells. By taking our 1 Na 1 Cl cell and copying it many times along each of the crystal axes, we will have simulated the same superstructure as our 8 Na 8 Cl cell but with much reduced computational cost.

Raw output

Only a few lists of information will be output from a calculation, in general. For the ions, the position, velocity and net force on each ion are recorded at each step. For electrons, the guess as to the electronic state function may be recorded as well. Finally, the total energy of the system is recorded. From these three types of information, we may deduce a number of properties.

Calculable properties

Unit cell parameters

Unit cell parameters (a,b,c,α,β,γ) can be computed from the final relaxed positions of the ions. [3] In a NaCl calculation, the final position of the Na ion might be (0,0,0) in picometer Cartesian coordinates and the final position of the Cl ion might be (282,282,282). From this, we see that the lattice constant would be 584 pm. For non-orthorhombic systems, the determination of cell parameters might be more complicated, but many ab-initio numerical packages have utilities to make this calculation simpler.

Once the lattice cell parameters are known, patterns for single crystal or powder diffraction can be readily predicted via Bragg's Law. [4]

Temperature and pressure

The temperature of the system can be estimated by use of the Equipartition Theorem, with three degrees of freedom for each ion. Since ionic velocities are generally recorded at each step in the numerical simulation, the average kinetic energy of each ion is easy to calculate. There exist schemes which attempt to control the temperature of the simulation by, e.g. enforcing each ion to have exactly the kinetic energy predicted by the Equipartition Theorem (Berendsen thermostat) or by allowing the system to exchange energy and momentum with a (more massive) fictitious enclosing system (Nose-Hoover thermostat).

The net force on each ion is generally calculated explicitly at each numerical step. From this, the stress tensor of the system can be calculated and usually is calculated by the numerical package. By varying the convergence criteria, one can either seek a lowest energy structure or a structure that produces a desired stress tensor. Thus, high pressures can be simulated as easily as ambient pressures. [5]

Elastic properties

The Young's modulus of a mineral can be predicted by varying one cell parameter at a time and observing the evolution of the stress tensor. [6] Because the raw output of a simulation includes energy and volume, the integrated version of the Birch-Murnaghan equation of state is often used to determine bulk modulus.

Electronic density of states

The electronic density functional is explicitly used in the calculation of the electronic ground state. Packages such as VASP have an option to calculate the electronic density of states per eV to facilitate the prediction of conduction bands and band gaps. [7]

Thermal transport properties

The Green-Kubo relations can be used to calculate the thermal transport properties of a mineral. Since the velocities of the ions are stored at each numerical step, one can calculate the time correlation of later velocities with earlier velocities. The integral of these correlations is related to the Fourier thermal coefficient.

Diffusion

By recording the ionic positions at each time step, one can observe how far, on average, each ion has moved from its original position. [8] The mean squared displacement of each ion type is related to the diffusion coefficient for a particle undergoing Brownian motion.

Related Research Articles

Chemical bond Lasting attraction between atoms that enables the formation of chemical compounds

A chemical bond is a lasting attraction between atoms, ions or molecules that enables the formation of chemical compounds. The bond may result from the electrostatic force between oppositely charged ions as in ionic bonds or through the sharing of electrons as in covalent bonds. The strength of chemical bonds varies considerably; there are "strong bonds" or "primary bonds" such as covalent, ionic and metallic bonds, and "weak bonds" or "secondary bonds" such as dipole–dipole interactions, the London dispersion force and hydrogen bonding.

Computational chemistry is a branch of chemistry that uses computer simulation to assist in solving chemical problems. It uses methods of theoretical chemistry, incorporated into computer programs, to calculate the structures and properties of molecules, groups of molecules, and solids. It is essential because, apart from relatively recent results concerning the hydrogen molecular ion, the quantum many-body problem cannot be solved analytically, much less in closed form. While computational results normally complement the information obtained by chemical experiments, it can in some cases predict hitherto unobserved chemical phenomena. It is widely used in the design of new drugs and materials.

Ionic bonding Chemical bonding involving attraction between ions

Ionic bonding is a type of chemical bonding that involves the electrostatic attraction between oppositely charged ions, or between two atoms with sharply different electronegativities, and is the primary interaction occurring in ionic compounds. It is one of the main types of bonding along with covalent bonding and metallic bonding. Ions are atoms with an electrostatic charge. Atoms that gain electrons make negatively charged ions. Atoms that lose electrons make positively charged ions. This transfer of electrons is known as electrovalence in contrast to covalence. In the simplest case, the cation is a metal atom and the anion is a nonmetal atom, but these ions can be of a more complex nature, e.g. molecular ions like NH+
4
or SO2−
4
. In simpler words, an ionic bond results from the transfer of electrons from a metal to a non-metal in order to obtain a full valence shell for both atoms.

Molecular dynamics Computer simulations to discover and understand chemical properties

Molecular dynamics (MD) is a computer simulation method for analyzing the physical movements of atoms and molecules. The atoms and molecules are allowed to interact for a fixed period of time, giving a view of the dynamic "evolution" of the system. In the most common version, the trajectories of atoms and molecules are determined by numerically solving Newton's equations of motion for a system of interacting particles, where forces between the particles and their potential energies are often calculated using interatomic potentials or molecular mechanics force fields. The method is applied mostly in chemical physics, materials science, and biophysics.

In solid-state physics, the electronic band structure of a solid describes the range of energy levels that electrons may have within it, as well as the ranges of energy that they may not have.

In chemistry, a hypervalent molecule is a molecule that contains one or more main group elements apparently bearing more than eight electrons in their valence shells. Phosphorus pentachloride, sulfur hexafluoride, chlorine trifluoride, the chlorite ion, and the triiodide ion are examples of hypervalent molecules.

The extended Hückel method is a semiempirical quantum chemistry method, developed by Roald Hoffmann since 1963. It is based on the Hückel method but, while the original Hückel method only considers pi orbitals, the extended method also includes the sigma orbitals.

The Jahn–Teller effect is an important mechanism of spontaneous symmetry breaking in molecular and solid-state systems which has far-reaching consequences in different fields, and is responsible for a variety of phenomena in spectroscopy, stereochemistry, crystal chemistry, molecular and solid-state physics, and materials science. The effect is named for Hermann Arthur Jahn and Edward Teller, who first reported studies about it in 1937. The Jahn–Teller effect, and the related Renner–Teller effect, are discussed in Section 13.4 of the spectroscopy textbook by Bunker and Jensen.

Ionic radius, rion, is the radius of a monatomic ion in an ionic crystal structure. Although neither atoms nor ions have sharp boundaries, they are treated as if they were hard spheres with radii such that the sum of ionic radii of the cation and anion gives the distance between the ions in a crystal lattice. Ionic radii are typically given in units of either picometers (pm) or angstroms (Å), with 1 Å = 100 pm. Typical values range from 31 pm (0.3 Å) to over 200 pm (2 Å).

SIESTA (computer program)

SIESTA is an original method and its computer program implementation, to perform efficient electronic structure calculations and ab initio molecular dynamics simulations of molecules and solids. SIESTA's efficiency stems from the use of strictly localized basis sets and from the implementation of linear-scaling algorithms which can be applied to suitable systems. A very important feature of the code is that its accuracy and cost can be tuned in a wide range, from quick exploratory calculations to highly accurate simulations matching the quality of other approaches, such as plane-wave and all-electron methods.

Semi-empirical quantum chemistry methods are based on the Hartree–Fock formalism, but make many approximations and obtain some parameters from empirical data. They are very important in computational chemistry for treating large molecules where the full Hartree–Fock method without the approximations is too expensive. The use of empirical parameters appears to allow some inclusion of electron correlation effects into the methods.

Ab initio quantum chemistry methods are computational chemistry methods based on quantum chemistry. The term ab initio was first used in quantum chemistry by Robert Parr and coworkers, including David Craig in a semiempirical study on the excited states of benzene. The background is described by Parr. Ab initio means "from first principles" or "from the beginning", implying that the only inputs into an ab initio calculation are physical constants. Ab initio quantum chemistry methods attempt to solve the electronic Schrödinger equation given the positions of the nuclei and the number of electrons in order to yield useful information such as electron densities, energies and other properties of the system. The ability to run these calculations has enabled theoretical chemists to solve a range of problems and their importance is highlighted by the awarding of the Nobel prize to John Pople and Walter Kohn.

Stopping power (particle radiation)

In nuclear and materials physics, stopping power is the retarding force acting on charged particles, typically alpha and beta particles, due to interaction with matter, resulting in loss of particle energy. Its application is important in areas such as radiation protection, ion implantation and nuclear medicine.

Car–Parrinello molecular dynamics or CPMD refers to either a method used in molecular dynamics or the computational chemistry software package used to implement this method.

Atomistix ToolKit (ATK) is a commercial software for atomic-scale modeling and simulation of nanosystems. The software was originally developed by Atomistix A/S, and was later acquired by QuantumWise following the Atomistix bankruptcy. QuantumWise was then acquired by Synopsys in 2017.

The hybrid QM/MM approach is a molecular simulation method that combines the strengths of ab initio QM calculations (accuracy) and MM (speed) approaches, thus allowing for the study of chemical processes in solution and in proteins. The QM/MM approach was introduced in the 1976 paper of Warshel and Levitt. They, along with Martin Karplus, won the 2013 Nobel Prize in Chemistry for "the development of multiscale models for complex chemical systems".

CP2K

CP2K is a freely available (GPL) quantum chemistry and solid state physics program package, written in Fortran 2008, to perform atomistic simulations of solid state, liquid, molecular, periodic, material, crystal, and biological systems. It provides a general framework for different methods: density functional theory (DFT) using a mixed Gaussian and plane waves approach (GPW) via LDA, GGA, MP2, or RPA levels of theory, classical pair and many-body potentials, semi-empirical and tight-binding Hamiltonians, as well as Quantum Mechanics/Molecular Mechanics (QM/MM) hybrid schemes relying on the Gaussian Expansion of the Electrostatic Potential (GEEP). The Gaussian and Augmented Plane Waves method (GAPW) as an extension of the GPW method allows for all-electron calculations. CP2K can do simulations of molecular dynamics, metadynamics, Monte Carlo, Ehrenfest dynamics, vibrational analysis, core level spectroscopy, energy minimization, and transition state optimization using NEB or dimer method.

Quantemol Ltd is based in University College London initiated by Professor Jonathan Tennyson FRS and Dr. Daniel Brown in 2004. The company initially developed a unique software tool, Quantemol-N, which provides full accessibility to the highly sophisticated UK molecular R-matrix codes, used to model electron polyatomic molecule interactions. Since then Quantemol has widened to further types of simulation, with plasmas and industrial plasma tools, in Quantemol-VT in 2013 and launched in 2016 a sustainable database Quantemol-DB, representing the chemical and radiative transport properties of a wide range of plasmas.

Computational materials science and engineering uses modeling, simulation, theory, and informatics to understand materials. The main goals include discovering new materials, determining material behavior and mechanisms, explaining experiments, and exploring materials theories. It is analogous to computational chemistry and computational biology as an increasingly important subfield of materials science.

FHI-aims Molecular dynamics modelling software

FHI-aims is a software package for computational molecular and materials science written in Fortran. It uses density functional theory and many-body perturbation theory to simulate chemical and physical properties of atoms, molecules, nanostructures, solids, and surfaces. Originally developed at the Fritz Haber Institute in Berlin the ongoing development of the FHI-aims source code is now driven by a world-wide community of collaborating research institutions.

References

  1. Gonze, Xavier; Finocchi, Fabio (2004). "Pseudopotentials Plane Waves–Projector Augmented Waves: A Primer". Physica Scripta. T109: 40. Bibcode:2004PhST..109...40G. doi: 10.1238/Physica.Topical.109a00040 .
  2. Kresse, G.; Furthmüller, J. (July 1996). "Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set". Computational Materials Science. 6 (1): 15–50. doi:10.1016/0927-0256(96)00008-0.
  3. Dag, Sefa; Wang, Lin-Wang (13 May 2010). "Packing Structure of Poly(3-hexylthiophene) Crystal: Ab Initio and Molecular Dynamics Studies". The Journal of Physical Chemistry B. 114 (18): 5997–6000. doi:10.1021/jp1008219.
  4. Pagliai, Marco; Muniz-Miranda, Francesco; Cardini, Gianni; Righini, Roberto; Schettino, Vincenzo (May 2011). "Spectroscopic properties with a combined approach of ab initio molecular dynamics and wavelet analysis". Journal of Molecular Structure. 993 (1–3): 438–442. Bibcode:2011JMoSt.993..438P. doi:10.1016/j.molstruc.2011.02.007.
  5. Wentzcovitch, Renata M.; Price, G.David (1996). "High pressure studies of Mantle minerals by ab initio variable cell shape molecular dynamics". Molecular Engineering. 6 (1–2). doi:10.1007/BF00161722.
  6. Ono, Shigeaki (10 October 2013). "Elastic Properties of CaSiO3 Perovskite from ab initio Molecular Dynamics". Entropy. 15 (10): 4300–4309. Bibcode:2013Entrp..15.4300O. doi: 10.3390/e15104300 .
  7. Feng, Min; Yang, Aria; Zuo, Xu; Vittoria, Carmine; Harris, Vincent G. (2010). "Ab initio study on copper ferrite". Journal of Applied Physics. 107 (9): 09A521. Bibcode:2010JAP...107iA521F. doi:10.1063/1.3338905.
  8. Zhang, Yi; Zhao, Yusheng; Chen, Changfeng (29 April 2013). "study of the stabilities of and mechanism of superionic transport in lithium-rich antiperovskites". Physical Review B. 87 (13). Bibcode:2013PhRvB..87m4303Z. doi: 10.1103/PhysRevB.87.134303 .