Atoms¶

Adding and Modifying Atoms¶

The following code provides examples for adding atoms to a ChemicalSystem and modifying their attributes.

from scm.libbase import ChemicalSystem
from scm.libbase import Atom

mol = ChemicalSystem()

# Three different ways of adding atoms
mol.add_atom("O", coords=[0, 0, 0])
mol.add_atom(Atom("H"), coords=[2, 0, 0])
mol.add_atom(1, coords=[0, 2, 0], unit="Angstrom")

# Change the coordinates of the Oxygen atom (units: Angstrom)
mol.atoms[0].coords = [1, 2, 3]

# This is equivalent to the line above:
mol.coords[0, :] = [1, 2, 3]

# Change the first hydrogen atom to deuterium by adjusting the mass (units: Dalton):
mol.atoms[1].mass = 2.014

Methods to interact with the atoms in your ChemicalSystem:

property ChemicalSystem.atoms: AtomList

A list of all Atom instances that are part of this ChemicalSystem.

property ChemicalSystem.coords: AngstromCoordsArray

The coordinates of the atoms in angstrom.

ChemicalSystem.contains_atom(atom: Atom) → bool

Checks if an Atom instance is part of a ChemicalSystem.

ChemicalSystem.atom_index(atom: Atom) → int

Given an Atom instance, returns its index in ChemicalSystem.atoms.

ChemicalSystem.add_atom(atom: Atom) → None

ChemicalSystem.add_atom(atom: Atom, coords: ArrayLike, unit: str = 'angstrom') → None

ChemicalSystem.add_atom(Z: int, coords: ArrayLike, unit: str = 'angstrom') → None

ChemicalSystem.add_atom(element: Element, coords: ArrayLike, unit: str = 'angstrom') → None

ChemicalSystem.add_atom(symbol: Atom.T_AtomSymbol, coords: ArrayLike, unit: str = 'angstrom') → None

Add a new Atom to the ChemicalSystem.

ChemicalSystem.remove_atom(atom_index: int) → None

Removes a single atom from the system, given its index.

ChemicalSystem.remove_atoms(atom_indices: ArrayLike) → None

Removes multiple atoms from the system, given their atom indices.

ChemicalSystem.set_atom(atom_index: int, atom: Atom) → None

Safe assignment of an Atom instance to atoms[atom_index].

Is essentially the same as self.atoms[atom_index] = atom, but makes sure that all the atomic attributes carried by the Atom instance are enabled on the target ChemicalSystem: the direct assignment will throw a ChemicalSystemError if the needed atomic attributes are not enabled, while calling this method will enable them automatically and will not throw.

If Atom is part of a ChemicalSystem, its will overwrite whatever coordinates atoms[atom_index] had before. If Atom is not part of a ChemicalSystem, the coordinates of atoms[atom_index] will remain unchanged.

Changing the order of atoms¶

ChemicalSystem.reorder_atoms(atom_order: List[int]) → None

Reorders the system’s atoms given a list of atom indices in the desired order.

The atom_order shall be a list of integers. The list shall be of length num_atoms and its values a permutation of range(0,num_atoms). It shall contain the indices of the atoms in the original system in the order in which the atoms should occur after reordering.

Example: for a 6 atom system the atom_order = [4, 5, 0, 1, 2, 3] would have the effect of moving the last two atoms to the front of the system. (One might also look at atom_order as a mapping from the atom index in the reordered system to the atom index in the original system. Using the above example atom_order[1] == 5, which is the index that the new atom at index 1 had in the original system.)

ChemicalSystem.sorted_atom_order(comp: atom_comparator_func_t) → List[int]

Returns the atom order based on a user defined comparator function as a list of atom indices.

The returned list can be passed into the reorder_atoms method to actually sort the atoms according to comp. See also documentation of the reorder_atoms and sort_atoms method.

ChemicalSystem.sort_atoms(comp: atom_comparator_func_t) → None

Sorts the system’s atoms according to a user defined comparator function.

Atoms are sorted in ascending order where the user defined comparator function implements the < operator.

Example for sorting atoms descending in their atomic number, i.e. heavy atoms first:

chemsys.sort_atoms(lambda A, B: A.Z > B.Z)

Note that the atom sorting is stable: two atoms for which comp(A,B) returns False will keep their relative order in the sorted system.

Atomic species¶

One often needs to group a system’s atoms into sets of indistinguishable atoms, essentially answering the question: which atoms are considered the same?

In AMS these groups of equivalent atoms are often called “atomic species”. Note that the assignment of atoms into species is application dependent: an electronic structure code may only care about the elements of the atoms, while a molecular dynamics driver will also need to distinguish the different isotopes of the same element. Engines implementing classical forcefields will even need to distinguish atoms based on their forcefield atom types and formal charges.

The ChemicalSystem class allows users to provide a binary predicate function that implements the comparison of two atoms, to allow for any definition of atomic species:

ChemicalSystem.determine_species(comp: atom_comparator_func_t) → Tuple[int, List[int], List[int]]

Determines the present atomic species, based on a user defined comparator function.

The comparator shall return True for two atoms if they are to be considered of the same species.

Returns a 3-tuple consisting of:

1. num_species is the total number of detected distinct species.

2. species_prototypes a list of size num_species with the index of the first atom in the system that was recognized as being of a distinct species. This atom is considered the species’ prototype atom and all atoms of the same species compare equal to it by the user defined comparator function comp.

3. atom_to_species a list of size num_atoms mapping each atom to its species.

See also the related get_species function, which returns the same data as a list of atom indices for each species.

In practice is usually easiest to implement the comparison as a lambda function right at the call site:

# Distinguish isotopes as separate species:
... = chemsys.determine_species(lambda A, B: A.Z == B.Z and A.mass == B.mass)

Note that you can access the ChemicalSystem instance from within the lambda, allowing access to things like regions or bonds, that are not part of the Atom class itself.

# Elements as species, but separately for region "substrate" and the rest:
... = chemsys.determine_species(
     lambda A, B: A.Z == B.Z and
     chemsys.is_atom_in_region(A, "substrate") == chemsys.is_atom_in_region(B, "substrate")
)

# Distinguish species based on the number of bonds
# e.g. all 4-bonded atoms are equivalent:
... = chemsys.determine_species(
lambda A, B: chemsys.bonds.num_bonds(chemsys.atom_index(A)) == chemsys.bonds.num_bonds(chemsys.atom_index(B))
)

If you primarily want to iterate over all atoms within a species, the get_species function may return the same information in a more convenient format:

ChemicalSystem.get_species(comp: atom_comparator_func_t) → List[List[int]]

Determines the present atomic species, based on a user defined comparator function.

The comparator shall return True for two atoms if they are to be considered of the same species.

Returns a list of atom indices for each species, allowing for easy iteration over all atoms belonging to a species. Atom indices within each species are in ascending order.

for species in chemsys.get_species():
    for iat in species:
        # do something with chemsys.atoms[iat]

See also the related determine_species function, which returns the same information in a different format.

Atom Attributes¶

Atom attributes provide additional, customizable information for each atom. They are organized into categories, each relevant to different computation engines or modules.

• Generic Attributes: These are common to all atoms, like mass.

• Engine/Module-Specific Attributes: Specific to certain computation engines or modules. For example, attributes relevant only to the ADF engine would be in the adf_attribs group and prefixed with ‘adf’ in the System Block.

You can toggle these attribute groups on or off. Here’s an example, where we initialize some atomic attributes in the System Block:

from scm.libbase import ChemicalSystem

mol = ChemicalSystem(
"""
System
     Atoms
          O 0.0 0.0 0.0 forcefield.type=O_water
          H 1.0 0.0 0.0 forcefield.type=H_water mass=2.014
          H 0.0 1.0 0.0 forcefield.type=H_water
     End
End
"""
)

# The Generic Attributes are always available (and have default values if not specified)
print(mol.atoms[1].mass) # Outputs 2.014
print(mol.atoms[2].mass) # Outputs 1.00798

# The forcefield attributes are enabled because there were defined in the system block
print(mol.atom_attributes_enabled('forcefield')) # Outputs True

# Accessing the forcefield atomic attributes
print(mol.atoms[0].forcefield.type) # Outputs O_water

# The adf attributes are not enabled
print(mol.atom_attributes_enabled('adf')) # Outputs False

# Enable the 'adf' attributes group before using it:
mol.enable_atom_attributes("adf")
mol.atoms[0].adf.f = "fragment_1"

Note: The Atom class also contains the atom attributes.

Atom attributes in the ChemicalSystem¶

property ChemicalSystem.adf_attribs: AtomAttributesADFList | None

An optional list of all AtomAttributesADF for all atoms.

property ChemicalSystem.band_attribs: AtomAttributesBANDList | None

An optional list of all AtomAttributesBAND for all atoms.

property ChemicalSystem.qe_attribs: AtomAttributesQEList | None

An optional list of all AtomAttributesQE for all atoms.

property ChemicalSystem.dftb_attribs: AtomAttributesDFTBList | None

An optional list of all AtomAttributesDFTB for all atoms.

property ChemicalSystem.reaxff_attribs: AtomAttributesReaxFFList | None

An optional list of all AtomAttributesReaxFF for all atoms.

property ChemicalSystem.forcefield_attribs: AtomAttributesForcefieldList | None

An optional list of all AtomAttributesForcefield for all atoms.

property ChemicalSystem.gui_attribs: AtomAttributesGUIList | None

An optional list of all AtomAttributesGUI for all atoms.

ChemicalSystem.enable_atom_attributes(group_prefix: Literal['gui', 'adf', 'band', 'forcefield', 'dftb', 'reaxff', 'qe']) → None

Enables the use of a group of atomic attributes.

ChemicalSystem.disable_atom_attributes(group_prefix: Literal['gui', 'adf', 'band', 'forcefield', 'dftb', 'reaxff', 'qe']) → None

Disables the use of a group of atomic attributes. Any set attributes within the group will be discarded.

ChemicalSystem.atom_attributes_enabled(group_prefix: Literal['gui', 'adf', 'band', 'forcefield', 'dftb', 'reaxff', 'qe']) → bool

Checks if a group of atomic attributes is enabled or not.

ChemicalSystem.has_ghost_atoms() → bool

Checks whether the ChemicalSystem contains any Ghost atoms.

Ghost atoms are a special construct in AMS for calculating basis set superposition errors. Please refer to the ADF and BAND manuals for details.

Engine/Module-Specific Atom Attributes Classes¶

Atom class¶

The Atom class contains all the relevant data for an individual atom:

class Atom

class Atom(Z: int)

class Atom(element: Element)

class Atom(symbol: T_AtomSymbol)

Class representing a single atom and its properties.

copy() → Atom

Creates a deep copy of the Atom.

has_identical_attributes(other: Atom) → bool

has_identical_attributes(other: Atom, properties: List[T_AtomAttributeGroup]) → bool

Helper for @overload to raise when called.

in_chemicalsystem() → bool

Check if this Atom instance is part of a ChemicalSystem.

transmute(Z: int)

transmute(element: Element)

transmute(symbol: T_AtomSymbol)

Helper for @overload to raise when called.

property Z: int

The atomic number of the atom.

property adf: AtomAttributesADF | None

Atomic attributes used by the ADF engine.

property band: AtomAttributesBAND | None

Atomic attributes used by the AMS BAND engine.

property coords: AngstromCoordsArray

The coordinates of the atom in angstrom.

property dftb: AtomAttributesDFTB | None

Atomic attributes used by the AMS DFTB engine.

property element: Element

Returns the atom’s element.

property forcefield: AtomAttributesForcefield | None

Atomic attributes used by the AMS Forcefield engine.

property gui: AtomAttributesGUI | None

Atomic attributes used for visualization in the GUI.

property is_ghost: bool

Whether the atom is a Ghost atom.

Ghost atoms are a special construct in AMS for calculating basis set superposition errors. Please refer to the ADF and BAND manuals for details.

property mass: float

The mass of the atom in dalton.

property qe: AtomAttributesQE | None

Atomic attributes used by the AMS QuantumESPRESSO engine.

property reaxff: AtomAttributesReaxFF | None

Atomic attributes used by the AMS ReaxFF engine.

property symbol: T_AtomSymbol

The atom’s symbol (e.g. ‘C’ or ‘Au’). May be ‘Gh.C’ for ghost atoms.