ONETEP

Introduction

ONETEP is a linear-scaling density functional theory code which exploit the near-sightness of the electronic density. It uses a set of atom-centered local orbitals (denoted NGWFs) which are optimised in situ to enable calculations with a minimal number of orbitals.

This interface makes it possible to use ONETEP as a calculator in ASE. You need to have a copy of the ONETEP code (and an appropriate license) to use this interface.

Environment variables

The environment variable ASE_ONETEP_COMMAND must hold the command to invoke the ONETEP calculation. The variable must be a string with a link to the ONETEP binary, and any other specific settings required for your environment (srun, mpirun, …)

You can setup this environment variable in your shell configuration file:

$ export ASE_ONETEP_COMMAND="export OMP_NUM_THREADS=4; mpirun -n 6 ~/onetep/bin/onetep.arch"

Or within python itself:

>>> environ["ASE_ONETEP_COMMAND"]="export OMP_NUM_THREADS=4; mpirun -n 6 ~/onetep/bin/onetep.arch"

ASE will automatically redirect stdout and stderr to the appropriate files, namely “$LABEL.out” and “$LABEL.err” where label is the name used for your ONETEP calculations

Pseudopotentials

ONETEP accepts PAW datasets and NCP pseudpotentials with formats USP and recpot. Pseudopotentials are passed directly to the Onetep calculator as a dictionnary definition. If no pseudopotentials are passed ASE will try to guess the files based on the element used and the pseudo_path variable.

# Explicitly providing each path
calc = Onetep(pseudopotentials = {'H': '/path/to/pseudos/H.usp', 'O': '/path/to/pseudos/O.usp'})
# Using pseudo_path
calc = Onetep(pseudo_path = '/path/to/pseudos', pseudopotentials = {'H': 'H.usp', 'O': 'O.usp'})
# ASE will try to guess them
calc = Onetep(pseudo_path = '/path/to/pseudos')

For ASE to correctly guess the pseudopotentials use a pseudo_path that contains only one kind of pseudopotential per element.

ONETEP Calculator

class ase.calculators.onetep.Onetep(*, profile=None, directory='.', parallel_info=None, parallel=True, **kwargs)[source]

Class for the ONETEP calculator, uses ase/io/onetep.py. Need the env variable “ASE_ONETEP_COMMAND” defined to properly work. All other options are passed in kwargs.

Parameters:
  • autorestart (Bool) – When activated, manages restart keywords automatically.

  • append (Bool) – Append to output instead of overwriting.

  • directory (str) – Directory where to run the calculation(s).

  • keywords (dict) – Dictionary with ONETEP keywords to write, keywords with lists as values will be treated like blocks, with each element of list being a different line.

  • xc (str) – DFT xc to use e.g (PBE, RPBE, …).

  • ngwfs_count (int|list|dict) –

    Behaviour depends on the type:

    int: every species will have this amount of ngwfs. list: list of int, will be attributed alphabetically to species: dict: keys are species name(s), value are their number:

  • ngwfs_radius (int|list|dict) –

    Behaviour depends on the type:

    float: every species will have this radius. list: list of float, will be attributed alphabetically to species: [10.0, 9.0] dict: keys are species name(s), value are their radius: {‘Na’: 9.0, ‘Cl’: 10.0}

  • pseudopotentials (list|dict) –

    Behaviour depends on the type:

    list: list of string(s), will be attributed alphabetically to specie(s): [‘Cl.usp’, ‘Na.usp’] dict: keys are species name(s) their value are the pseudopotential file to use: {‘Na’: ‘Na.usp’, ‘Cl’: ‘Cl.usp’}

  • pseudo_path (str) –

    Where to look for pseudopotential, correspond to the pseudo_path keyword of ONETEP.

    Note

    write_forces is always turned on by default when using this interface.

    Note

    Little to no check is performed on the keywords provided by the user via the keyword dictionary, it is the user responsibility that they are valid ONETEP keywords.

Simple calculations can be setup calling the Onetep calculator without any parameters, in this case ONETEP default parameters will be used. For more complex cases using the keywords parameters is necessary. The ‘keywords’ parameters is a dictionnary in which each key is a string that represent a ONETEP keywords.

Examples

Here is an example of setting up a calculation on a water molecule:

from ase.build import molecule
from ase.calculators.onetep import Onetep
from os import environ

# water molecule from ASE database, centered in a ~ 24 Å box
wat = molecule('H2O')
wat.center(12)
environ["ASE_ONETEP_COMMAND"]="export OMP_NUM_THREADS=8; mpirun -n 2 ~/onetep/bin/onetep.arch"
# Ouput will be in "water.out"
calc = Onetep(label = 'water', xc = 'PBE', paw = True, pseudo_path = '/path/to/pseudos')
wat.calc = calc
wat.get_potential_energy()

Another more complex example on Pt13:

from os import environ

import numpy as np

from ase.build import molecule
from ase.calculators.onetep import Onetep
from ase.cluster import Octahedron
from ase.optimize.sciopt import SciPyFminBFGS
# Pt13 from ase.cluster
nano = Octahedron('Pt', 3, 1)
nano.set_cell(np.eye(3)*12)
nano.center()

label = 'pt13'

environ["ASE_ONETEP_COMMAND"]="export OMP_NUM_THREADS=8; mpirun -n 8 ~/onetep/bin/onetep.arch"

# ONETEP default are atomic units, one can specify 'cutoff_energy' : '600 eV' if needed.
keywords = {
    'xc' : 'rpbe',
    'do_properties' : True,
    'cutoff_energy' : 35,
    'output_detail': 'verbose',
    'elec_energy_tol': 1.0e-5/len(atoms),
}

# Ouput will be in "pt13.out",
# append = True will not overwrite file at each step
calc = Onetep(
    label = label,
    edft = True,
    append = True,
    pseudo_path = '/path/to/pseudos',
    keywords = keywords)

nanoparticle.calc = calc

opt = SciPyFminBFGS(atoms = nano, trajectory = label + ".traj", logfile = label + ".log")
opt.run(fmax=0.01)

Here is an example of setting up an EELS and LDOS calculations on an N-substituted graphene sheet, demonstrating several more advanced functionalities (eg tags, species groups, and overrides to pseudopotentials and atomic solver strings):

import numpy as np

from ase.build import graphene_nanoribbon
from ase.calculators.onetep import Onetep
from ase.io import write
from numpy.linalg import norm
from numpy.random import choice

sheet = graphene_nanoribbon(10, 10, type='zigzag', vacuum = 10)

# Get all distances to center of mass
com = sheet.get_center_of_mass()
distances_to_com = norm(sheet.positions - com, axis = 1)

# Find atoms close to com and change one randomly to N
p, = np.where(distances_to_com < 5)
to_nitro = choice(p)
sheet[to_nitro].symbol = 'N'

shell_rad = np.array([1.5, 2.5, 3.0, 4.0, 4.5])

tags = np.zeros(len(sheet), dtype=np.int32)

# We want to tag atoms that are close to the introduced nitrogen
for idx, rad in enumerate(reversed(shell_rad)):
    # All distances N-C
    dist = norm(sheet[to_nitro].position - sheet.get_positions(), axis = 1)
    # Which ones are closest to rad?
    p, = np.where(dist < rad)
    # Cannot be the nitrogen itself
    p = p[p != to_nitro]
    # Tags them
    tags[p] = len(shell_rad) - idx

sheet.set_tags(tags)

tags = ['' if i == 0 else i for i in tags]

species = np.unique(np.char.add(sheet.get_chemical_symbols(), tags))

keywords = {
    'species_core_wf' : ['N /path/to/pseudo/corehole.abinit'],
    'species_solver' : ['N SOLVE conf=1s1 2p4'],
    'pseudo_path': '/Users/tomdm/PseudoPotentials/SSSP_1.2.1',
    'xc' : 'PBE',
    'paw': True,
    'do_properties': True,
    'cutoff_energy' : '500 eV',
    'species_ldos_groups': species,
    'task' : 'GeometryOptimization'
}

calc = Onetep(
    label = 'N_doped_graphene_001',
    keywords = keywords
)

# Checking the input before running the calculation
write('to_check.dat', sheet, format='onetep-in', keywords = keywords)

sheet.calc = calc
# Will actually run the geometry optimisation
# using ONETEP internal BFGS
sheet.get_potential_energy()

Quickly restart with solvation effect using the soft sphere model

from ase.io import read
from ase.io.onetep import get_onetep_keywords

# Read from the previous run...
optimized_sheet = read("N_doped_graphene_001.out")

# Function to retrieve keywords dict from input file...
keywords = get_onetep_keywords('N_doped_graphene_001.dat')

# We add solvation keywords
keywords.update(
    {
    'is_implicit_solvent': True,
    'is_include_apolar': True,
    'is_smeared_ion_rep': True,
    'is_dielectric_model': 'fix_cavity',
    'is_dielectric_function' : 'soft_sphere'
    }
)

optimized_sheet.calc = Onetep(...)

...