Elk Code FAQ
0. General questions
0.1 How can I obtain the code?
Download the latest version from
0.2 How can I be informed of new releases?
You can join the
for notification of new releases, features and bugs in the code.
0.3 Does the name Elk mean anything?
It could stand for "El
ectrons in k
-space", or it could just refer
to the animal Alces alces
which is also called a moose in North America.
0.4 What are the units used by the code?
Unless otherwise specified, all the units are
. In particular,
the electron rest mass, me
; the electronic charge,
; Planck's constant over two pi, h-bar; and the electrostatic constant
all equal one. The unit of energy is the Hartree
(Ha) and the unit of length is the Bohr. Rydbergs, electron volts (eV) and
Ångströms are never used.
0.5 How does Elk compare with other FP-LAPW codes?
Elk has similar capabilities as other codes, but is specifically designed
to make the implementation of new developments in DFT as easy as possible. To
this end, the code is cleanly written, fully documented and open source.
0.6 How is the documentation of the code handled?
Every subroutine and function has a LaTeX header which can be read by the
and compiled into a single document. Use the command make doc
directory to produce the file elk.pdf
0.7 Is there a graphical user interface (GUI) available?
No. The code requires a single user input file, elk.in
, which for
most cases should be easier to understand than a GUI.
0.8 If I modify the Elk code do I have to release my modifications?
No. The GPL
does not oblige
you to release your modified version, it simply requires that if you do, it
should be under the GPL.
0.9 Is the code suitable for production use?
Yes, for many purposes. It is always good practice to check the consistency of
the results with other codes. Features which are marked as
should not be used for production.
0.10 How can I cite the code in an article?
Simply cite the webpage http://elk.sourceforge.net/
. Citation of Elk
is encouraged but not mandatory.
1. Compiling the code
1.1 How do I compile the code?
Run the setup
script in the elk
directory. You can then
edit the file make.inc
and set the compiler options for maximum
performance on your system. Finally run the command make all
. This will
compile the main code as well as several auxiliary programs. The executable
is located in the src
directory and can be copied to
a location in the executable search path.
1.2 What do the auxilliary programs do?
The auxiliary programs include Spacegroup for producing crystal geometries from
spacegroup data and EOS for fitting equations of state to energy-volume data.
They can be found in their respective directories.
1.3 Why does the code not compile on my computer?
Elk is written in strict Fortran 90 syntax and should compile without
warnings on all F90 compilers. Some compilers may have difficulty with the BLAS
and LAPACK libraries. In this case you will have to refer to your compiler
manual and adjust the F77 options in the file make.inc
. If the main
part of the code fails to compile then please contact one of the authors.
1.4 Which version of BLAS/LAPACK does the code require?
Version 3. We strongly recommend the use of machine-optimised libraries as
Elk uses them heavily. If these are not available then
can be used instead.
1.5 Which fast Fourier transform (FFT) package does Elk use?
It uses a modified version of
is a highly efficient implementation of an in-place N-dimensional complex FFT.
However you can easily replace this with your own machine-optimised FFT package.
See the routine zfftifc
2. Running the code
2.1 How do I restart a calculation?
=1 and restart. You can change almost any parameter in the
input files except the number of atoms.
2.2 How can a running calculation be stopped cleanly?
Create an empty file STOP
in the running directory (use
). The code will then complete the current self-consistent
loop, write the file STATE.OUT
2.3 How can I force a running calculation to write STATE.OUT?
Create an empty file WRITE
in the running directory (use
). You will know STATE.OUT
has been written
2.4 How can I run Elk on a parallel computer?
Elk uses OpenMP
to run in parallel and
your Fortran compiler must support this. Simply set the appropriate compiler
options in make.inc
and the environment variable OMP_NUM_THREADS
to the number of threads you require. OpenMP cannot be used on a distributed
cluster at this stage (although there is a new Intel compiler which supports
this). Note that the parallelisation is coarse-grained and performed mostly over
-point set. Parallel BLAS, LAPACK and FFT libraries can also be used
for fine-grained parallelism. Phonon calculations can be spread natively over a
distributed system, so long as each node sees the same directory.
2.5 Why do my OpenMP runs keep crashing?
We think this is because local variables for each thread are allocated on the
stack rather than the heap, and some compilers don't allocate enough stack
space. Set the compiler options in make.inc
to increase the available
stack space. A second possibility is that your LAPACK library is not
thread-safe: try compiling the code using the native BLAS/LAPACK libraries.
2.6 How can I specify a path for use as scratch space?
Use the variable scrpath
. The eigenvalue, eigenvector and occupancy
files will be written to this directory. Note that you will need to retain
these files in order to perform further tasks.
2.7 Why do I keep getting the warning error
Warning(occupy): minimum eigenvalue less than minimum linearisation
If this warning goes away after the first few iterations then the results may
be trusted. If it persists then this indicates that
• The muffin-tins are highly mismatched i.e. there are both large and
small muffin-tins in the unit cell. This can be rectified by either setting
=.true. in elk.in
; or increasing the smallest and
decreasing the largest rmt
in the species files
is too small: increase it in small steps until the warning
is too large: decrease it in small steps until the
2.8 How can I make Elk run faster?
• Adjust the compiler options in make.inc
for optimal performance
• Use machine optimised versions of BLAS, LAPACK and the FFT library
• If you are running on a network file system, then set the scratch path
) to a local directory
• Reduce rgkmax
(but make sure the results are still converged)
• Run the code in parallel using OpenMP
2.9 How can I run the code without updating the density and potential in
Set maxscl to 1, this will run one self-consistent loop but wont write the
density and potential at the end of it. The eigenvalues and eigenvectors are
2.10 How can I make the code write STATE.OUT after every
Normally, the density and potentials are written to the file STATE.OUT
only after completion of the self-consistent loop. By setting nwrite
a positive integer the file will be written during the loop every
3. Crystal structure
3.1 The input file requires the lattice vectors and the atomic positions. How
can I set up the geometry using only the lattice parameters, spacegroup number
and Wyckoff coordinates?
This is not incorporated directly in the code. Instead, you can use the
auxiliary program Spacegroup to generate the vectors and positions in this way.
This will also produce an output file which can be viewed with
The Spacegroup utility also allows for the construction of supercells. See the
for more details.
3.2 I have used the conventional unit cell instead of the primitive cell. Is
Elk smart enough to use the smaller cell?
Not by default. You have to set the variable primcell
and the conventional cell will then be automatically reduced to
the primitive cell.
3.3 My muffin-tins overlap, what should I do?
You can either set the muffin-tin radii in the relevant species files, or
to be .true.
and the optimal radii will be
3.4 How can I fit an energy-volume equation of state?
Use the EOS utility provided.
3.5 How can I set up a molecule?
Set the input variable molecule
The atomic positions
in the atoms
block are then assumed to be in Cartesian coordinates. See
the Users' Guide and examples for more details.
4. Plotting output data
4.1 How can I plot the 1D output of the code?
The 1D plotting output is always in the form of a simple text-based list. We use
the open source program
for producing 1D
plots. The plot files can be imported into Grace without modification.
4.2 How can I plot the 2D/3D scalar or vector field output of the code?
The 2- and 3-dimensional output of the code always consists of a text-based list
of the coordinates and the scalar or vector value of the field. We have
sucessfully used the open source package
for all types of 2 and 3D plots.
4.3 Where is the Fermi level in the DOS, band structure and Fermi surface
Always at zero.
4.4 How can I plot a band structure for particular high-symmetry directions?
Firstly, you should define reciprocal space lattice coordinates for vertices of
the directions you are interested in. To do this consult the images of Brillouin
zones for various crystal structures (the files are located in the
directory). These vertices should be entered in the
block. After performing the calculation you should have the
which contain the band
structure and lines indicating the location of the vertices in the plot, respectively.
4.5 Can Elk produce a band structure plot with line thickness
proportional to the s, p, d or f character of a
Yes, use tasks
=21. Band structures are then plotted to files
for atom aaaa
of species ss
which include the band characters for each l
component of that atom in
columns 4 onwards. Column 3 contains the sum over l
of the characters.
4.6 Why is my energy vs. volume curve jagged?
This is most often due to poor convergence with respect to the number of
augmented plane waves and k
-points, see J. Phys. Conden. Matt. 2
4395. Increase rgkmax
, or both.
4.7 How do I plot a total or partial density of states (DOS)?
to 10 and run the code. This will produce a total and
partial spin-resolved density of states. Note that to improve the quality of
the DOS, you can first run one self-consistent loop (set maxscl
to 1) with a larger k
-point set. This will not update
the density in STATE.OUT
but will produce a new set of eigenvalues and
Fermi energy. Compute the DOS as usual afterwards. See also the dos
block in the Users' Guide.
4.8 What is the structure of a PDOS_Sss_Aaaaa.OUT file?
projections of the atomic site resolved DOS are
arranged in logical order in the PDOS
) = (0,0), (1,-1), (1,0), (1,1), (2,-2), (2,-1), (2,0), (2,1),
(2,2), etc., and separated by blank lines. For example, the d
DOS is found in blocks 5-9.
4.9 How can the partial DOS be resolved into irreducible representations (IR)
like eg and t2g?
The partial DOS file PDOS_Sss_Aaaaa.OUT
projections which correspond to the IR of the site symmetry for atom
of species ss
. To find out which projections belong to
the same IR, consult the file ELMIREP.OUT
. This file lists the
eigenvalues of a random matrix in the spherical harmonic basis which has been
symmetrised with the site symmetries. Degenerate eigenvalues indicate which
projections belonging to the same IR. This feature can be switched off by
4.10 Can I get a DOS summed over m quantum numbers?
Yes. By default, the partial density of states is resolved over
) quantum numbers. If you want the partial DOS to be summed
(and thus depend on l only) set dosmsum
to be .true.
4.11 I have performed a calculation for a system which definitely doesn't
have d, f,... electrons whereas there are blocks of DOS
corresponding to such electrons in the PDOS file. Is that normal?
Yes. This is because in a solid spherical symmetry is lost and the wave functions
have non-zero coefficients for all angular momenta.
4.12 Why does my plot not match well at the muffin-tin boundary?
You can improve the matching by increasing lmaxapw
, or the number of points in the muffin-tin
which is set by nrmt
in the species files. If you change any of these
values, you must reconverge your calculation before plotting. ELF plots are
also sensitive to lradstp
, which sets the step size for the coarse
radial mesh. Setting lradstp
=1 will improve the continuity of the ELF
4.13 How can I obtain orbital occupation numbers after performing a LDA+U
Look in the files DMATLU.OUT
. It will contain the density (occupancy)
matrix in the basis of spherical harmonics. Note: non-diagonal elements of the
occupation matrix can have any sign, diagonal elements are real and positive.
Total shell occupancy is the sum over all diagonal elements. For more details
5.1 How is magnetism implemented in the code?
HB = ge/(4c) σ·Bext,
Elk handles magnetism in a unique way thanks to an idea of Lars Nordström.
In the first-variational step, the Hamiltonian matrix is set up in the APW basis
using only the scalar effective potential, vs. The eigenvalues
and wavefunctions obtained from this are then used as a basis for the
second-variational Hamiltonian in which the effective magnetic field is included
as σ⋅Bs. This Hamiltonian is also diagonalised and
full spinor wavefunctions obtained which are then used to determine the density,
ρ, and magnetisation vector field, m.
5.2 What is the purpose of the external magnetic fields?
By default, the global field bfieldc and the individual muffin-tin
fields bfcmt are added to the second variational Hamiltonian only as a
means of breaking spin symmetry. These fields are considered infinitesimal
and their associated energy is not added to the total. Therefore they should be
made only as large as is required to break spin symmetry.
5.3 Supposing I want my external magnetic fields to be finite and have their
energy added to the total?
Just set bfinite to .true.
5.4 How do I include spin-orbit coupling?
Set spinorb to .true. The σ⋅L term is then
included in the second-variational step.
5.5 How do I perform an anti-ferromagnetic magnetic (AFM) calculation?
You need at least two atoms in your unit cell. Point the magnetic field on the
first atom in the +z direction and on the second atom in the -z
5.6 How do I perform a fixed spin moment (FSM) calculation?
Set fixspin to be .true. and specify the moment you want with
the variable momfix. Note that you must specify all three components of
momfix and in the same direction as the external magnetic field.
5.7 Can I perform a FSM calculation with spin-orbit coupling and
non-collinear magnetism switched on?
5.8 Does non-collinear magnetism work with GGA functionals?
Yes, although strictly speaking it should only be used with LSDA functionals
because a different spin-rotation is applied to the magnetisation at every
point in space.
5.9 How can I run a spin-spiral calculation
Set spinsprl to be .true. and select the particular
q-vector you want with vqlss in reciprocal lattice coordinates.
5.10 Why is the degeneracy slightly lifted of states I know are degenerate in
the spin-polarised case?
Assuming the states in question are supposed to be degenerate, the lifting is
due to the exchange-correlation effective field being coverted from spherical
harmonics to spherical coordinates. This always entails some small breaking of
symmetry. The effect can be made negligible by increasing lmaxvr.
5.11 Is it important to set the number of empty states to a large number
for magnetic calculations?
Yes! The variable nempty detemines the number of unoccupied
first-variational states. Since the second-variational step is done using the
first-variational states as the basis, it is crucial to have a sufficient
number of them. For difficult cases you may have set nempty to 30 or
5.12 What are the units of the external magnetic field used by the code?
The relevant term in the Hamiltonian is
where ge is the electron g-factor. In terms of SI units,
one unit of bfieldc is equal to 1715.255397557 Tesla.
6. Forces, phonons and structural optimisation
6.1 Do forces work with spin-orbit coupling and non-collinear magnetism?
6.2 How can I constrain the positions of particular atoms during a structural
Set the atomic mass to a negative value in the species files. The code considers
these to have infinite mass and does not move them.
6.3 How can I restart a structural optimisation run?
Replace the geometry information in elk.in with that in
GEOMETRY_OPT.OUT. You can simply append this to the file with
cat GEOMETRY_OPT.OUT >> elk.in
because the code will always use the last instance of any input variable.
Then set tasks=3 and rerun.
6.4 The incomplete basis set (IBS) corrections take a long time to
compute, is it possible to switch them off?
Yes, just set tfibs to .false.
6.5 I'm only interested in certain phonon modes, how can I avoid calculating
all the others?
As with the constrained structural optimisation, make the atomic mass in the
species file negative. The code will then assume these masses are infinite and
the corresponding rows and columns in the dynamical matrices are zero.
6.6 Is the acoustic sum rule enforced when computing phonon modes?
6.7 Why are phonon calculations so slow?
This is because an automated supercell approach to computing the dynamical
matrices is used. The idea is to compute the phonon modes on a relatively
coarse grid (eg. 4x4x4) and interpolate to arbitrary q-points. Always
choose the grid numbers to be multiples of each other i.e. 4x4x8 is efficient
but 2x3x5 will be very inefficient. By far the most time consuming step is the
diagonalisation of the supercell Hamiltonian. Therefore the planewave cut-off
variable rgkmax should be made as small as possible (but no smaller!).
Dynamical matrix calculations can also be spread across a cluster of nodes which
share the same directory.
6.8 Is the LO-TO splitting in polar semiconductors correctly calculated?
No. This will be introduced in the near future.
7. Basis set and convergence
7.1 Which paramemters should I check for convergence?
rgkmax, gmaxvr, lmaxapw, lmaxvr, nempty,
lradstp are the most important ones. Checking for convergence with respect to
the k-point mesh is mandatory, as in any DFT calculation. Depending on the mesh,
swidth has to be chosen appropriately. A small swidth is necessary for
accurate gaps and magnetic moments.
7.2 Are the default cutoffs optimal?
No. They are reasonable defaults for many systems, but they are usually a compromise
between accuracy and speed.
7.3 Should I try to converge the total energy?
No. Only energy differences are meaningful and should be converged. Energy differences
converge much faster than the total energy itself.
7.4 How can I check the completeness of the augmentation?
Run several calculations of the same system with different muffin-tin radii. Change
these in the species files. If the augmentation is sufficiently complete, your results
(bands, gaps, moments, etc.) will essentially not depend on the MT radii within some