ZBLPotential

class ZBLPotential(particleType1, particleType2, A, b, r_i=None, r_cut=None)

Constructor of the potential.

Parameters:
  • particleType1 (ParticleType or ParticleIdentifier) – Identifier of the first particle type.
  • particleType2 (ParticleType or ParticleIdentifier) – Identifier of the second particle type.
  • A (PhysicalQuantity of type energy * length) – Potential parameter.
  • b (PhysicalQuantity of type length) – Potential parameter.
  • r_i (PhysicalQuantity of type length) – The inner cutoff radius where the smoothing of the potential starts. r_i must be smaller than r_cut.
  • r_cut (PhysicalQuantity of type length) – The cutoff radius of this potential.
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()

Get the default parameters of this potential and return them in form of a dictionary of <parameter name, default value> key-value pairs.

getParameter(parameterName)

Get the current value of the parameter parameterName.

setA(A)

Set the parameter A.

Parameters:A (PhysicalQuantity of type energy * length) – Potential parameter.
setCutoff(r_cut)

Set the cutoff radius for this potential.

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

Set the inner cutoff radius for this potential.

Parameters:r_i (PhysicalQuantity of type length) – The inner cutoff radius (distance where the spline-smoothing starts).
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.
setb(b)

Set the parameter b.

Parameters:b (PhysicalQuantity of type length) – Potential parameter.

Usage Examples

Define a Buckingham-type potential for titanium oxide with additional ZBL repulsive term between titanium and oxygen atoms by adding particle types and interaction functions to the TremoloXPotentialSet.

# Set up lattice
lattice = SimpleTetragonal(4.593*Angstrom, 2.959*Angstrom)

# Define elements
elements = [Titanium, Titanium, Oxygen, Oxygen, Oxygen, Oxygen]

# Define coordinates
fractional_coordinates = [[ 0.    ,  0.    ,  0.    ],
                          [ 0.5   ,  0.5   ,  0.5   ],
                          [ 0.3051,  0.3051,  0.    ],
                          [ 0.6949,  0.6949,  0.    ],
                          [ 0.8051,  0.1949,  0.5   ],
                          [ 0.1949,  0.8051,  0.5   ]]

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

# -------------------------------------------------------------
# Calculator
# -------------------------------------------------------------

potentialSet = TremoloXPotentialSet(name='Trinastic_HfOSiTaTi_2013')
potentialSet.addParticleType(ParticleType(
    symbol='O',
    mass=15.9994*atomic_mass_unit,
    charge=-1.2
))
potentialSet.addParticleType(ParticleType(
    symbol='Ti',
    mass=47.867*atomic_mass_unit,
    charge=2.4
))

potential = BuckinghamPotential(
    particleType1='O',
    particleType2='O',
    A=1388.77*eV,
    rho=0.3623*Angstrom,
    r_i=6.0*Angstrom,
    r_cut=7.5*Angstrom
)
potentialSet.addPotential(potential)
potential = LennardJonesMNPotential(
    particleType1='O',
    particleType2='O',
    A=0.0*Angstrom**12*eV,
    B=175.0*eV*Angstrom**6,
    m=12.0,
    n=6.0,
    r_cut=10.0*Angstrom
)
potentialSet.addPotential(potential)
potential = BuckinghamPotential(
    particleType1='O',
    particleType2='Ti',
    A=5105.12*eV,
    rho=0.2253*Angstrom,
    r_i=6.0*Angstrom,
    r_cut=7.5*Angstrom
)
potentialSet.addPotential(potential)
potential = LennardJonesMNPotential(
    particleType1='O',
    particleType2='Ti',
    A=0.0*Angstrom**12*eV,
    B=20.0*eV*Angstrom**6,
    m=12.0,
    n=6.0,
    r_cut=10.0*Angstrom
)
potentialSet.addPotential(potential)
potential = MorsePotential(
    particleType1='O',
    particleType2='Ti',
    r_0=1.96*Angstrom,
    k=1.9*1/Angstrom,
    E_0=0.3478*eV,
    r_i=6.5*Angstrom,
    r_cut=7.0*Angstrom
)
potentialSet.addPotential(potential)

# Add the ZBL-potential
# Calculate the ZBL-prefactor A=e**2/(4*pi*eps_0)*Zi*Zj
ZBL_prefactor = elementary_charge**2/(4.0*numpy.pi*vacuum_permitivity)*Titanium.atomicNumber()*Oxygen.atomicNumber()
ZBL_length = 0.8854*Bohr/(Titanium.atomicNumber()**0.23 + Oxygen.atomicNumber()**0.23)
potential = ZBLPotential(
    particleType1='Ti',
    particleType2='O',
    A=ZBL_prefactor,
    b=ZBL_length,
    r_i=0.7*Angstrom,
    r_cut=2.0*Angstrom
)
potentialSet.addPotential(potential)

potentialSet.setCoulombSolver(CoulombDSF(r_cut=9.0*Angstrom, alpha=0.2))
calculator = TremoloXCalculator(parameters=potentialSet)
calculator.setInternalOrdering("default")
calculator.setVerletListsDelta(0.25*Angstrom)

bulk_configuration.setCalculator(calculator)
bulk_configuration.update()

Notes

This potential implements the universal Ziegler-Biersack-Littmark (ZBL) repulsive potential [ZB85]

\[\begin{split}V^{\rm ZBL} = \sum_i \sum_{j > i} \frac{1}{4\pi\epsilon_0} \frac{Z_i Z_j}{r_{ij}} \phi\left(\frac{r_{ij}}{a_{ij}}\right)\end{split}\]

where

\[\phi(x) = 0.1818e^{-3.2 x}+0.5099e^{-0.9423 x}+0.2802e^{-0.4029 x}+0.02817e^{-0.2016 x}\]

and

\[a_{ij} = \frac{0.8854a_0}{Z_i^{0.23} + Z_j^{0.23}} \, .\]

In TremoloX, the ZBL-potential between a pair of particle types i and j is implemented as follows

\[V^{\rm ZBL}_{ij}(r_{ij}) = \frac{A_{ij}}{r_{ij}} \phi\left(\frac{r_{ij}}{b_{ij}}\right)\]

where the constructor arguments A and b correspond to \(A_{ij}\) and \(b_{ij}\).

To obtain the original ZBL-potential, set

\[A_{ij} = \frac{Z_i Z_j}{4\pi\epsilon_0}\]

and

\[b_{ij} = a_{ij} = \frac{0.8854a_0}{Z_i^{0.23} + Z_j^{0.23}} \, .\]

At long distances this potential is brought smoothly to zero between the inner cutoff r_i and the outer cutoff r_cut, using a 5th order spline function. This ensures that there are no discontinuites in the forces as atoms are brought closer together, and that the energy is properly conserved. The applied potential \(U(r)\) is given as:

\[U(r) = V_{ij}(r) \times S(r)\]

where \(V_{ij}(r)\) is the pair potential and \(S(r)\) is the spline function. The values of the spline function are:

  • 1 when \(r \le r_i\)
  • In the range \([0,1]\) when \(r_i < r < r_{cut}\)
  • 0 when \(r \ge r_{cut}\)