The Pymatgen python library allows to setup
solid-state calculations using a flexible set of classes as well as an API
to an online data base of structures. Its Structure
and Lattice
objects are directly supported by the DFTK load_atoms
and load_lattice
functions, such that DFTK may be readily used to run calculation on systems
defined in pymatgen. Using the pymatgen_structure
function a conversion
from DFTK to pymatgen structures is also possible. In the following we
use this to create a silicon supercell and find its LDA ground state
using direct minimisation. To run this example Julia's PyCall
package
needs to be able to find an installation of pymatgen
.
First we setup the silicon lattice in DFTK.
using DFTK
a = 10.263141334305942 # Lattice constant in Bohr
lattice = a / 2 .* [[0 1 1.]; [1 0 1.]; [1 1 0.]]
Si = ElementPsp(:Si, psp=load_psp("hgh/lda/Si-q4"))
atoms = [Si => [ones(3)/8, -ones(3)/8]];
Next we make a [2, 2, 2]
supercell using pymatgen
pystruct = pymatgen_structure(lattice, atoms)
pystruct.make_supercell([2, 2, 2])
lattice = load_lattice(pystruct)
atoms = [Si => [s.frac_coords for s in pystruct.sites]];
┌ Warning: `vendor()` is deprecated, use `BLAS.get_config()` and inspect the output instead │ caller = npyinitialize() at numpy.jl:67 └ @ PyCall /home/runner/.julia/packages/PyCall/L0fLP/src/numpy.jl:67
Setup an LDA model and discretize using
a single k-point and a small Ecut
of 5 Hartree.
model = model_LDA(lattice, atoms)
basis = PlaneWaveBasis(model; Ecut=5, kgrid=(1, 1, 1))
PlaneWaveBasis discretization: Ecut : 5.0 Ha fft_size : (32, 32, 32) kgrid type : Monkhorst-Pack kgrid : [1, 1, 1] num. irred. kpoints : 1 Discretized Model(lda_xc_teter93, 3D): lattice (in Bohr) : [0 , 10.2631 , 10.2631 ] [10.2631 , 0 , 10.2631 ] [10.2631 , 10.2631 , 0 ] unit cell volume : 2162.1 Bohr³ atoms : Si₁₆ atom potentials : ElementPsp(Si, psp=hgh/lda/si-q4) num. electrons : 64 spin polarization : none temperature : 0 Ha terms : Kinetic() AtomicLocal() AtomicNonlocal() Ewald() PspCorrection() Hartree() Xc(:lda_xc_teter93)
Find the ground state using direct minimisation (always using SCF is boring ...)
scfres = direct_minimization(basis, tol=1e-5);
Iter Function value Gradient norm 0 1.116062e+02 1.517882e+00 * time: 0.8374879360198975 1 1.074000e+01 8.582361e-01 * time: 2.374027967453003 2 -1.134891e+01 9.564115e-01 * time: 3.0860068798065186 3 -3.392864e+01 7.674659e-01 * time: 4.179577827453613 4 -4.744996e+01 5.731013e-01 * time: 5.210864782333374 5 -5.688425e+01 2.348844e-01 * time: 6.2739949226379395 6 -5.983825e+01 1.758409e-01 * time: 6.982828855514526 7 -6.094845e+01 5.582174e-02 * time: 7.708448886871338 8 -6.133676e+01 9.125249e-02 * time: 8.413647890090942 9 -6.160948e+01 3.923540e-02 * time: 9.125579833984375 10 -6.180592e+01 2.782618e-02 * time: 9.793018817901611 11 -6.195302e+01 2.462728e-02 * time: 10.50530195236206 12 -6.201198e+01 1.780274e-02 * time: 11.215266942977905 13 -6.207398e+01 1.724701e-02 * time: 11.910501956939697 14 -6.210939e+01 1.915338e-02 * time: 12.623705863952637 15 -6.215032e+01 1.554781e-02 * time: 13.3208909034729 16 -6.218157e+01 1.261166e-02 * time: 14.052361965179443 17 -6.220570e+01 1.017819e-02 * time: 14.727994918823242 18 -6.222428e+01 7.246597e-03 * time: 15.435925006866455 19 -6.223493e+01 7.214679e-03 * time: 16.158989906311035 20 -6.224150e+01 6.314786e-03 * time: 16.858760833740234 21 -6.224625e+01 4.292858e-03 * time: 17.566049814224243 22 -6.225074e+01 4.058197e-03 * time: 18.272429943084717 23 -6.225429e+01 3.602924e-03 * time: 18.96968388557434 24 -6.225701e+01 2.957667e-03 * time: 19.662516832351685 25 -6.225901e+01 2.542853e-03 * time: 20.36894178390503 26 -6.226025e+01 1.819695e-03 * time: 21.061124801635742 27 -6.226095e+01 1.054844e-03 * time: 21.76589798927307 28 -6.226130e+01 1.066389e-03 * time: 22.4663348197937 29 -6.226145e+01 6.920819e-04 * time: 23.182384967803955 30 -6.226153e+01 5.939634e-04 * time: 23.892266988754272 31 -6.226158e+01 4.624261e-04 * time: 24.60138487815857 32 -6.226161e+01 4.021984e-04 * time: 25.302117824554443 33 -6.226163e+01 3.950532e-04 * time: 25.98001480102539 34 -6.226164e+01 2.510007e-04 * time: 26.66957998275757 35 -6.226165e+01 2.122974e-04 * time: 27.343291997909546 36 -6.226166e+01 1.651182e-04 * time: 28.031521797180176 37 -6.226166e+01 9.719340e-05 * time: 28.719368934631348 38 -6.226166e+01 8.177096e-05 * time: 29.411049842834473 39 -6.226166e+01 7.317111e-05 * time: 30.11129593849182 40 -6.226167e+01 4.129686e-05 * time: 30.806140899658203 41 -6.226167e+01 3.298540e-05 * time: 31.497408866882324 42 -6.226167e+01 3.020884e-05 * time: 32.1739239692688 43 -6.226167e+01 1.806679e-05 * time: 32.90002989768982 44 -6.226167e+01 1.691626e-05 * time: 33.58825898170471 45 -6.226167e+01 1.223037e-05 * time: 34.267826795578 46 -6.226167e+01 1.068305e-05 * time: 34.949291944503784 47 -6.226167e+01 1.002592e-05 * time: 35.65868401527405 48 -6.226167e+01 7.064774e-06 * time: 36.35013699531555 49 -6.226167e+01 5.377785e-06 * time: 37.05099892616272 50 -6.226167e+01 3.751058e-06 * time: 37.77816987037659 51 -6.226167e+01 2.727381e-06 * time: 38.58723497390747
scfres.energies
Energy breakdown (in Ha): Kinetic 25.7671060 AtomicLocal -18.8557714 AtomicNonlocal 14.8522689 Ewald -67.1831486 PspCorrection -2.3569765 Hartree 4.8485372 Xc -19.3336820 total -62.261666460506