Embedded Atom Method (EAM) Tabulation¶
An EAM model is defined by constructing instances of atsim.potentials.EAMPotential
describing each species within the model. EAMPotential
encapsulates the density and embedding functions specific to each species’ many bodied interactions. In addition the purely pairwise interactions within the EAM are defined using a list of atsim.potentials.Potential
objects.
Once the EAM model has been described in terms of EAMPotential
and Potential
objects it can be tabulated for specific simulation codes. This is done by using the EAM_Tabulation
objects from the atsim.potentials.eam_tabulation
module:
Class | Format | Simulation Code | Example |
---|---|---|---|
SetFL_EAMTabulation |
setfl , |
LAMMPS | |
SetFL_FS_EAMTabulation |
pair_style eam/fs | LAMMPS | Example 3a: Tabulate Al-Fe Finnis-Sinclair Potentials Using SetFL_FS_EAMTabulation for LAMMPS |
TABEAM_EAMTabulation |
TABEAM |
DL_POLY | Example 2b: Tabulate Al-Cu Alloy Potentials Using TABEAM_EAMTabulation for DL_POLY |
TABEAM_FinnisSinclair_EAMTabulation |
EEAM TABEAM |
DL_POLY | Example 3b: Tabulate Al-Fe Finnis-Sinclair Potentials Using TABEAM_FinnisSinclair_EAMTabulation for DL_POLY |
Excel_EAMTabulation |
.xlsx |
||
Excel_FinnisSinclair_EAMTabulation |
.xlsx |
Even though the use of EAM_Tabulation
objects is preferred a legacy procedural interface is also provided. This is described here: Procedural EAM Tabulation.
Examples¶
Example 1: Using SetFL_EAMTabulation
to Tabulate Ag Potential for LAMMPS¶
This example shows how to use the SetFL_EAMTabulation
class to tabulate an EAM model for the simulation of Ag metal. How to use this tabulation within LAMMPS will then be demonstrated. The final tabulation script can be found in eam_example1.py
.
See also
- A
potable
version of this example is given here: Sutton Ag EAM Example.
Model Description¶
Within this example the Ag potential of Sutton will be tabulated [1]. Within the EAM the energy (\(E_i\)) of an atom \(i\) whose species is \(\alpha\) is given by:
Note
- \(\rho_\beta(r_{ij})\) is the density function which gives the electron density for atom \(j\) with species \(\beta\) as a function of its separation from atom \(i\), \(r_{ij}\).
- The electron density for atom \(i\) is obtained by summing over the density (\(\rho_\beta (r_{ij}\)) contributions due to its neighbours.
- The embedding function \(F_\alpha(\rho)\) is used to calculate the many-bodied energy contribution from this summed electron density.
- The sum \(\frac{1}{2} \sum_{j \neq i} \phi_{\alpha \beta} (r_{ij})\) gives the pair-potential contribution to atom \(i\)’s energy.
- \(\phi_{\alpha \beta} (r_{ij})\) are simply pair potentials that describe the energy between two atoms as a function of their separation.
The embedding function used by Sutton is:
and the density function is:
whilst pair interactions are given by:
The model parameters are given as:
Parameter | Value |
---|---|
\(m\) | 6 |
\(n\) | 12 |
\(a\) | \(2.928323832 \text{Å} \text{eV}^{\frac{1}{3}}\) |
\(b\) | \(2.485883762 \text{eV}^\frac{1}{12} \text{Å}\) |
Define the Model¶
It is now necessary to describe the model in python code. Hard-coding the model parameters from the previous table, embed()
and density()
functions can be defined for \(F_{\text{Ag}} (\rho)\) and \(\rho_{\text{Ag}}\) respectively:
import math
from atsim.potentials import EAMPotential, Potential
from atsim.potentials.eam_tabulation import SetFL_EAMTabulation
def embed(rho):
return -math.sqrt(rho)
def density(rij):
if rij == 0:
return 0.0
return (2.928323832 / rij) ** 6.0
The embedding and density functions should then be wrapped in an EAMPotential
object to create a single item list:
eam_potentials = [EAMPotential("Ag", 47, 107.8682, embed, density)]
Similarly the pair potential component, \(\phi_{\text{Ag}-\text{Ag}} (r_{ij})\), of the model can easily be defined as:
def pair_AgAg(rij):
if rij == 0:
return 0.0
return (2.485883762/rij) ** 12
This can then be wrapped in a atsim.potentials.Potential
object to create a list of pair potentials.
pair_potentials = [Potential('Ag', 'Ag', pair_AgAg)]
Now all the components of the model have been defined a table file can be created in the setfl
format. Before doing this, it is necessary to choose appropriate density and separation cut-offs together with the number of rows in the density/pair functions (nr
) and embedding function (nrho
) respectively:
- Here 50000 density values will be tabulated to a cutoff of 50.0.
- The pair-potential cut-off and the maximum \(r_{ij}\) value for the density function is 12 Å and both will have 12000 rows.
An instance of atsim.potentials.eam_tabulation.SetFL_EAMTabulation
is created with the EAMPotential
and Potential
objects. This object is then used to tabulate the Ag potential by calling the write()
method with the Ag.eam.alloy
file object:
pair_potentials = [Potential('Ag', 'Ag', pair_AgAg)]
cutoff_rho = 50.0
nrho = 50000
cutoff = 12.0
nr = 12000
tabulation = SetFL_EAMTabulation(
pair_potentials,
eam_potentials,
cutoff, nr,
cutoff_rho, nrho)
with open("Ag.eam.alloy", 'w') as outfile:
tabulation.write(outfile)
Putting this together the following script is obtained (this script can also be downloaded eam_example1.py
:
#! /usr/bin/env python
import math
from atsim.potentials import EAMPotential, Potential
from atsim.potentials.eam_tabulation import SetFL_EAMTabulation
def embed(rho):
return -math.sqrt(rho)
def density(rij):
if rij == 0:
return 0.0
return (2.928323832 / rij) ** 6.0
def pair_AgAg(rij):
if rij == 0:
return 0.0
return (2.485883762/rij) ** 12
def main():
# Create EAMPotential
eam_potentials = [EAMPotential("Ag", 47, 107.8682, embed, density)]
pair_potentials = [Potential('Ag', 'Ag', pair_AgAg)]
cutoff_rho = 50.0
nrho = 50000
cutoff = 12.0
nr = 12000
tabulation = SetFL_EAMTabulation(
pair_potentials,
eam_potentials,
cutoff, nr,
cutoff_rho, nrho)
with open("Ag.eam.alloy", 'w') as outfile:
tabulation.write(outfile)
if __name__ == "__main__":
main()
Running this script will produce a table file named Ag.eam.alloy
in the same directory as the script:
python eam_example1.py
Using the Ag.eam.alloy
file within LAMMPS¶
This section of the example will now demonstrate how the table file can be used used to perform a static energy minimisation of an FCC Ag structure in LAMMPS.
Place the following in a file called fcc.lmpstruct
in the same directory as the Ag.eam.alloy
file you created previously. This describes a single FCC cell with a wildly inaccurate lattice parameter:
Title
4 atoms
1 atom types
0.0 5.000000 xlo xhi
0.0 5.000000 ylo yhi
0.0 5.000000 zlo zhi
0.000000 0.000000 0.000000 xy xz yz
Masses
1 107.86820000000000163709 #Ag
Atoms
1 0 1 0.000000 0.000000 0.000000 0.000000
2 0 1 0.000000 2.500000 2.500000 0.000000
3 0 1 0.000000 0.000000 2.500000 2.500000
4 0 1 0.000000 2.500000 0.000000 2.500000
The following LAMMPS input file describes a minimisation run. The lines describing potentials are highlighted. Put its contents in a file called example_eam_alloy_minimize.lmpin
:
units metal
boundary p p p
atom_style full
read_data fcc.lmpstruct
pair_style eam/alloy
pair_coeff * * Ag.eam.alloy Ag
fix 1 all box/relax x 0.0 y 0.0 z 0.0
minimize 0.0 1.0e-8 1000 100000
The pair_style eam/alloy
command tells LAMMPS to use the EAM and expect pair_coeff
commands mapping atom types to particular table files:
pair_style eam/alloy
The following pair_coeff
directive indicates that the interaction between atom-type 1 (Ag) with itself should use the setfl
formatted file contained within Ag.eam.alloy
. The Ag
label at the end of the line indicates that atom-type 1 should be associated with this label in the table file:
pair_coeff * * Ag.eam.alloy Ag
The example can then be run by invoking LAMMPS:
lammps -in example_eam_alloy_minimize.lmpin
Example 2a: Tabulate Al-Cu Alloy Potentials Using SetFL_EAMTabulation
for LAMMPS¶
Within the following example the process required to generate and use a setfl
file that tabulates the Al-Cu alloy model of Zhou et al [2]. In comparison to the previous example this example contains density and embedding functions for multiple elements and includes pair-potentials specific to pairs of interacting species. The eam_example2a.py
file gives a complete example of how the Zhou model can be tabulated.
Model Description¶
The model makes use of the EAM as described above (see Example 1 Model Description) . The density function, \(\rho_\beta (r_{ij})\) for an atom \(j\) of species \(\beta\) separated from atom \(i\) by \(r_{ij}\) is:
where \(f_e\), \(r_e\), \(\omega\) and \(\lambda\) are parameters specific to species \(\beta\). The pair-potential function acting between species \(\alpha\)–\(\beta\) is obtained by combining the density functions of the interacting species:
The homogeneous pair-interactions, \(\phi_{\alpha\alpha}(r_{ij})\) and \(\phi_{\beta\beta}(r_{ij})\) have the form:
again, \(A\), \(B\), \(\gamma\), \(\omega\), \(\kappa\) and \(\omega\) are parameters specific to the species \(\alpha\).
The embedding function for each species, \(F_\alpha (\rho)\), is defined over three density ranges using the following:
The model parameters for Cu and Al are given in the following table:
Parameter Cu Al \(r_e\) 2.556162 2.863924 \(f_e\) 1.554485 1.403115 \(\rho_e\) 21.175871 20.418205 \(\rho_s\) 21.175395 23.195740 \(\gamma\) 8.127620 6.613165 \(\omega\) 4.334731 3.527021 \(A\) 0.396620 0.314873 \(B\) 0.548085 0.365551 \(\kappa\) 0.308782 0.379846 \(\lambda\) 0.756515 0.759692 \(F_{n0}\) -2.170269 -2.807602 \(F_{n1}\) -0.263788 -0.301435 \(F_{n2}\) 1.088878 1.258562 \(F_{n3}\) -0.817603 -1.247604 \(F_0\) -2.19 -2.83 \(F_1\) 0 0 \(F_2\) 0.561830 0.622245 \(F_3\) -2.100595 -2.488244 \(\eta\) 0.310490 0.785902 \(F_e\) -2.186568 -2.824528
Note
The Al \(A\) value is given as 0.134873 in Zhou’s original Phys. Rev. B paper. However parameter file provided by Zhou for this model, at http://www.ctcms.nist.gov/potentials/Zhou04.html gives the parameter as 0.314873. It is this latter value that is used here.
In addition the final term of the embedding function has been modified to match that used in fortran tabulation code also provided at http://www.ctcms.nist.gov/potentials/Zhou04.html
Define the Model¶
A series of python functions are defined to describe the embedding, density and pair interaction functions. To encourage code re-use a number of function factories are defined. Using the parameters passed to them they return specialised functions appropriate for the parameters. The given factory functions make use of python’s support for closures in their implementation.
The makeFunc()
factory function is used to define density functions. As this functional form is also used as a component of the pair-potentials makeFunc()
is re-used within the makePairPotAA()
factory function.
def makeFunc(a, b, r_e, c):
# Creates functions of the form used for density function.
# Functional form also forms components of pair potential.
def func(r):
return (a * math.exp(-b*(r/r_e - 1)))/(1+(r/r_e - c)**20.0)
return func
The following factory returns the functions used to describe the homogeneous Al-Al and Cu-Cu pair-potential interactions:
def makePairPotAA(A, gamma, r_e, kappa,
B, omega, lamda):
# Function factory that returns functions parameterised for homogeneous pair interactions
f1 = makeFunc(A, gamma, r_e, kappa)
f2 = makeFunc(B, omega, r_e, lamda)
def func(r):
return f1(r) - f2(r)
return func
Whilst makePairPotAB()
describes the Al-Cu pair-potential:
def makePairPotAB(dens_a, phi_aa, dens_b, phi_bb):
# Function factory that returns functions parameterised for heterogeneous pair interactions
def func(r):
return 0.5 * ((dens_b(r)/dens_a(r) * phi_aa(r)) + (dens_a(r)/dens_b(r) * phi_bb(r)))
return func
The makeEmbed()
function describes the embedding function:
def makeEmbed(rho_e, rho_s, F_ni, F_i, F_e, eta):
# Function factory returning parameterised embedding function.
rho_n = 0.85*rho_e
rho_0 = 1.15*rho_e
def e1(rho):
return sum([F_ni[i] * (rho/rho_n - 1)**float(i) for i in range(4)])
def e2(rho):
return sum([F_i[i] * (rho/rho_e - 1)**float(i) for i in range(4)])
def e3(rho):
return F_e * (1.0 - eta*math.log(rho/rho_s)) * (rho/rho_s)**eta
def func(rho):
if rho < rho_n:
return e1(rho)
elif rho_n <= rho < rho_0:
return e2(rho)
return e3(rho)
return func
Lists of EAMPotential
and Potential
objects are created and returned as a tuple by the makePotentialObjects()
function within eam_example2a.py
. Before invoking the factory functions we just defined, the model parameters are assigned to easily identifiable variables within this function:
def makePotentialObjects():
# Potential parameters
r_eCu = 2.556162
f_eCu = 1.554485
gamma_Cu = 8.127620
omega_Cu = 4.334731
A_Cu = 0.396620
B_Cu = 0.548085
kappa_Cu = 0.308782
lambda_Cu = 0.756515
rho_e_Cu = 21.175871
rho_s_Cu = 21.175395
F_ni_Cu = [-2.170269, -0.263788, 1.088878, -0.817603]
F_i_Cu = [-2.19, 0.0, 0.561830, -2.100595]
eta_Cu = 0.310490
F_e_Cu = -2.186568
r_eAl = 2.863924
f_eAl = 1.403115
gamma_Al = 6.613165
omega_Al = 3.527021
# A_Al = 0.134873
A_Al = 0.314873
B_Al = 0.365551
kappa_Al = 0.379846
lambda_Al = 0.759692
rho_e_Al = 20.418205
rho_s_Al = 23.195740
F_ni_Al = [-2.807602, -0.301435, 1.258562, -1.247604]
F_i_Al = [-2.83, 0.0, 0.622245, -2.488244]
eta_Al = 0.785902
F_e_Al = -2.824528
# Define the density functions
dens_Cu = makeFunc(f_eCu, omega_Cu, r_eCu, lambda_Cu)
dens_Al = makeFunc(f_eAl, omega_Al, r_eAl, lambda_Al)
# Finally, define embedding functions for each species
embed_Cu = makeEmbed(rho_e_Cu, rho_s_Cu, F_ni_Cu, F_i_Cu, F_e_Cu, eta_Cu)
embed_Al = makeEmbed(rho_e_Al, rho_s_Al, F_ni_Al, F_i_Al, F_e_Al, eta_Al)
# Wrap them in EAMPotential objects
eam_potentials = [
EAMPotential("Al", 13, 26.98, embed_Al, dens_Al),
EAMPotential("Cu", 29, 63.55, embed_Cu, dens_Cu)]
# Define pair functions
pair_CuCu = makePairPotAA(A_Cu, gamma_Cu, r_eCu, kappa_Cu,
B_Cu, omega_Cu, lambda_Cu)
pair_AlAl = makePairPotAA(A_Al, gamma_Al, r_eAl, kappa_Al,
B_Al, omega_Al, lambda_Al)
pair_AlCu = makePairPotAB(dens_Cu, pair_CuCu, dens_Al, pair_AlAl)
# Wrap them in Potential objects
pair_potentials = [
Potential('Al', 'Al', pair_AlAl),
Potential('Cu', 'Cu', pair_CuCu),
Potential('Al', 'Cu', pair_AlCu)]
return eam_potentials, pair_potentials
Now the functions required by the EAMPotential
instances for Al and Cu can be created:
# Define the density functions
dens_Cu = makeFunc(f_eCu, omega_Cu, r_eCu, lambda_Cu)
dens_Al = makeFunc(f_eAl, omega_Al, r_eAl, lambda_Al)
# Finally, define embedding functions for each species
embed_Cu = makeEmbed(rho_e_Cu, rho_s_Cu, F_ni_Cu, F_i_Cu, F_e_Cu, eta_Cu)
embed_Al = makeEmbed(rho_e_Al, rho_s_Al, F_ni_Al, F_i_Al, F_e_Al, eta_Al)
Now these are wrapped up in EAMPotential
objects to give the eamPotentials
list:
eam_potentials = [
EAMPotential("Al", 13, 26.98, embed_Al, dens_Al),
EAMPotential("Cu", 29, 63.55, embed_Cu, dens_Cu)]
Similarly, using the makePairPotAA()
and makePairPotAB()
function factories the Potential
objects required for the tabulation are defined:
# Define pair functions
pair_CuCu = makePairPotAA(A_Cu, gamma_Cu, r_eCu, kappa_Cu,
B_Cu, omega_Cu, lambda_Cu)
pair_AlAl = makePairPotAA(A_Al, gamma_Al, r_eAl, kappa_Al,
B_Al, omega_Al, lambda_Al)
pair_AlCu = makePairPotAB(dens_Cu, pair_CuCu, dens_Al, pair_AlAl)
# Wrap them in Potential objects
pair_potentials = [
Potential('Al', 'Al', pair_AlAl),
Potential('Cu', 'Cu', pair_CuCu),
Potential('Al', 'Cu', pair_AlCu)]
Now we have all the objects required for SetFL_EAMTabulation
. The next excerpt calls makeObjects()
to get the EAM and pair-potential objects before creating the tabulation object, and invoking its write()
method to write the data into a file called Zhou_AlCu.eam.alloy
:
def main():
eam_potentials, pair_potentials = makePotentialObjects()
# Perform tabulation
# Make tabulation
cutoff_rho = 100.0
nrho = 2000
cutoff = 6.0
nr = 2000
tabulation = SetFL_EAMTabulation(
pair_potentials,
eam_potentials,
cutoff, nr,
cutoff_rho, nrho
)
with open("Zhou_AlCu.eam.alloy", 'w') as outfile:
tabulation.write(outfile)
Putting this all together gives the following script (which can also be downloaded using the following link eam_example2a.py
:). Running this (python eam_example2a.py
) produces the Zhou_AlCu.eam.alloy
file in current working directory.
#! /usr/bin/env python
import math
from atsim.potentials import EAMPotential, Potential
from atsim.potentials.eam_tabulation import SetFL_EAMTabulation
def makeFunc(a, b, r_e, c):
# Creates functions of the form used for density function.
# Functional form also forms components of pair potential.
def func(r):
return (a * math.exp(-b*(r/r_e - 1)))/(1+(r/r_e - c)**20.0)
return func
def makePairPotAA(A, gamma, r_e, kappa,
B, omega, lamda):
# Function factory that returns functions parameterised for homogeneous pair interactions
f1 = makeFunc(A, gamma, r_e, kappa)
f2 = makeFunc(B, omega, r_e, lamda)
def func(r):
return f1(r) - f2(r)
return func
def makePairPotAB(dens_a, phi_aa, dens_b, phi_bb):
# Function factory that returns functions parameterised for heterogeneous pair interactions
def func(r):
return 0.5 * ((dens_b(r)/dens_a(r) * phi_aa(r)) + (dens_a(r)/dens_b(r) * phi_bb(r)))
return func
def makeEmbed(rho_e, rho_s, F_ni, F_i, F_e, eta):
# Function factory returning parameterised embedding function.
rho_n = 0.85*rho_e
rho_0 = 1.15*rho_e
def e1(rho):
return sum([F_ni[i] * (rho/rho_n - 1)**float(i) for i in range(4)])
def e2(rho):
return sum([F_i[i] * (rho/rho_e - 1)**float(i) for i in range(4)])
def e3(rho):
return F_e * (1.0 - eta*math.log(rho/rho_s)) * (rho/rho_s)**eta
def func(rho):
if rho < rho_n:
return e1(rho)
elif rho_n <= rho < rho_0:
return e2(rho)
return e3(rho)
return func
def makePotentialObjects():
# Potential parameters
r_eCu = 2.556162
f_eCu = 1.554485
gamma_Cu = 8.127620
omega_Cu = 4.334731
A_Cu = 0.396620
B_Cu = 0.548085
kappa_Cu = 0.308782
lambda_Cu = 0.756515
rho_e_Cu = 21.175871
rho_s_Cu = 21.175395
F_ni_Cu = [-2.170269, -0.263788, 1.088878, -0.817603]
F_i_Cu = [-2.19, 0.0, 0.561830, -2.100595]
eta_Cu = 0.310490
F_e_Cu = -2.186568
r_eAl = 2.863924
f_eAl = 1.403115
gamma_Al = 6.613165
omega_Al = 3.527021
# A_Al = 0.134873
A_Al = 0.314873
B_Al = 0.365551
kappa_Al = 0.379846
lambda_Al = 0.759692
rho_e_Al = 20.418205
rho_s_Al = 23.195740
F_ni_Al = [-2.807602, -0.301435, 1.258562, -1.247604]
F_i_Al = [-2.83, 0.0, 0.622245, -2.488244]
eta_Al = 0.785902
F_e_Al = -2.824528
# Define the density functions
dens_Cu = makeFunc(f_eCu, omega_Cu, r_eCu, lambda_Cu)
dens_Al = makeFunc(f_eAl, omega_Al, r_eAl, lambda_Al)
# Finally, define embedding functions for each species
embed_Cu = makeEmbed(rho_e_Cu, rho_s_Cu, F_ni_Cu, F_i_Cu, F_e_Cu, eta_Cu)
embed_Al = makeEmbed(rho_e_Al, rho_s_Al, F_ni_Al, F_i_Al, F_e_Al, eta_Al)
# Wrap them in EAMPotential objects
eam_potentials = [
EAMPotential("Al", 13, 26.98, embed_Al, dens_Al),
EAMPotential("Cu", 29, 63.55, embed_Cu, dens_Cu)]
# Define pair functions
pair_CuCu = makePairPotAA(A_Cu, gamma_Cu, r_eCu, kappa_Cu,
B_Cu, omega_Cu, lambda_Cu)
pair_AlAl = makePairPotAA(A_Al, gamma_Al, r_eAl, kappa_Al,
B_Al, omega_Al, lambda_Al)
pair_AlCu = makePairPotAB(dens_Cu, pair_CuCu, dens_Al, pair_AlAl)
# Wrap them in Potential objects
pair_potentials = [
Potential('Al', 'Al', pair_AlAl),
Potential('Cu', 'Cu', pair_CuCu),
Potential('Al', 'Cu', pair_AlCu)]
return eam_potentials, pair_potentials
def main():
eam_potentials, pair_potentials = makePotentialObjects()
# Perform tabulation
# Make tabulation
cutoff_rho = 100.0
nrho = 2000
cutoff = 6.0
nr = 2000
tabulation = SetFL_EAMTabulation(
pair_potentials,
eam_potentials,
cutoff, nr,
cutoff_rho, nrho
)
with open("Zhou_AlCu.eam.alloy", 'w') as outfile:
tabulation.write(outfile)
if __name__ == '__main__':
main()
Using the Zhou_AlCu.eam.alloy
file within LAMMPS¶
Within LAMMPS the setfl
files generated by SetFL_EAMTabulation
are used with the eam/alloy pair_style. The pair_coeff
directive used with this pair_style
effectively maps LAMMPS species numbers to the element names within the table file.
Single Element Systems¶
Assuming a LAMMPS system containing only Al (i.e. Al is species 1) then the pair_style
and pair_coeff
directives would be given as:
pair_style eam/alloy
pair_coeff * * Zhou_AlCu.eam.alloy Al
Likewise if a copper system was being simulated:
pair_style eam/alloy
pair_coeff * * Zhou_AlCu.eam.alloy Cu
Mixed Al-Cu System¶
For an Al-Cu system where Al is species 1 and Cu species 2 then the directives would be:
pair_style eam/alloy
pair_coeff * * Zhou_AlCu.eam.alloy Al Cu
Or if Cu was 1 and Al 2:
pair_style eam/alloy
pair_coeff * * Zhou_AlCu.eam.alloy Cu Al
Example 2b: Tabulate Al-Cu Alloy Potentials Using TABEAM_EAMTabulation
for DL_POLY¶
The tabulation script used with Example 2a can be easily modified to produce the TABEAM
format expected by the DL_POLY simulation code. See the tabulation script for this example: eam_example2b.py
.
The EAMPotential
and Potential
lists are created in exactly the same way as Example 2a, however rather than creating an instance of SetFL_EAMTabulation
in the main()
function it is modified to use the DL_POLY specific TABEAM_EAMTabulation
class instead and to write into a file named TABEAM
. The main()
function of eam_example2b.py
is now given:
def main():
eamPotentials, pairPotentials = makePotentialObjects()
# Perform tabulation
# Make tabulation
cutoff_rho = 100.0
nrho = 2000
cutoff = 6
nr = 2000
tabulation = TABEAM_EAMTabulation(
pairPotentials, eamPotentials, cutoff, nr, cutoff_rho, nrho)
with open("TABEAM", 'w') as outfile:
tabulation.write(outfile)
Using the TABEAM
file with DL_POLY¶
Running eam_example2b.py
will create a file named TABEAM
in the working directory. This should be copied into the simulation directory containing the DL_POLY input files (CONTROL
, CONFIG
and FIELD
).
The following should be added at the bottom of the FIELD
file:
metal 3
Al Al eam
Cu Cu eam
Al Cu eam
Example 3a: Tabulate Al-Fe Finnis-Sinclair Potentials Using SetFL_FS_EAMTabulation
for LAMMPS¶
This example will show how to reproduce the EAM model described by Mendelev et al. for Fe segregation at grain boundaries within Al [3]. As a result this example effectively shows how to reproduce the AlFe_mm.eam.fs
file provided with the LAMMPS source distribution using the SetFL_FS_EAMTabulation
class.
See also
The file format created by atsim.potentials.eam_tabulation.SetFL_FS_EAMTabulation
is supported by the LAMMPS pair_style eam/fs
command. This adds an additional level of flexibility in comparison to the eam/alloy
style; when calculating the density surrounding an atom with species \(\alpha\), each neighbouring atom’s contribution to the density is calculated as a function of its separation from the central atom using \(\rho_{\alpha\beta}(r_{ij})\). This means that the density function is now specific to both the central atom species, \(\alpha\) and that of the surrounding atom, \(\beta\). By comparison when using eam/alloy
tabulations the same \(\rho_\beta(r_{ij})\) function is used, no matter the type of the central atom. This means that the equation describing eam/fs
style models becomes:
Here a binary Al, Fe, model is being described and the resultant eam/fs
file should contain definitions for the following:
- Pair-Potentials: \(\phi_{\text{Al}\text{Al}}(r_{ij})\), \(\phi_{\text{Fe}\text{Fe}}(r_{ij})\) and \(\phi_{\text{Al}\text{Fe}}(r_{ij})\).
- Embedding-Functions: \(F_{\text{Al}}(\rho)\) and \(F_{\text{Fe}}(\rho)\).
- Density-Functions: \(\rho_{\text{Al}\text{Al}}(r_{ij})\), \(\rho_{\text{Fe}\text{Fe}}(r_{ij})\), \(\rho_{\text{Al}\text{Fe}}(r_{ij})\) and \(\rho_{\text{Fe}\text{Al}}(r_{ij})\).
From this it can be seen that, when using eam/fs
style potentials, the density functions must have both the \(\alpha\beta\) and \(\beta\alpha\) interactions specified.
Although both the \(\alpha\beta\) and \(\beta\alpha\) can be described using eam/fs
files, the Mendelev model used in this example uses the same density function for both Al-Fe and Fe-Al cross density functions [3].
Using SetFL_FS_EAMTabulation
to Tabulate the Model¶
As in previous examples it is necessary to define pair, density and embedding functions in python code that are then wrapped in EAMPotential
and Potential
objects to be passed to the tabulation function. For brevity only the names of the functions, as defined in the attached example file (eam_example3a.py
) are now given:
- Pair-Potentials:
def ppfuncAlAl(r):
- Al-Al pair-potential \(\phi_{\text{Al}\text{Al}}(r_{ij})\).def ppfuncAlFe(r):
- Al-Fe pair-potential \(\phi_{\text{Al}\text{Fe}}(r_{ij})\).def ppfuncFeFe(r):
- Fe-Fe pair-potential \(\phi_{\text{Fe}\text{Fe}}(r_{ij})\).
- Embedding-Functions:
def AlEmbedFunction(rho):
- Al embedding function \(F_{\text{Al}}(\rho)\).def FeEmbedFunction(rho):
- Fe embedding function \(F_{\text{Fe}}(\rho)\).
- Density-Functions:
def AlAlDensityFunction(r):
- Al density function \(\rho_{\text{Al}\text{Al}}(r_{ij})\).def FeFeDensityFunction(r):
- Fe density function \(\rho_{\text{Al}\text{Al}}(r_{ij})\).def FeAlDensityFunction(r):
- Al-Fe density function \(\rho_{\text{Al}\text{Fe}}(r_{ij})\).
Note
The functional forms used within the Mendelev paper [3] are somewhat long, and including their implementations here would detract from the readability of this example. However, they are included in the attached python file: eam_example3a.py
.
These functions are used within the main()
function of the eam_example3a.py
file which is now shown:
def main():
# Define list of pair potentials
pairPotentials = [
Potential('Al', 'Al', ppfuncAlAl),
Potential('Al', 'Fe', ppfuncAlFe),
Potential('Fe', 'Fe', ppfuncFeFe)]
# Assemble the EAMPotential objects
eamPotentials = [
# Al
EAMPotential('Al', 13, 26.98154, AlEmbedFunction,
{'Al': AlAlDensityFunction,
'Fe': FeAlDensityFunction},
latticeConstant=4.04527,
latticeType='fcc'),
# Fe
EAMPotential('Fe', 26, 55.845, FeEmbedFunction,
{'Al': FeAlDensityFunction,
'Fe': FeFeDensityFunction},
latticeConstant=2.855312,
latticeType='bcc')]
# Number of grid points and cut-offs for tabulation.
cutoff_rho = 300.0
nrho = 10000
cutoff = 6.5
nr = 10000
tabulation = SetFL_FS_EAMTabulation(
pairPotentials, eamPotentials, cutoff, nr, cutoff_rho, nrho)
with open("Mendelev_Al_Fe.eam.fs", "w") as outfile:
tabulation.write(outfile)
- Breaking
main()
into its components, first a list ofPotential
objects is created, this is common with the other tabulation methods already discussed:
# Define list of pair potentials
pairPotentials = [
Potential('Al', 'Al', ppfuncAlAl),
Potential('Al', 'Fe', ppfuncAlFe),
Potential('Fe', 'Fe', ppfuncFeFe)]
- Next, the
EAMPotential
objects for Al and Fe are instantiated. This is where a Finnis-Sinclair model differs from the standard EAM seen earlier. Instead of a single density callable being passed to theEAMPotential
constructor a dictionary of density functions is passed instead (see highlighted lines):
eamPotentials = [
# Al
EAMPotential('Al', 13, 26.98154, AlEmbedFunction,
{'Al': AlAlDensityFunction,
'Fe': FeAlDensityFunction},
latticeConstant=4.04527,
latticeType='fcc'),
# Fe
EAMPotential('Fe', 26, 55.845, FeEmbedFunction,
{'Al': FeAlDensityFunction,
'Fe': FeFeDensityFunction},
latticeConstant=2.855312,
latticeType='bcc')]
- The density function dictionary keys refer to the \(\beta\) species in each \(\alpha\beta\) pair. This means that:
- for the Al
EAMPotential
instance:
- \(\rho_{\text{Al}\text{Al}}\) =
AlAlDensityFunction()
,- \(\rho_{\text{Al}\text{Fe}}\) =
FeAlDensityFunction()
.- for the Fe
EAMPotential
instance:
- \(\rho_{\text{Fe}\text{Al}}\) =
FeAlDensityFunction()
,- \(\rho_{\text{Fe}\text{Fe}}\) =
FeFeDensityFunction()
.
- Finally, having defined the list of
EAMPotential
instances these are passed to the constructor ofSetFL_FS_EAMTabulation
, in this case writing the data toMendelev_Al_Fe.eam.fs
in the current directory:
tabulation = SetFL_FS_EAMTabulation(
pairPotentials, eamPotentials, cutoff, nr, cutoff_rho, nrho)
with open("Mendelev_Al_Fe.eam.fs", "w") as outfile:
tabulation.write(outfile)
Using the Mendelev_Al_Fe.eam.fs
file within LAMMPS¶
For a binary system where Al and Fe have IDs of 1 and 2 the Mendelev_Al_Fe.eam.fs
file is specified to LAMMPS as follows:
pair_style eam/fs
pair_coeff * * Mendelev_Al_Fe.eam.fs Al Fe
Example 3b: Tabulate Al-Fe Finnis-Sinclair Potentials Using TABEAM_FinnisSinclair_EAMTabulation
for DL_POLY¶
Using exactly the same model definition as for Example 3a, the Al-Fe model can be re-tabulated for DL_POLY with minimal modification to the main()
function. The modified version of the tabulation script can be found in eam_example3b.py
.
The main()
function is given below:
def main():
# Define list of pair potentials
pairPotentials = [
Potential('Al', 'Al', ppfuncAlAl),
Potential('Al', 'Fe', ppfuncAlFe),
Potential('Fe', 'Fe', ppfuncFeFe)]
# Assemble the EAMPotential objects
eamPotentials = [
# Al
EAMPotential('Al', 13, 26.98154, AlEmbedFunction,
{'Al': AlAlDensityFunction,
'Fe': FeAlDensityFunction},
latticeConstant=4.04527,
latticeType='fcc'),
# Fe
EAMPotential('Fe', 26, 55.845, FeEmbedFunction,
{'Al': FeAlDensityFunction,
'Fe': FeFeDensityFunction},
latticeConstant=2.855312,
latticeType='bcc')]
# Number of grid points and cut-offs for tabulation.
cutoff_rho = 300.0
nrho = 10000
cutoff = 6.5
nr = 10000
tabulation = TABEAM_FinnisSinclair_EAMTabulation(
pairPotentials, eamPotentials, cutoff, nr, cutoff_rho, nrho)
with open("TABEAM", "w") as outfile:
tabulation.write(outfile)
Excluding the import statement at the top of the file, only two lines have been changed (highlighted). The first changes the filename to TABEAM
whilst the second tells python to create an object of TABEAM_FinnisSinclair_EAMTabulation
instead of atsim.potentials.eam_tabulation.SetFL_FS_EAMTabulation
:
tabulation = TABEAM_FinnisSinclair_EAMTabulation(
pairPotentials, eamPotentials, cutoff, nr, cutoff_rho, nrho)
with open("TABEAM", "w") as outfile:
tabulation.write(outfile)
That’s it, nothing else has changed.
Using the TABEAM
file with DL_POLY¶
Running eam_example3b.py
produces a file names TABEAM
within the working directory. This should be placed in the same directory as the other DL_POLY input files (CONTROL
, CONFIG
and FIELD
). Then the following should be added to the end of the FIELD
file:
metal 3
Al Al eeam
Fe Fe eeam
Al Fe eeam
Note
The Extended EAM (eeam) variant of the TABEAM
file generated here is only supported in DL_POLY versions >= 4.05.
Footnotes
[1] | A.P. Sutton, and J. Chen, “Long-range Finnis-Sinclair potentials”, Philos. Mag. Lett. 61 (1990) 139 doi:10.1080/09500839008206493. |
[2] |
|
[3] | (1, 2, 3) M.I. Mendelev, D.J. Srolovitz, G.J. Ackland, and S. Han, “Effect of Fe Segregation on the Migration of a Non-Symmetric Σ5 Tilt Grain Boundary in Al”, J. Mater. Res. 20 (2011) 208. |