isolated_atom calculation style

Lucas M. Hale, lucas.hale@nist.gov, Materials Science and Engineering Division, NIST.

Introduction

The isolated_atom calculation style evaluates the base energies of all atomic models associated with an interatomic potential. For some potentials, the isolated energy values are necessary to properly compute the cohesive energy of crystal structures. This also provides a simple test whether a potential implementation is compatible with a version of LAMMPS.

Version notes

  • 2020-09-22: Notebook first added.

Additional dependencies

Disclaimers

  • NIST disclaimers

  • Some potentials have two cutoffs with atomic energies outside the first being the “isolated” energy while outside the second have zero energy. The first isolated energy values for those potentials can be found using the diatom_scan calculation instead.

Method and Theory

The calculation loops over all symbol models of the potential and creates a system with a single particle inside a system with non-periodic boundary conditions. The potential energy of each unique isolated atom is evaluated without relaxation/integration.

The cohesive energy, \(E_{coh}\), of a crystal structure is given as the per-atom potential energy of the crystal structure at equilibrium \(E_{crystal}/N\) relative to the potential energy of the same atoms infinitely far apart, \(E_i^{\infty}\)

\[E_{coh} = \frac{E_{crystal} - \sum{N_i E_{i}^{\infty}}}{N},\]

Where the \(N_i\) values are the number of each species \(i\) and \(\sum{N_i} = N\).

For most potentials, \(E_i^{\infty}=0\) meaning that the measured potential energy directly corresponds to the cohesive energy. However, this is not the case for all potentials as some have offsets either due to model artifacts or because it allowed for a better fitted model.

Demonstration

1. Setup

1.1. Library imports

Import libraries needed by the Notebook. The external libraries used are:

[1]:
# Standard library imports
import os
from pathlib import Path
import datetime

# http://www.numpy.org/
import numpy as np

# https://github.com/usnistgov/atomman
import atomman as am
import atomman.lammps as lmp
import atomman.unitconvert as uc

# https://github.com/usnistgov/iprPy
import iprPy

print('Notebook last executed on', datetime.date.today(), 'using iprPy version', iprPy.__version__)
Notebook last executed on 2020-09-22 using iprPy version 0.10.2

1.2. Default calculation setup

[2]:
# Specify calculation style
calc_style = 'isolated_atom'

# If workingdir is already set, then do nothing (already in correct folder)
try:
    workingdir = workingdir

# Change to workingdir if not already there
except:
    workingdir = Path('calculationfiles', calc_style)
    if not workingdir.is_dir():
        workingdir.mkdir(parents=True)
    os.chdir(workingdir)

# Initialize connection to library
library = iprPy.Library(load=['lammps_potentials'])

2. Assign values for the calculation’s run parameters

2.1. Specify system-specific paths

  • lammps_command is the LAMMPS command to use (required).

  • mpi_command MPI command for running LAMMPS in parallel. A value of None will run simulations serially.

[3]:
lammps_command = 'lmp_serial'
mpi_command = None

2.2. Load interatomic potential

  • potential_name gives the name of the potential_LAMMPS reference record in the iprPy library to use for the calculation.

  • potential is an atomman.lammps.Potential object (required).

[4]:
potential_name = '1999--Mishin-Y--Ni--LAMMPS--ipr1'

# Retrieve potential and parameter file(s)
potential = library.get_lammps_potential(id=potential_name, getfiles=True)

3. Define calculation function(s) and generate template LAMMPS script(s)

3.1 run0.template

[5]:
with open('run0.template', 'w') as f:
    f.write("""#LAMMPS input script that evaluates a system's energy without relaxing

<atomman_system_pair_info>

thermo_style custom step pe
thermo_modify format float %.13e

run 0""")

3.2 isolated_atom()

[6]:
def isolated_atom(lammps_command, potential, mpi_command=None):
    """
    Evaluates the isolated atom energy for each elemental model of a potential.

    Parameters
    ----------
    lammps_command :str
        Command for running LAMMPS.
    potential : atomman.lammps.Potential
        The LAMMPS implemented potential to use.
    mpi_command : str, optional
        The MPI command for running LAMMPS in parallel.  If not given, LAMMPS
        will run serially.

    Returns
    -------
    dict
        Dictionary of results consisting of keys:

        - **'energy'** (*dict*) - The computed potential energies for each
          symbol.
    """
    # Build filedict if function was called from iprPy
    try:
        assert __name__ == pkg_name
        calc = iprPy.load_calculation(calculation_style)
        filedict = calc.filedict
    except:
        filedict = {}

    # Initialize dictionary
    energydict = {}

    # Initialize single atom system
    box = am.Box.cubic(a=1)
    atoms = am.Atoms(atype=1, pos=[[0.5, 0.5, 0.5]])
    system = am.System(atoms=atoms, box=box, pbc=[False, False, False])

    # Get lammps units
    lammps_units = lmp.style.unit(potential.units)

    # Define lammps variables
    lammps_variables = {}

    # Loop over symbols
    for symbol in potential.symbols:
        system.symbols = symbol

        # Add charges if required
        if potential.atom_style == 'charge':
            system.atoms.prop_atype('charge', potential.charges(system.symbols))

        # Save configuration
        system_info = system.dump('atom_data', f='isolated.dat',
                                  potential=potential,
                                  return_pair_info=True)
        lammps_variables['atomman_system_pair_info'] = system_info

        # Write lammps input script
        template_file = 'run0.template'
        lammps_script = 'run0.in'
        template = iprPy.tools.read_calc_file(template_file, filedict)
        with open(lammps_script, 'w') as f:
            f.write(iprPy.tools.filltemplate(template, lammps_variables,
                                             '<', '>'))

        # Run lammps and extract data
        output = lmp.run(lammps_command, lammps_script, mpi_command)
        energy = output.simulations[0]['thermo'].PotEng.values[-1]
        energydict[symbol] = uc.set_in_units(energy, lammps_units['energy'])

    # Collect results
    results_dict = {}
    results_dict['energy'] = energydict

    return results_dict

4. Run calculation function(s)

[7]:
results_dict = isolated_atom(lammps_command, potential, mpi_command=mpi_command)
[8]:
results_dict.keys()
[8]:
dict_keys(['energy'])

5. Report results

5.1. Define units for outputting values

  • energy_unit is the unit of energy to display values in.

[9]:
energy_unit = 'eV'

5.2 Display isolated energies

[10]:
for symbol, energy in results_dict['energy'].items():
    print(symbol, uc.get_in_units(energy, energy_unit), energy_unit)
Ni -3.0970029925997e-11 eV