AtomsBase integration
AtomsBase.jl is a common interface for representing atomic structures in Julia. DFTK directly supports using such structures to run a calculation as is demonstrated here.
using DFTK
Feeding an AtomsBase AbstractSystem to DFTK
In this example we construct a silicon system using the ase.build.bulk
routine from the atomistic simulation environment (ASE), which is exposed by ASEconvert as an AtomsBase AbstractSystem
.
# Construct bulk system and convert to an AbstractSystem
using ASEconvert
system_ase = ase.build.bulk("Si")
system = pyconvert(AbstractSystem, system_ase)
FlexibleSystem(Si₂, periodic = TTT):
bounding_box : [ 0 2.715 2.715;
2.715 0 2.715;
2.715 2.715 0]u"Å"
Atom(Si, [ 0, 0, 0]u"Å")
Atom(Si, [ 1.3575, 1.3575, 1.3575]u"Å")
Si
Si
To use an AbstractSystem in DFTK, we attach pseudopotentials, construct a DFT model, discretise and solve:
system = attach_psp(system; Si="hgh/lda/si-q4")
model = model_LDA(system; temperature=1e-3)
basis = PlaneWaveBasis(model; Ecut=15, kgrid=[4, 4, 4])
scfres = self_consistent_field(basis, tol=1e-8);
n Energy log10(ΔE) log10(Δρ) Diag Δtime
--- --------------- --------- --------- ---- ------
1 -7.921684723798 -0.69 5.6
2 -7.926167418854 -2.35 -1.22 1.0 215ms
3 -7.926837141798 -3.17 -2.37 1.8 237ms
4 -7.926861512551 -4.61 -3.02 2.8 256ms
5 -7.926861636120 -6.91 -3.36 2.0 250ms
6 -7.926861666509 -7.52 -3.72 1.8 218ms
7 -7.926861680518 -7.85 -4.37 1.1 202ms
8 -7.926861681793 -8.89 -4.89 2.1 233ms
9 -7.926861681862 -10.16 -5.29 1.8 241ms
10 -7.926861681872 -11.02 -5.96 2.0 228ms
11 -7.926861681873 -12.04 -6.77 2.0 234ms
12 -7.926861681873 -13.39 -7.42 2.4 264ms
13 -7.926861681873 + -14.45 -8.22 2.2 234ms
If we did not want to use ASE we could of course use any other package which yields an AbstractSystem object. This includes:
Reading a system using AtomsIO
using AtomsIO
# Read a file using [AtomsIO](https://github.com/mfherbst/AtomsIO.jl),
# which directly yields an AbstractSystem.
system = load_system("Si.extxyz")
# Now run the LDA calculation:
system = attach_psp(system; Si="hgh/lda/si-q4")
model = model_LDA(system; temperature=1e-3)
basis = PlaneWaveBasis(model; Ecut=15, kgrid=[4, 4, 4])
scfres = self_consistent_field(basis, tol=1e-8);
n Energy log10(ΔE) log10(Δρ) Diag Δtime
--- --------------- --------- --------- ---- ------
1 -7.921698101427 -0.69 5.8
2 -7.926165146536 -2.35 -1.22 1.0 215ms
3 -7.926837719298 -3.17 -2.37 1.8 238ms
4 -7.926861518543 -4.62 -3.03 2.8 258ms
5 -7.926861643788 -6.90 -3.41 1.6 255ms
6 -7.926861670183 -7.58 -3.81 1.9 222ms
7 -7.926861680419 -7.99 -4.33 1.6 215ms
8 -7.926861681776 -8.87 -4.87 1.9 224ms
9 -7.926861681864 -10.06 -5.35 1.6 245ms
10 -7.926861681871 -11.11 -5.89 1.6 219ms
11 -7.926861681873 -11.95 -7.09 1.9 228ms
12 -7.926861681873 -13.57 -7.61 3.6 300ms
13 -7.926861681873 -15.05 -8.46 2.0 230ms
The same could be achieved using ExtXYZ by system = Atoms(read_frame("Si.extxyz"))
, since the ExtXYZ.Atoms
object is directly AtomsBase-compatible.
Directly setting up a system in AtomsBase
using AtomsBase
using Unitful
using UnitfulAtomic
# Construct a system in the AtomsBase world
a = 10.26u"bohr" # Silicon lattice constant
lattice = a / 2 * [[0, 1, 1.], # Lattice as vector of vectors
[1, 0, 1.],
[1, 1, 0.]]
atoms = [:Si => ones(3)/8, :Si => -ones(3)/8]
system = periodic_system(atoms, lattice; fractional=true)
# Now run the LDA calculation:
system = attach_psp(system; Si="hgh/lda/si-q4")
model = model_LDA(system; temperature=1e-3)
basis = PlaneWaveBasis(model; Ecut=15, kgrid=[4, 4, 4])
scfres = self_consistent_field(basis, tol=1e-4);
n Energy log10(ΔE) log10(Δρ) Diag Δtime
--- --------------- --------- --------- ---- ------
1 -7.921691762045 -0.69 5.9
2 -7.926170013481 -2.35 -1.22 1.0 203ms
3 -7.926840264288 -3.17 -2.37 1.9 236ms
4 -7.926864880882 -4.61 -3.01 2.5 253ms
5 -7.926865037564 -6.80 -3.33 1.9 222ms
6 -7.926865076638 -7.41 -3.71 1.5 223ms
7 -7.926865091700 -7.82 -4.38 1.2 195ms
Obtaining an AbstractSystem from DFTK data
At any point we can also get back the DFTK model as an AtomsBase-compatible AbstractSystem
:
second_system = atomic_system(model)
FlexibleSystem(Si₂, periodic = TTT):
bounding_box : [ 0 5.13 5.13;
5.13 0 5.13;
5.13 5.13 0]u"a₀"
Atom(Si, [ 1.2825, 1.2825, 1.2825]u"a₀")
Atom(Si, [ -1.2825, -1.2825, -1.2825]u"a₀")
Si
Si
Similarly DFTK offers a method to the atomic_system
and periodic_system
functions (from AtomsBase), which enable a seamless conversion of the usual data structures for setting up DFTK calculations into an AbstractSystem
:
lattice = 5.431u"Å" / 2 * [[0 1 1.];
[1 0 1.];
[1 1 0.]];
Si = ElementPsp(:Si, psp=load_psp("hgh/lda/Si-q4"))
atoms = [Si, Si]
positions = [ones(3)/8, -ones(3)/8]
third_system = atomic_system(lattice, atoms, positions)
FlexibleSystem(Si₂, periodic = TTT):
bounding_box : [ 0 5.13155 5.13155;
5.13155 0 5.13155;
5.13155 5.13155 0]u"a₀"
Atom(Si, [ 1.28289, 1.28289, 1.28289]u"a₀")
Atom(Si, [-1.28289, -1.28289, -1.28289]u"a₀")
Si
Si