MeamPairPotential

class MeamPairPotential(particleType1, particleType2, latticeType, nearestNeighbors, alpha, referenceDistance, referenceEnergy, attrac=0.0, repuls=0.0, nn2=False, zbl=True)

Constructor of the potential.

Parameters:
  • particleType1 (ParticleType or ParticleIdentifier) – Identifier of the particle type.
  • particleType2 (ParticleType or ParticleIdentifier) – Identifier of the particle type.
  • latticeType (str) – Type of the reference structure.
  • nearestNeighbors (int) – Number of nearest neighbors in the reference structure.
  • alpha (float) – Exponential decay factor for pair potential.
  • referenceDistance (PhysicalQuantity of type length) – The equilibrium nearest neighbor distance..
  • referenceEnergy (PhysicalQuantity of type energy) – The cohesive energy.
  • attrac (float) – The attractive factor in the pair potential.
  • repuls (float) – The repulsive factor in the pair potential.
  • nn2 (bool) – Triggers second nearest neighbor formulation for symmetric pair.
  • zbl (bool) – Toggles zbl blending for small distances.
getAllParameterNames()

Return the names of all used parameters as a list.

getAllParameters()

Return all parameters of this potential and their current values as a <parameterName / parameterValue> dictionary.

static getDefaults()

This potential does not really have useful defaults, so I’ll just make some up for testing only

getParameter(parameterName)

Get the current value of the parameter parameterName.

setAlpha(alpha)

Set the parameter alpha.

Parameters:alpha (float) – Exponential decay factor for pair potential.
setAttrac(attrac)

Set the parameter attrac.

Parameters:attrac (float) – The attractive factor in the pair potential.
setCutoff(r_cut)

Set the cutoff radius for this potential.

Parameters:r_cut (PhysicalQuantity of type length) – The cutoff radius of this potential.
setLatticeType(latticeType)

Set the reference lattice type.

Parameters:latticeType (str) – Type of the reference structure.
setNN2(nn2)

Set the parameter nn2.

Parameters:nn2 (bool) – Triggers second nearest neighbor formulation for symmetric pair.
setNearestNeighbors(nearestNeighbors)

Set the parameter nearestNeighbors.

Parameters:nearestNeighbors (int) – Number of nearest neighbors in the reference structure.
setParameter(parameterName, value)

Set the parameter parameterName to the given value.

Parameters:
  • parameterName (str) – The name of the parameter that will be modified.
  • value – The new value that will be assigned to the parameter parameterName.
setReferenceDistance(distance)

Set the parameter distance.

Parameters:distance (PhysicalQuantity of type length) – The equilibrium nearest neighbor distance.
setReferenceEnergy(energy)

Set the parameter energy.

Parameters:energy (PhysicalQuantity of type energy) – The cohesive energy.
setRepuls(repuls)

Set the parameter repuls.

Parameters:repuls (float) – The repulsive factor in the pair potential.
setZBL(zbl)

Set the parameter ZBL.

Parameters:zbl (bool) – Toggles zbl blending for small distances.

Usage Examples

Define a MEAM potential for a SiO2 quartz crystal by adding particle types and interaction functions to the TremoloXPotentialSet.

# Set up lattice
vector_a = [5.4306, 0.0, 0.0]*Angstrom
vector_b = [0.0, 5.4306, 0.0]*Angstrom
vector_c = [0.0, 0.0, 5.4306]*Angstrom
lattice = UnitCell(vector_a, vector_b, vector_c)

# Define elements
elements = [Silicon, Silicon, Silicon, Silicon, Silicon, Silicon, Silicon,
            Silicon]

# Define coordinates
fractional_coordinates = [[ 0.  ,  0.  ,  0.  ],
                          [ 0.25,  0.25,  0.25],
                          [ 0.5 ,  0.5 ,  0.  ],
                          [ 0.75,  0.75,  0.25],
                          [ 0.5 ,  0.  ,  0.5 ],
                          [ 0.75,  0.25,  0.75],
                          [ 0.  ,  0.5 ,  0.5 ],
                          [ 0.25,  0.75,  0.75]]

# Set up configuration
bulk_configuration = BulkConfiguration(
    bravais_lattice=lattice,
    elements=elements,
    fractional_coordinates=fractional_coordinates
    )

# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------
potentialSet = TremoloXPotentialSet(name='MEAM_Si_2007')
potentialSet.addParticleType(ParticleType(
    symbol='Si',
    mass=28.0855*atomic_mass_unit,
    charge=None,
    sigma=None,
    sigma14=None,
    epsilon=None,
    epsilon14=None,
    atomicNumber=14,
    tags=[],
))


_option = MeamGlobalOption(
    delr=0.1*Angstrom,
    erose=2,
    wf_mixing=2,
    r_cut=4.5*Angstrom,
    augment_1st=0.0,
    embedding_negative=False,
    density_scaling=False,
)
potentialSet.addOption(_option)
_potential = MeamElementPotential(
    particleType=ParticleIdentifier('Si', []),
    latticeType='dia',
    nearestNeighbors=4,
    alpha=4.89890486934,
    beta=[3.55, 2.5, 0.0, 7.5, ],
    referenceDistance=2.35*Angstrom,
    referenceEnergy=4.63*eV,
    scalingFactor=0.58,
    weightingFactors=[1.8, 5.25, -2.61, ],
    rho=1.0,
    gamma=3,
    attrac=0.0,
    repuls=0.0,
    nn2=True,
    zbl=False,
)
potentialSet.addPotential(_potential)
_potential = MeamScreeningPotential(
    particleType1=ParticleIdentifier('Si', []),
    particleType2=ParticleIdentifier('Si', []),
    particleType3=ParticleIdentifier('Si', []),
    Cmin=1.41,
    Cmax=2.8,
)
potentialSet.addPotential(_potential)

calculator = TremoloXCalculator(parameters=potentialSet)
calculator.setVerletListsDelta(0.25*Angstrom)

bulk_configuration.setCalculator(calculator)
bulk_configuration.update()

Notes

The modified embedded-atom-method (MEAM) potential extends the conventional embedded-atom-method by including directional bonding [Bas92]. This allows the application to both metallic and covalent systems, as well as mixed systems.

The energy function for the Modified Embedded Atom Method as presented e.g. in [GWS03] is given as the sum of some pair term and some embedding term as follows:

\[E= \sum_i \left\{ F_i(\bar{\rho}_i) + \frac{1}{2} \sum_{j\ne i} \Phi_{ij}(r_{ij})\right\}\]

The embedding function is of the form

\[F_i(\bar{\rho}_i) = A_i E_i^0\bar{\rho}_i \log \bar{\rho}_i\]

for \(\bar{\rho}_i > 0\) and for \(\bar{\rho}_i \le 0\)

\[F_i(\bar{\rho}_i) = 0\]

if embedding_negative is 0 or

\[F_i(\bar{\rho}_i) = -A_i E_i^0 \bar{\rho}_i\]

if it is 1, where \(A_i\) and \(E_i^0\) are element dependent parameters. The background electron density \(\bar{\rho}_i\) is given by:

\[\bar{\rho}_i = \frac{\bar{\rho}_i^{(0)}}{\rho_i^0} G_i (\Gamma_i)\]

where

\[\Gamma_i = \sum_{k=1}^3 t_i^{(k)}\left(\frac{\bar{\rho}_i^{(k)}}{\bar{\rho}_i^{(0)}}\right)^2\]

and \(G_i\) is one of

(1)\[\begin{split}G_i(\Gamma) = \begin{cases} \sqrt{1+\Gamma} & gamma=0\\ \exp(\Gamma / 2) & gamma=1\\ \text{sign}(1+\Gamma) \sqrt{|1+\Gamma|} & gamma=2\\ \frac{2}{1+\exp(-\Gamma)} &gamma=3\\ \sqrt{1+\Gamma} & gamma=4 \end{cases}\end{split}\]

With the element dependent constant \(\rho_i^0\), the number of nearest neighbors in the reference structure \(Z_{i0}\) and \(\Gamma_i^\text{ref}\) \(\Gamma\) evaluated in the reference structure, the composition dependent density scaling can be chosen to be either

(2)\[\rho_i^0 = \rho_{i0}Z_{i0}G_i(\Gamma_i^\text{ref})\]

or

(3)\[\rho_i^0 = \rho_{i0} Z_{i0}\]

Note that for gamma=0 and gamma=2 \(G_i(\Gamma_i^\text{ref})\) is always set to 1. The partial electron densities are defined as

\[\begin{split}\begin{eqnarray*} \bar{\rho}_i^{(0)} &=& \sum_{j\ne i} \rho_j^{a(0)}(r_{ij})S_{ij} \\ (\bar{\rho}_i^{(1)})^2 &=& \sum_{\alpha = 1}^3 \left[\sum_{j \ne i} \rho_j^{a(1)} \frac{r_{ij\alpha}}{r_{ij}} S_{ij} \right]^2 \\ (\bar{\rho}_i^{(2)})^2 &=& \sum_{\alpha =1}^3 \sum_{\beta =1}^3 \left[\sum_{j\ne i} \rho_j^{a(2)} \frac{r_{ij\alpha}r_{ij\beta}}{r_{ij}^2} S_{ij} \right]^2 - \frac{1}{3} \left[ \sum_{j\ne i} \rho_j^{a(2)}S_{ij}\right]^2\\ (\bar{\rho}_i^{(3)})^2 &=& \sum_{\alpha =1}^3 \sum_{\beta =1}^3 \sum_{\gamma =1}^3 \left[\sum_{j\ne i} \rho_j^{a(3)} \frac{r_{ij\alpha}r_{ij\beta}r_{ij\gamma}}{r_{ij}^3} S_{ij}\right]^2 - \frac{3}{5}\sum_{\alpha =1}^3\left[\sum_{j \ne i} \rho_j^{a(3)} \frac{r_{ij\alpha}}{r_{ij}} S_{ij}\right]^2 \end{eqnarray*}\end{split}\]

with the atomic electron densities given by

\[\rho_i^{a(k)}(r_{ij}) = \rho_{i0} \exp \left[-\beta_i^{(k)} \left(\frac{r_{ij}}{r_i^0} - 1 \right)\right]\]

Finally, the average weighting factors are given by

(4)\[\begin{split}t_i^{(k)} = \begin{cases} \frac{1}{\bar{\rho}_i^{(0)}} \sum_{j\ne i} t_{0,j}^{(k)} \rho_j^{a(0)} S_{ij} & \text{wf-mixing=0} \\ \frac{ \sum_{j\ne i} t_{0,j}^{(k)} \rho_j^{a(0)} S_{ij}} { \sum_{j\ne i} \left(t_{0,j}^{(k)}\right)^2 \rho_j^{a(0)} S_{ij}}& \text{wf-mixing=1} \\ t_{0,i}^{(k)} & \text{wf-mixing=2} \end{cases}\end{split}\]

with element dependent parameters \(t_{0,i}^{(k)}\).

The pair term \(\Phi_{ij}\) is defined as

(5)\[\begin{split}\begin{eqnarray} \Phi_{ij}(r_{ij}) &=& \bar{\Phi}_{ij}(r_{ij}) S_{ij}\\ \bar{\Phi}_{ij}(r_{ij}) &=& \frac{1}{Z_{ij0}} [2E_i^u(r_{ij}) - F_i(\hat{\rho}_i (r_{ij})) - F_j(\hat{\rho}_j(r_{ij})] \end{eqnarray}\end{split}\]

where \(Z_{ij0}\) is the number of nearest neighbors in the reference structure, \(F_i(\hat{\rho}_{ij})\) the embedding function evaluated in the reference structure and

(6)\[\begin{split}E_i^u(r_{ij}) = \begin{cases} -E_{ij}^c \left(1+a_{ij}^{*} + a_{ij}^{(3)} {a_{ij}^{*}}^3 \frac{r^0_{ij}}{r_{ij}}\right) \exp(-a_{ij}^*) & \text{erose=0} \\ -E_{ij}^c \left(1+a_{ij}^* + \left(-\text{attrac}+\frac{\text{repuls}}{r_{ij}}\right){a_{ij}^*}^3\right)\exp(-a_{ij}^*) & \text{erose=1}\\ -E_{ij}^c\left(1+a_{ij}^* + a_{ij}^{(3)}{a_{ij}^*}^3\right)\exp(-a_{ij}^*) & \text{erose=2} \end{cases}\end{split}\]

with

\[a_{ij}^*(r_{ij}) = \alpha_{ij} \left(\frac{r_{ij}}{r^0_{ij}} - 1\right)\]

and

\[\begin{split}a_{ij}^{(3)} = \begin{cases} \text{repuls} & a_{ij}^* < 0\\ \text{attrac} & a_{ij}^* \ge 0 \end{cases}\end{split}\]

Some reference structures use a more complex form of (5). For small \(r_{ij}\), this pair term is blended with a ZBL potential by default.newline Last, the screening function \(S_{ij}\) is defined by

\[\begin{split}\begin{eqnarray*} S_{ij}&=&\bar{S}_{ij} f_c\left(\frac{r_c - r_{ij}}{\Delta r}\right)\\ \bar{S}_{ij} &=& \prod_{k\ne i,j} S_{ikj}\\ S_{ikj} &=& f_c\left(\frac{C_{ikj}-C_{\text{min},ikj}}{C_{\text{max},ikj} -C_{\text{min},ikj}}\right)\\ C_{ikj}&=& 1+2\frac{r_{ij}^2 r_{ik}^2 + r_{ij}^2r_{jk}^2-r_{ij}^4}{r_{ij}^4-(r_{ik}^2-r_{jk}^2)^2}\\ f_c(x)&=&\begin{cases} 1 & x \ge 1 \\ [1-(1-x)^4]^2 & 0<x<1\\ 0 & x \le 0 \end{cases} \end{eqnarray*}\end{split}\]

where \(r_c\) is the global cutoff radius, \(\Delta r\) is the length of a smoothing region for \(r\) near \(r_c\) and \(C_\text{min}\) and \(C_\text{max}\) are element-triple dependent parameters.

Note that the cutoff radius for this term needs to be larger than \(r_c\), depending on \(C_{\text{max}}\), because even triples with some particle distance \(> r_c\) may have a screening effect.

The parameter lattice defines the reference structure for the given particle type. The parameters alpha, beta_k, re, Ec, scaling_factor, weighting_factor_k, rho, gamma, attrac and repuls correspond to the variables \(\alpha\), \(\beta^{(k)}\), \(r^0\), \(E^c\), \(A\), \(t_0^{(k)}\), \(\rho_{0}\), attrac and repuls from the previous section respectively. nn2 enables the second nearest neighbor formulation as described in cite{lee2000second}, zbl enables ZBL blending for small distances if set to 1.

re and Ec correspond to \(r^0\) and \(E^c\). Note that all element_data entries must be given before the first element_pair entry.newline If augment_1st is set to 1, we set

\[t_{0,i}^{(1)} = t_{0,i}^{(1)} + \frac{3}{5} t_{0,i}^{(3)}\]

for all \(i\). This is usually only needed for older parameter sets, modern ones already include this correction. The default density_scaling uses (2), setting it to 1 uses (3) instead.

The latticeType parameter must be one of the following strings:

  • dia = diamond
  • fcc = face centered cubic
  • bcc = body centered cubic
  • dim = dimer
  • b1 = rock salt
  • hcp = hexagonal close-packed
  • c11 = MoSi2 structure
  • l12 = Cu3Au structure
  • b2 = CsCl structure

[Bas92]M. I. Baskes. Modified embedded-atom potentials for cubic materials and impurities. Phys. Rev. B, 46:2727–2742, Aug 1992. URL: http://link.aps.org/doi/10.1103/PhysRevB.46.2727, doi:10.1103/PhysRevB.46.2727.
[GWS03]PM Gullet, G Wagner, and A Slepoy. Numerical tools for atomistic simulations. SANDIA Report, 8782:2003, 2003.