Application Notes

This application note provides an overview of the forcefield based simulation of crystalline C₆₀ (Buckminsterfullerene) using the LAMMPS molecular dynamics simulation package. The emphasis is on the overall philosophy of LAMMPS calculations in the MedeA® environment.

MEDEA's MT-Elastic Properties automates the calculation of elastic properties from first-principles.
For a given input system, MT applies symmetry-relevant strains and computes the resulting stress tensor for each deformed structure. The elastic properties are derived by a multi-dimensional least-square fit of the strain-stress data.

 MEDEA MT in Depth: Forsterite Mg₂SiO₄

The introduction of vacancies in a crystalline solid causes a local rearrangement of atoms around
the defect. The purpose of this case study is the prediction of such changes in α-quartz.

Structure around Oxygen Vacancy in α-Quartz

The electronic structures of the metallic and insulating phases of VO₂ are calculated using density functional theory and modern hybrid functionals as recently implemented in MedeA-VASP. Strongly contrasting previous calculations as based on local or semilocal approximations, which missed the insulating behavior of the low-temperature phase, these new calculations accurately reproduce the optical band gap and thus bring to an end a fifty-year old controversy on the origin of the metal-insulator transition of VO₂.

Accurate Band Gaps of Correlated Transition-Metal Oxides from MedeA-VASP with Hybrid-Functionals

We demonstrate the capabilities of MedeA®with selected examples, focusing on the lattice thermal
conductivity using forcefield methods as implemented in the software environment. The thermal conductivity is of high interest in different fields. In thermoelectrics, materials are sought with a high electrical conductivity combined with a low thermal conductivity as can be found in doped semiconductors with a high density of states near the band edges. In the present paper, we investigate the thermal conductivity of Si-Ge alloys and discuss the influence of defects, and disorder. All the computations are done using MedeA's LAMMPS and Thermal Conductivity modules, with the Reverse Non-Equilibrium Molecular Dynamics (RNEMD) approach.

Lattice thermal conductivity calculations in Si-Ge

In this note, we compare computed structural properties of clay minerals with experimental data, wherever such data is available. Calculated properties include cell parameters, atomic positions (in particular H positions) and internal surface areas.

Structural Properties And Surface Area of Clay Minerals

Patrick Soukiassian, Erich Wimmer, Edvige Celasco, Claudia Giallombardo, Simon Bonanni, Luca Vattuone, Letizia Savio, Anontio Tejeda, Mathieu Silly, Marie D’angelo, Fausto Sirotti, Mario Rocca

Nanostructuring a surface is a key and mandatory engineering step toward advances in nanotechnology. A team of french and italian scientists and of a franco-american company has just shown that hydrogen/deuterium (H/D) induces the opening of nanotunnels below the surface of an advanced semiconductor, silicon carbide (SiC). Such a finding is an especially interesting one, particularly in views of the specific properties of SiC. These investigations have been performed using advanced experimental tools such as synchrotron radiation and vibrational spectroscopy techniques, and state-of-the-art theoretical simulations. Depending on the H/D SiC surface exposures, these nanotunnels undergo through a sequence of semiconducting/metallic/semiconducting transitions. Such nanotunnels open very promising prospects toward applications in electronic, chemistry, storage, sensors and biotechnology.

Faujasitic zeolites are employed in a range of industrial contexts including separation processes, gas purification and dehydration, and shape selective catalysis. The simulation of the interaction of small molecules with faujasite based zeolites is increasingly practical and the resulting information provides the basis for efficient materials selection and design. When experimental data are scarce and empirical models insufficiently accurate, the MedeA®software environment can be employed to provide quantitative property data for faujasitic zeolite systems.

The colloidal behavior of asphaltenes explains the high viscosity of heavy oils and the occurrence of solid deposits in reservoir rocks, production wells or transport lines.
Asphaltenes make refining more difficult because they cause coke deposits and tend to lower the yield in high quality fuels or chemicals. Due to their significant content of sulfur and metals (Ni, V) and high aromaticity, energy-intensive conversion and hydrotreatment processes are required.

Colloidal behavior of confined model asphaltenes using molecular dynamics

The formation of micelles by surfactants was followed by molecular dynamics calculations performed with MedeA-LAMMPS and using the PCFF+ forcefield. An initial model with a random distribution of C9TAC surfactant molecules was built using the MedeA®building capabilities, such as the Molecular Builder and the Amorphous Materials Builder. Results are in agreement with previous simulation studies and available experimental data.

Typical HP-HT conditions may be defined by fluid pressures in excess of 50 MPa or temperatures above 150°C, as encountered in deep reservoirs below the North Sea, the Caspian Sea, the Gulf of Mexico and offshore Brazil, among others. Under these conditions natural gas may contain hydrocarbons with chain length of as much as 30 carbon atoms and a methane content higher than 60% (molar concentration). Because of the high temperature, water content may be significant, and further, H2S or CO2 contents may be elevated in some HP-HT reservoirs.

Properties of Natural Gases in Classical and High Pressure-High Temperature Conditions

The contact resistance between metals and semiconductors in nanoelectronic devices is mainly determined by the Schottky barrier. Controlling the Schottky barrier height (SBH) hence means being able to manipulate the contact resistance. and thereby to reduce the energy consumption as well as the heat production of electronic devices. While so far phenomenological considerations were able to determine the SBH only at a qualitative level, a parameter-free quantitative evaluation is possible via atomistic simulations using the MedeA®software platform.. This application note illustrates the calculation and modification of the SBH for a system of direct technological importance, namely a NiSi/Si contact. Furthermore, the effect of dopant atoms on the SBH, which are needed to tune the SBH for minimal contact resistance of n- and p-doped semiconductors, is investigated. Reduction of the Schottky barrier height is achieved by doping with Ba. This note also demonstrates how the preferred positions of dopant elements can be determined with S as example.

Prediction of Schottky Barrier in Electronic Devices

The effect of resin molecular architecture on the small strain elastic constants of diamine-cured epoxy thermosets has been studied using classical all-atom simulations conducted within the MedeA®simulation environment. Batches of thermoset systems have been created using chemically similar di, tri and tetrafunctional resins, followed by calculation of stiffness and compliance matrices for each individual model. Analysis of the batches of topologically and geometrically distinct structures using the Hill-Walpole approach yields upper and lower bounds estimates of the moduli differing by typically 1%, enabling critical comparison with experimentally-measured values.

Effect of Resin Molecular Architecture on Epoxy Thermoset Mechanical Properties

The energy band structure of InAs is computed with the HSE06 hybrid functional using MedeA-VASP. The calculations give a direct band gap of 0.35 eV, which coincides with the experimental value. The Γ-L separation in the lowest conduction band is computed to be 1.11 eV. Screened-exchange FLAPW calculations reported a value of 1.21 eV. The earlier handbook value was 0.74 eV. More recent experiments reported 1.1±0.05 eV, which is perfectly consistent with the present calculations.

Energy Band Structure of InAs with MedeA-VASP

The stress-strain behavior of a Cu nanowire is simulated with the MedeA®environment using a quasi-classical embedded atom potential and the LAMMPS molecular dynamics code. The monocrystalline wire has a diameter of 3.3 nm and is stretched in the [100] direction. The simulations show an initial elastic region with a linear increase in stress, which reaches a maximum just before the onset of slip in (111) planes. Upon further strain the model reveals the formation of more slip planes and necking until the break point is reached.

Stress-Strain Behavior of a Cu Nanowire Simulated with MedeA-LAMMPS

This application note describes the calculation of densities, cohesive energy densities (CED), and enthalpies of vaporization for a range of straight chain hydrocarbon fluids. The construction, simulation, and analysis methodologies employed are reviewed; and the accuracy of atomistic simulation for organic materials and polymers illustrated. Mean absolute errors are 0.23% and 0.28% for densities and heats of vaporization respectively for the PCFF+ forcefield.

Organic Materials Properties: Densities,  Cohesive Energies, and Heats of Vaporization

This application note provides an overview of the use of EAM based forcefield simulations in the MedeA®environment using the MedeA LAMMPS interface. The emphasis is on the background to EAM based simulation and the properties which may be obtained using the method.

Embedded Atom Method (EAM) Simulation in the MedeA Environment

In the design of substrates for the epitaxial growth of InxGa1-xN alloys it is useful to know the elastic properties of the semiconductor as a function of composition. This application note shows the use of MedeA®with VASP 5.2 and the mechanical-thermal (MT) module in computing these properties. Judging by the results for the binaries GaN and InN, the level of accuracy is comparable with that achieved in experiments.

Structure and Elastic Properties of (In,Ga)N

First-principles calculations reveal that the magnetic moments of atoms on an Fe(001) surface are 30% larger than in the bulk. This enhancement decays within about three layers towards the bulk, which demonstrates the highly localized character of enhanced surface magnetism in transition metals such as iron.

Magnetism of Fe Surface

Ab Initio calculations (VASP) correctly predict the monoclinic phase to be the most stable at low temperatures, followed by the tetragonal, and the high-temperature cubic phases. The structural information and heat of formation obtained provide a sound starting point for calculations of thermodynamic properties.

Despite the enormous progress in experimental surface science, notably with spectroscopic methods exploiting synchrotron radiation and scanning tunneling microscopy, computations remain one of the most useful sources for accurate data on surface structures. In fact, quite often it is the combination of experimental and computed results, which gives the most reliable data of surface structures. As an example, let us apply MEDEA to the Si(001) surface. This surface is the typical substrate in the manufacturing of semiconducting devices.

Surface Reconstruction of Si(001)

The temperature-induced phase transition from monoclinic to tetragonal ZrO₂ is predicted from first principles calculations using a quasi-harmonic approach for the vibrational enthalpy and entropy. The computed transition temperature is within 15% of the experimental value. Relative trends due to vacancies, alloying elements, and mechanical stress can be expected to have a higher accuracy. The present results show the importance of thermal expansion, which is here also obtained from first principles.

Temperature-Dependent Phase Transitions of ZrO₂

Interfaces are present in most materials and have a large impact on mechanical properties such as stiffness and yield strength. Given that the properties of an interface can radically change by the presence of even minute amounts of impurities, it is of great interest to predict the effect of segregated atoms at interfaces.

As systematic experimental information on the impact of specific defect types on the grain boundary strength is hard to obtain, computational modelling is of great help.

Strength of Ni Grain Boundary and the Effect of Boron

Accurate measurements of diffusion coefficients of atoms in solids are difficult and deviations between different experiments can be several orders of magnitude. For the benchmark case of hydrogen diffusion in nickel first-principles calculations give a remarkable agreement with available experimental data especially near room temperature. Thus, computations of diffusion coefficients can be comparable in reliability with measured data. Simulations are possible for situations such as high strain, or slow processes where measurements are difficult or impractical.

Diffusion of Hydrogen in Nickel

The structure of the Fe₂O₃ (0001) surface as a function of oxygen partial pressure and temperature is computed using first-principles thermodynamics. The results reveal an active four-fold coordinated surface Fe atom which releases oxygen atoms at approximately 850 K at ambient oxygen partial pressure. This property is likely to be related to the catalytic activity of hematite for selective oxidation reactions such as the oxidation from ethylbenzene to styrene.

Structure of an iron oxide (Fe₂O₃) surface as function of temperature and O₂ pressure

First-principles calculations reveal the atomistic structures of the active phases of CoMoS and NiMoS hydrodesulfurization catalysts. The reliable determination of the catalyst surface is critical, as it represents the starting point for subsequent adsorption and reaction path simulations. The predicted dominant structures are consistent with experimental STM, EXAFS and magnetic susceptibility measurements.

Atomic Structure of Hydrodesulfurization (HDS) Catalysts

First-principles calculations reveal a three-fold increase in the Young’s modulus of graphite as it is lithiated (C→LiC₆). A linear expression is determined that describes the approximate stiffness of Li intercalated graphite as a function of loading which may lead to greatly improved continuum models of electrode deformation and failure.

Graphite Electrode Elastic Properties upon Li Intercalation

First-principles calculations reveal the dominant acid sites on the amorphous silica-alumina (ASA) catalyst surface. Based on the strength of interaction with Lewis base probe molecules, the bridging and pseudo-bridging silanol (PBS) species are determined to be the most acidic and therefore catalytically active groups. Evident in the calculations is the interplay between the Lewis acid and the Brønsted acid sites on the ASA catalyst surface, giving rise to the enhanced acidity of the PBS species. This work provides mechanistic insight to inform efforts at rationally engineering enhanced ASA-based solid acids.

Acidity of Amorphous Silica-Alumina (ASA) Catalysts

First-principles calculations using VASP reveal the lowest energy reaction pathway for the catalyzed skeletal isomerization of 2-pentene by the acidic zeolite H-ZSM. Three potential mechanisms were evaluated: an ethyl-shift pathway, a dimethylcyclopropane (DMCP) intermediate pathway and a pathway involving an edge-protonated DMCP species. The results indicate that the DMCP intermediate pathway is the kinetically preferred pathway with a classical barrier height of 98 kJ/mol. Evident in the calculations is the influence of the transient intermediate stability along the reaction path; with secondary carbenium ions leading to energetically favored mechanisms.

Catalytic Isomerization of 2-pentene in H-ZSM-22

The compressibility, tensile strength, and mechanical resistance to shear of a solid are related to the corresponding moduli (bulk, Young’s, and shear modulus), which are derived from the coefficients of elasticity. First-principles calculations of these fundamental mechanical properties give values of the same quality as experimental data, but at a substantially smaller effort and cost. This is demonstrated here for cubic silicon carbide, β-SiC, corundum, α-Al₂O₃, and a tourmaline with a fairly complex crystal structure. First-principles calculations are a valuable source for these fundamental materials property data.

Elastic coefficients and moduli for cubic silicon carbide (β-SiC), corundum (α-Al₂O₃), and a tourmaline crystal (Schorl)

The crystal structure of a purely organic, hydrogen-bonded molecular crystal is very well described by density functional theory with a gradient-corrected Perdew-Becke-Ernzerhof potential. The computations were preformed with the VASP program using the projector augmented wave method with a plane wave basis set. The agreement between computed and experimental lattice parameters is better than 2% with a tendency of the calculations to overestimate the bond lengths. The calculations provide equilibrium positions for the hydrogen atoms, which are difficult to place based on x-ray diffraction data.

Crystal Structure of Glucose: Placing Hydrogen Atoms by Computations

In many metal-ceramic composites the interface between the metallic and ceramic phases determines the mechanical properties of the material. A prominent example is the WC-Co composite, where the combination of high WC hardness and Co ductility results in advantageous mechanical properties for applications in the tool manufacturing industry. The outstanding performance of WC-Co can be explained by the low interface energies (high stability) of the contacting surfaces of WC and Co.

Interface Energy of Metal-Ceramic Interface Co/WC Using ab initio  Thermodynamics

The physics and chemistry of materials containing defects is of great interest in areas such as semiconductors, metal alloys and compounds, magnetic systems and optical materials.
Destabilizing or stabilizing crystalline bulk systems, surfaces and interfaces through additives is a common technique, for example recent research into complex hydrides aims to tune the kinetics of the adsorption/desorption changing the stoichiometry and composition of the base materials.

Building and Analyzing Indium defects in GaAs

The energy of adsorption and dissociation of molecules on surfaces plays a critical role in technological processes such as chemical vapor deposition, catalysis, and corrosion. The present case shows the calculation of the energy of the dissociative chemisorption of a silane molecule on a Si (001) surface.

Energy of Dissociative Chemisorption of SiH₄ on Si (001) Surface

This application shows the interaction of carbon monoxide with rutile.
An answer is given to the question whether CO binds with the carbon or the oxygen molecule to the surface.

CO Adsorption on a TiO₂ Surface

The molecular builder (Molecular Builder) is part of the MEDEA standard suite of building tools. This tutorial provides an overview of the Molecular Builder’s basic functionality.

In this application note, we show how to build the isomers ethyl alcohol and di-ethyl ether and calculate their respective heats of formation.

Heat of Formation of Ethyl Alcohol and Dimethyl Ether

A key process in the semiconductor manufacturing is the reactive adsorption of molecules such as silane (SiH4) and dichlorosilane (SiCl2H2) on the surfaces of silicon wafers. This case study demonstrates the calculation of the geometry of a silane molecule on a reconstructed Si(001) surface.

Dissociation of SiH₄ on Si(001) Surface

First-principles computations correctly describe the ferroelectric distortions and the macroscopic polarization of BaTiO₃ in agreement with experiment. Computations of the vibrational properties (phonons) reveal that a cubic perovskite structure of BaTiO₃ becomes stable under compression of the lattice. This demonstrates the usefulness of first-principles calculations in the design and optimization of ferroelectric materials.

Ferroelectric Properties of BaTiO₃

In addition to the lattice parameters, the rutile structure has an additional degree of freedom, namely the position of the oxygen atom. This case study shows the simultaneous calculation of all these structural degrees of freedom.

Structure and Bond Lengths in Titanium Dioxide (TiO₂)

This application shows the calculation of the elastic constants of TiB₂, a hexagonal structure.

Elastic Constants of TiB₂

The surface energy of a material is defined as the energy required to create a surface (h k l) from the bulk material. Surface energies are usually given in units of J/m2.

Surface Energy of Molybdenum

The work function of the metal gate in a CMOS stack depends on the composition and structure of the interfaces. This is demonstrated here for the case of a Si-HfO2-W stack by introducing a Hf vacancy at the Si/HfO2 interface. At a concentration of 1.2 vacancies per nm2 the work functions is increased by 500 meV.

Modeling work function changes in CMOS stacks containing HfO₂ high-k dielectrics

Hydrides containing alkaline-earth metals are prototypes for hydrogen storage materials. For this purpose the heat of formation and the mechanical properties are of fundamental interest. First- principles electronic structure methods provide systematic values for these materials properties. The agreement with experimental data for the heat of formation is good. Presently, no experimental data for the elastic coefficients of these metal hydrides are available thus leaving the computed data as the sole source.

Alkaline-earth hydrides

Increasingly stringent environmental regulations require a lowering of sulfur in Diesel fuels. This is accomplished by a catalytic process transforming sulfur-containing molecules into H2S, which is removed from the liquid phase. Larger sulfur-organic molecules are more difficult to attack and new catalytic materials are needed. The present screening study demonstrates how the combination of experimental activity data, crystallographic information from structural databases and first- principles computation of binding energies are used to identify potential new candidates.

Screening of desulfurization catalysts

Elemental germanium is a semiconductor with a measured indirect band gap of 0.66 eV. Using a hybrid functional as implemented in VASP 5.2, the computed value is 0.66 eV while standard density functional approaches incorrectly predict Ge to have no band gap. Other features of the band structure such as the direct gap at Γ are also well reproduced by the current level of theory, namely 0.8 eV (measured) and 0.73 eV (computed), thus demonstrating the reliability of this level of approach in predicting energy band structures. This sets the stage for using computations to modify the band structure for example by uniaxial strain to meet specific design criteria.

Energy band structure of germanium

Iodine is a fission product of uranium. It can attack the inner side of zircaloy cladding in nuclear power reactors leading to cracking and fracture. Computations show that iodine molecules adsorb and dissociate on a zirconium surface without an energy barrier. The binding energy of iodine on this surface is large (nearly 300 kJ/mol per iodine atom), but the barriers for surface diffusion is only 6.8 kJ/mol. This gives rise to rapid surface diffusion allowing iodine atoms to reach the crack tips faster than the propagation of cracks.

Adsorption and Dissociation of Iodine Molecules on a Zr Surface

The correct lattice parameters of a crystalline structure are of fundamental importance for any reliable computational predictions of materials properties. This case study shows the calculation of the lattice parameters of titanium carbide, TiC.

Structure of bulk Titanium Carbide (TiC)

The performance of catalytic materials depends on complex phenomena linked to chemical
composition, preparation, activation procedures, and surface conditions under operational
conditions. This complexity requires a comprehensive arsenal of R&D approaches including
theoretical and computational methods. While many fundamental research efforts are
currently directed at a detailed understanding of surface reaction mechanisms, PREDIBOND™
focuses on bond strength and local chemical environment as central descriptors of chemical

PREDIBOND™ in Heterogeneous Catalysis: Predicting Activity Patterns

This application note deals with positioning molecules on surfaces. As an example we will investigate the adsorption of ethyl alcohol (ethanol) on a Cu (111) surface. In doing so we will consider two possible configurations for the adsorbed molecule:
1. adsorbed ethanol
2. dehydrogenated ethanol, i.e. an ethoxygen.
Running structure relaxations using VASP produces a first estimate of the relative stability of these two systems.

Dehydrogenation Energy of Ethyl Alcohol (Ethanol) on a Cu (111) Surface

Calculation of elastic properties is straightforward with MedeA's Mechanical and Thermal Properties module.

Elastic Properties of Diamond

This case study covers the practical use of MEDEA to calculate thermochemical functions for solids, molecules and atoms. We will use VASP and PHONON for this, but the current document focuses on the thermochemistry and not the details of the calculations.

Practical Thermochemistry: Sodium Metal, Chlorine Gas, and Solid Sodium  Chloride

Despite the enormous progress in experimental surface science, notably with spectroscopic
methods exploiting synchrotron radiation and scanning tunneling microscopy, computations
remain one of the most useful sources for accurate data on surface structures. In fact, quite often
it is the combination of experimental and computed results, which gives the most reliable data of
surface structures. As an example, let us apply MEDEA to the Si(001) surface. This surface is the
typical substrate in the manufacturing of semiconducting devices.

Surface Reconstruction of Si(001)

The insertion of interstitial impurities in a host lattice causes local deformations of the lattice. The
purpose of this case study is the comparison of such deformations caused by boron and fluorine
impurities in a silver lattice.

Deformation of Silver Lattice by Interstitial Boron and Fluorine Impurities

Ferromagnetism has its quantum mechanic origin in the difference of spin-up and spin-down
electron densities. It is driven by a balance between a gain in exchange energy due to larger spin-
polarization and a loss in Coulomb repulsion and kinetic energy of the electron system.

The purpose of this study is the computation of the cleavage energy of a material, i.e. the energy
required to split a material into two parts. This could be a bulk material, a grain boundary, or an
interface. To this end, one needs to compute the total energy of the bulk solid and the material
with a free surface.

Cleavage Energy of TiN

In Chromium and Chromium alloys antiferromagnetic ordering and spin-density-waves (SDW) states are at the origin of many physical properties like thermal expansion, elastic constants, and electrical resistivity among others.
This document summarizes structural and elastic properties of Chromium, computed from
first principles.

Chromium: Structure and Elastic Properties

The convergence of the total electronic energy as computed with VASP is determined by two key
computational parameters, namely the number of basis functions (plane wave cutoff) and the
number of k-points (k-spacing). In addition the integration of the states near the Fermi level is
influenced by a smearing parameter.

Convergence of Total Energy with Plane Wave Cutoff and k-Mesh:  Mo, Al, and LiF

The cohesive energy of a solid is defined as the energy required for separating the condensed material into isolated free atoms. Cohesive energies range from about 0.1 eV or 10 kJ/mol for inert gases up to about 8 eV or 800 kJ/mol per atom for strongly bound materials such as diamond or tungsten.

Cohesive Energy of Diamond