ICME GRP1 HW1
The concept of embedding energy, that is, the amount of energy needed to add or remove one atom from a lattice, was originally proposed by Friedel (1952) and developed by Stott and Zaremba (1980). Daw and Baskes (1984) produced the Embedded Atom Method (hereafter EAM). This was a numerical method for calculating the embedding energy in metals. Baskes (1993) followed this with the Modified Embedded Atom Method (hereafter MEAM), which extended the application of EAM to bodies of more than one phase. This was validated by using it to address the two phase bonding of Aluminum and Silicon.
The potentials developed by MEAM are not unique in that many different potentials have been generated using this same method, each tailored to address the simulation of specific physical properties of a given system. To address this, Tschopp, Solanki, Baskes, Gao, Sun, and Horstemeyer (2011), provided a heuristic, similar to an Ashby property map, for the formulation of new potentials, and assessing the problem domain of their applicability.
MEAM, and therefore also EAM, establishes a functional relationship between the cohesive energy of an atom and the electron density of volume in which said atom is placed or removed. The function ρ(r) is used to represent the contribution of neighboring atoms in the lattice to the electron density. The total energy E, is computing by the following relation:
Where F is the embedding energy of the ith material, r is the distance between atoms of the ith and jth materials, and ∅ represents the pair potential of the ith and jth materials. Please note that the superscript in the above formula is not indicative of indicial summation notation unless otherwise noted. This energy is frequently passed along to Molecular Dynamics (MD) code to calculate the forces at work on each atom.
In calculating these atomic forces, the dipole force tensor β_ij represents the equivalent to a stress tensor at a given point at the atomistic length scale, if, hypothetically, stress could be defined and, therefore, applied at that scale. Its form is given by:
where ω^i represents the atomic volume of phase I, fki is the force vector between atoms, and rmij is the displacement vector between two atoms. This is used, in turn, to calculate the actual continuum stress tensor, σ, given by:
where N* is the total number of atoms in the volume of material of interest.
Gall, Horstemeyer, Schilfgaarde, and Baskes (1999) further defined uniaxial strain as:
where v1 is the atomic velocity, t is time, and L is length.
Potential Calibration, Optimization, & Validation
Of critical note is the fact that the inter-atomic potential, ∅, in the formulation of total energy E, is not unique. It is, in fact, a function of a number of parameters, whose inclusion in ∅’s functional form is based upon the intended use of the pair potential. This limits the use of a given ∅ to the domain of problems for which it was originally formulated to help describe. Tschopp, et. al (2011) proposed a “Design Map” methodology, similar to an Ashby property map, in order to address the assessment of a given potentials domain of applicability, as well as develop new potentials tailored to simulated specific physical system characteristics. This “design map” methodology is decomposed into two phases, a global approach and a local approach.
The global approach may be thought of as a coarse refinement of potential parameters, whereby said parameters are selected from a set of all possible parameters affecting the inter atomic potential of a system. It may be further decomposed into three sub-phases; those being Initialization, Potential Space Evaluation, and Sensitivity Analysis.
During Initialization, all potential parameters, appropriate formation energies and their corresponding ab initio and/or experimental values are identified. These parameters and energies are passed on to the Potential Space Evaluation phase, where the ranges of potential parameters are defined in an effort to calibrate the potential. Each parameter is defined such that the potential is brought into line with the formation energies. Additionally simulated ab initio values are compared with experimental values during this phase. During the Sensitivity Analysis, the number of potential parameters is reduced as each is evaluated in order to determine its influence on the formation energies, and only those parameters with significant effect are retained. In this phase, the potential parameters are normalized by the highest value among them, and each is varied in turn, to determine which of the parameters have the greatest effect on the total energy, E.
The local approach may be decomposed, much like the global approach, save that it is broken down into five sub phases rather than three. These are Potential Space evaluation, Potential Space Sampling, Analytical Model Generation, Potential Design Map Development, and Potential Design Map Validation.
During the Potential Space Evaluation, the range of the parameters selected during the global approach is again calibrated, this time, more precisely. Each of the parameters is normalized individually while the remaining parameters are held fixed. This necessarily bounds the parameter range.
The Potential Space Sampling phase is the most computationally demanding, as a sampling of the range of each potential parameter is taken in order to calculate a sampling of formation energies. Monte Carlo, and Latin Hypercube Sampling have been used to generate a sampling distribution in the past, the latter by Tschopp et al. Parallelization of this process is possible, as the sample space for each parameter may be generated independently.
The Analytical Model Generation phase sees the use of regression analysis techniques to construct a polynomial approximation of the relationship between the potential parameters and the formation energies. This approximation should be checked against nano-scale calculations in order to ascertain its fitness before moving on to the next sub phase. Tschopp et. al, use the R2 correlation coefficient to make such comparisons.
In the Potential Design Map Development phase, the parameters used to generate each potential are organized based upon their relationship with the physical phenomena to be simulated with said potential. This constitutes the Design Map. This map is used to generate an optimized interatomic potential by combining the difference between the analytical formation energies and their ab initio/experimental values and assigning weights to the differences to obtain a global analytical expression, fobj. These weights are designed based on their relevance to the particular properties being simulated with the optimized potential. Tschopp, e. al make use of a constrained, nonlinear optimization technique to generate the optimized potential from the previously generated fobj. This process may be repeated to generate further potentials optimized for different uses. It is important to acknowledge that this phase requires active thought on part of the modeler in terms of deciding which parameters are most critical to the accurate simulation of the physical phenomena under investigation and assign weights to the parameters accordingly.
The Potential Design Map Validation phase is the final phase of the Local Approach. In this phase, the potential is validated by using it to generate physical properties of the system in question that were not used to aid in the calibration of its parameters.
Intrinsic Extrinsic Stacking Fault Energy Convergence
The Density functional theory (DFT) code VASP is used to formulate a description of the ground state properties of BCC, FCC and HCP aluminum. DFT works by describing a system of interacting sub-atomic particles (i.e. electrons) through the use of its density in place of its many-body wave function. Kohn and Sham pointed out that solids can often be considered as close to the limit of the homogeneous electron gas. Thus the Local Density Approximation treats solid materials as having the same density and properties as a homogeneous electron gas.
Results and Discussion
Energy Volume Curves
After running several convergence studies in search of optimal K-point grids, energy cutoffs, and global break conditions, energy volume curves for all three crystal lattice structures were obtained.
Constants and Material Properties
|Cohesive Energy (eV)||4.17177|
|Lattice Parameter (Ang)||3.9787|
|Bulk Modulus (GPa)||82|
|(100) Surface Energy (eV)||890|
|(110) Surface Energy (eV)||960|
|(111) Surface Energy (eV)||780|
|Extrinsic Stacking Fault Energy (eV)||133|
|Intrinsic Stacking Fault Energy||133|
|Vacancy Formation Energy (eV)||0.5|
|Octahedral Interstitial Formation Energy (eV)||2.8|
|Tetrahedral Interstitial Formation Energy (eV)||3.68|
- ↑ M. Horstemeyer, Integrated Computation Materials Engineering (ICME) for Metals, John Wiley & Sons, Inc., 2012.