Tutorial

This page gives an overview of using AtomsBase in practice and introduces the conventions followed across the AtomsBase ecosystem. It serves as a reference for both users interested in doing something with an AbstractSystem object as well as developers wishing to integrate their code with AtomsBase. See Interface for a precise specification of the AtomsBase interface.

For the examples we will mostly draw on the case of atomistic systems using the FlexibleSystem data structure. The Atom and FlexibleSystem data structure we focus on here should provide good defaults for many purposes.

High-level introduction

The main purpose of AtomsBase is to conveniently pass atomistic data between Julia packages. For example the following snippet loads an extxyz file using AtomsIO and returns it as an AtomsBase-compatible system (in data):

using AtomsIO
data = load_system("Si.extxyz")

Next we use ASEconvert to convert this system to python, such that we can make use of the atomistic simulation environment (ASE) to form a (2, 1, 1) supercell, which is afterwards converted back to Julia (by forming another AtomsBase-compatible system).

using ASEconvert
supercell = pyconvert(AbstractSystem, data * pytuple((2, 1, 1)))

Finally the supercell is passed to DFTK, where we attach pseudopotentials and run a PBE calculation:

using DFTK
cell_with_pseudos = attach_psp(supercell, Si="hgh/pbe/si-q4")
model = model_PBE(cell_with_pseudos)
basis = PlaneWaveBasis(model, Ecut=30, kgrid=(5, 5, 5)
self_consistent_field(basis).energy.total

For more high-level examples see also:

Atom interface and conventions

An Atom object can be constructed just by passing an identifier (e.g. symbol like :C, atomic number like 6) and a vector of positions as

atom = Atom(:C, [0, 1, 2.]u"bohr")

Atom(C, Z = 6, m = 12.011 u):
    position          : [0,1,2]u"a₀"
    species           : C

This automatically fills the atom with standard data such as the atomic mass. See also the reference of the Atom function for more ways to construct an atom.

Such data can be accessed using the AtomsBase interface functions, such as:

atomic_mass(atom)

12.011 u

atomic_symbol(atom)

atomic_number(atom)

position(atom)

3-element StaticArraysCore.SVector{3, Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(a₀,), 𝐋, nothing}}} with indices SOneTo(3):
 0.0 a₀
 1.0 a₀
 2.0 a₀

velocity(atom)

3-element StaticArraysCore.SVector{3, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(nm, s^-1), 𝐋 𝐓^-1, nothing}}} with indices SOneTo(3):
 0.0 nm s^-1
 0.0 nm s^-1
 0.0 nm s^-1

See atom.jl and the respective API documentation for details. Notice in particular that atomic_number specifies the element whereas species and also atomic_symbol may be more specific and may e.g. specify isotopes such as deuterium

deuterium = Atom(:D, [0, 1, 2.]u"bohr")
atomic_number(deuterium) == 1
atomic_symbol(deuterium) == :D

true

An equivalent dict-like interface based on keys, haskey, get and pairs is also available. For example

keys(atom)

(:position, :velocity, :species, :mass)

atom[:species]

pairs(atom)

Base.Generator{NTuple{4, Symbol}, AtomsBase.var"#60#61"{Atom{3, Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(a₀,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(nm, s^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐌, Unitful.FreeUnits{(u,), 𝐌, nothing}}}}}(AtomsBase.var"#60#61"{Atom{3, Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(a₀,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(nm, s^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐌, Unitful.FreeUnits{(u,), 𝐌, nothing}}}}(Atom(C,  [       0,        1,        2]u"a₀")), (:position, :velocity, :species, :mass))

This interface seamlessly generalises to working with user-specific atomic properties as will be discussed next.

Optional atomic properties

Custom properties can be attached to an Atom by supplying arbitrary keyword arguments upon construction. For example to attach a pseudopotential for using the structure with DFTK, construct the atom as

atom = Atom(:C, [0, 1, 2.]u"bohr", pseudopotential="hgh/lda/c-q4")

Atom(C, Z = 6, m = 12.011 u):
    position          : [0,1,2]u"a₀"
    species           : C
    pseudopotential   : hgh/lda/c-q4

which will make the pseudopotential identifier available as

atom[:pseudopotential]

"hgh/lda/c-q4"

Notice that such custom properties are fully integrated with the standard atomic properties, e.g. automatically available from the keys, haskey and pairs functions, e.g.:

@show haskey(atom, :pseudopotential)
pairs(atom)

Base.Generator{NTuple{5, Symbol}, AtomsBase.var"#60#61"{Atom{3, Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(a₀,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(nm, s^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐌, Unitful.FreeUnits{(u,), 𝐌, nothing}}}}}(AtomsBase.var"#60#61"{Atom{3, Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(a₀,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(nm, s^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐌, Unitful.FreeUnits{(u,), 𝐌, nothing}}}}(Atom(C,  [       0,        1,        2]u"a₀")), (:position, :velocity, :species, :mass, :pseudopotential))

Updating an atomic property proceeds similarly. E.g.

newatom = Atom(atom; atomic_mass=13u"u")

Atom(C, Z = 6, m = 12.011 u):
    position          : [0,1,2]u"a₀"
    species           : C
    atomic_mass       : 13 u
    pseudopotential   : hgh/lda/c-q4

makes a new carbon atom with all properties identical to atom (including custom ones), but setting the atomic_mass to 13 units.

To simplify interoperability some optional properties are reserved. For these:

Throughout the AtomsBase ecosystem these property keys carry a well-defined meaning. I.e. if they are supported by a code, the code expects them to hold the data specified below.
If a consuming code cannot make use of these properties, it should at least warn the user about it. For example if a library or simulation code does not support such a feature and drops the respective information it should @warn or (even better) interrupt execution with an error.

Property name	Unit / Type	Description
`:charge`	`Charge`	Net charge of the atom
`:covalent_radius`	`Length`	Covalent radius tabulated for the atom
`:vdw_radius`	`Length`	VdW radius tabulated for the atom
`:magnetic_moments`	`Union{Float64,Vector{Float64}}`	Initial magnetic moment
`:pseudopotential`	`String`	Pseudopotential or PAW keyword or `""` if Coulomb potential employed

A convenient way to iterate over all data stored in an atom offers the pairs function:

for (k, v) in pairs(atom)
    println("$k  =  $v")
end

position  =  Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(a₀,), 𝐋, nothing}}[0.0 a₀, 1.0 a₀, 2.0 a₀]
velocity  =  Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(nm, s^-1), 𝐋 𝐓^-1, nothing}}[0.0 nm s^-1, 0.0 nm s^-1, 0.0 nm s^-1]
species  =  C
mass  =  12.011 u
pseudopotential  =  hgh/lda/c-q4

System interface and conventions

Once the atoms are constructed these can be assembled into a system. For example to place a hydrogen molecule into a cubic box of 10Å and periodic boundary conditions, use:

box = ([10.0, 0.0, 0.0]u"Å", [0.0, 10.0, 0.0]u"Å", [0.0, 0.0, 10.0]u"Å")
hydrogen = periodic_system([Atom(:H, [0, 0, 1.]u"Å"),
                            Atom(:H, [0, 0, 3.]u"Å")],
                            box)

FlexibleSystem(H₂, pbc = TTT):
    bounding_box      : [      10        0        0;
                                0       10        0;
                                0        0       10]u"Å"

    Atom(H,  [       0,        0,        1]u"Å")
    Atom(H,  [       0,        0,        3]u"Å")

An update constructor for systems is supported as well (see AbstractSystem). For example

using AtomsBase: AbstractSystem
AbstractSystem(hydrogen;
               bounding_box=([5.0, 0.0, 0.0]u"Å",
                             [0.0, 5.0, 0.0]u"Å",
                             [0.0, 0.0, 5.0]u"Å"))

FlexibleSystem(H₂, pbc = TTT):
    bounding_box      : [       5        0        0;
                                0        5        0;
                                0        0        5]u"Å"

    Atom(H,  [       0,        0,        1]u"Å")
    Atom(H,  [       0,        0,        3]u"Å")

To update the atomic composition of the system, this function supports an atoms (or particles) keyword argument to supply the new set of atoms to be contained in the system.

Note that in this example FlexibleSystem( ... ) would have worked as well (since we are updating a FlexibleSystem). However, using the AbstractSystem constructor to update the system is more general as it allows for type-specific dispatching when updating other data structures implementing the AbstractSystem interface.

Similar to the atoms, system objects similarly support a functional-style access to system properties as well as a dict-style access:

bounding_box(hydrogen)

(Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(Å,), 𝐋, nothing}}[10.0 Å, 0.0 Å, 0.0 Å], Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(Å,), 𝐋, nothing}}[0.0 Å, 10.0 Å, 0.0 Å], Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(Å,), 𝐋, nothing}}[0.0 Å, 0.0 Å, 10.0 Å])

hydrogen[:periodicity]

(true, true, true)

pairs(hydrogen)

Base.Generator{Tuple{Symbol, Symbol}, AtomsBase.var"#7#8"{FlexibleSystem{3, Atom{3, Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(Å,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(Å, fs^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐌, Unitful.FreeUnits{(u,), 𝐌, nothing}}}, PeriodicCell{3, Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(Å,), 𝐋, nothing}}}}}}(AtomsBase.var"#7#8"{FlexibleSystem{3, Atom{3, Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(Å,), 𝐋, nothing}}, Unitful.Quantity{Float64, 𝐋 𝐓^-1, Unitful.FreeUnits{(Å, fs^-1), 𝐋 𝐓^-1, nothing}}, Unitful.Quantity{Float64, 𝐌, Unitful.FreeUnits{(u,), 𝐌, nothing}}}, PeriodicCell{3, Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(Å,), 𝐋, nothing}}}}}(FlexibleSystem(H₂, pbc = TTT)), (:bounding_box, :periodicity))

Moreover atomic properties of a specific atom or all atoms can be directly queried using the indexing notation:

hydrogen[1, :position]  # Position of first atom

3-element StaticArraysCore.SVector{3, Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(Å,), 𝐋, nothing}}} with indices SOneTo(3):
 0.0 Å
 0.0 Å
 1.0 Å

hydrogen[:, :position]  # All atomic symbols

2-element Vector{StaticArraysCore.SVector{3, Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(Å,), 𝐋, nothing}}}}:
 [0.0 Å, 0.0 Å, 1.0 Å]
 [0.0 Å, 0.0 Å, 3.0 Å]

Finally, supported keys of atomic properties can be directly queried at the system level using atomkeys and hasatomkey. Note that these functions only apply to atomic properties which are supported by all atoms of a system. In other words if a custom atomic property is only set in a few of the contained atoms, these functions will not consider it.

atomkeys(hydrogen)

(:position, :velocity, :species, :mass)

For constructing atomic systems the functions atomic_system, isolated_system, periodic_system are oftentimes more convenient as they provide specialisations for some standard atomic system setups. For example to setup a hydrogen system with periodic BCs, we can issue

bounding_box = ([10.0, 0.0, 0.0]u"Å", [0.0, 10.0, 0.0]u"Å", [0.0, 0.0, 10.0]u"Å")
hydrogen = periodic_system([:H => [0, 0, 1.]u"Å",
                            :H => [0, 0, 3.]u"Å"],
                           bounding_box)

FlexibleSystem(H₂, pbc = TTT):
    bounding_box      : [      10        0        0;
                                0       10        0;
                                0        0       10]u"Å"

    Atom(H,  [       0,        0,        1]u"Å")
    Atom(H,  [       0,        0,        3]u"Å")

To setup a silicon unit cell we can use fractional coordinates (which is common for solid-state simulations):

bounding_box = 10.26 / 2 .*  ([0, 0, 1]u"bohr", [1, 0, 1]u"bohr", [1, 1, 0]u"bohr")
silicon = periodic_system([:Si =>  ones(3)/8,
                           :Si => -ones(3)/8],
                           bounding_box, fractional=true)

FlexibleSystem(Si₂, pbc = TTT):
    bounding_box      : [       0        0     5.13;
                             5.13        0     5.13;
                             5.13     5.13        0]u"a₀"

    Atom(Si, [  1.2825,  0.64125,   1.2825]u"a₀")
    Atom(Si, [ -1.2825, -0.64125,  -1.2825]u"a₀")

Alternatively we can also place an isolated H2 molecule in vacuum, which is the standard setup for molecular simulations:

hydrogen = isolated_system([:H => [0, 0, 1.]u"bohr",
                            :H => [0, 0, 3.]u"bohr"])

FlexibleSystem(H₂, pbc = FFF):
    Atom(H,  [       0,        0,        1]u"a₀")
    Atom(H,  [       0,        0,        3]u"a₀")

Optional system properties

Similar to atoms, FlexibleSystem and FastSystem also support storing arbitrary data, for example

system = isolated_system([:H => [0, 0, 1.]u"bohr", :H => [0, 0, 3.]u"bohr"]; extra_data=42)

FlexibleSystem(H₂, pbc = FFF):
    extra_data        : 42

    Atom(H,  [       0,        0,        1]u"a₀")
    Atom(H,  [       0,        0,        3]u"a₀")

Again these custom properties are fully integrated with keys, haskey, pairs and get.

@show keys(system)

(:bounding_box, :periodicity, :extra_data)