ELk Code FAQ

0.1 How can I obtain the code?
Download the latest version from SourceForge.net.

0.2 How can I be informed of new releases?
You can join the mailing list for notification of new releases, features and bugs in the code.

0.3 Does the name Elk mean anything?
It could stand for "Electrons 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 atomic. In particular, the electron rest mass, m_e ; the electronic charge, e ; Planck's constant over two pi, h-bar; and the electrostatic constant 4πε₀ 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 Protex perl script and compiled into a single document. Use the command make doc in the elk/src 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 experimental 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.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 elk 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 GotoBLAS or ATLAS can be used instead.

1.5 Which fast Fourier transform (FFT) package does Elk use?
It uses a modified version of FFTPACK5. This 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 for details.

2.1 How do I restart a calculation?
Set tasks=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 touch STOP). The code will then complete the current self-consistent loop, write the file STATE.OUT and stop.

2.3 How can I force a running calculation to write STATE.OUT?
Create an empty file WRITE in the running directory (use touch WRITE). You will know STATE.OUT has been written when WRITE disappears.

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 the k-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 energy
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 autormt=.true. in elk.in; or increasing the smallest and decreasing the largest rmt in the species files
• gmaxvr is too small: increase it in small steps until the warning disappears.
• rgkmax is too large: decrease it in small steps until the warning disappears.

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 (scrpath) 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 STATE.OUT?
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 updated however.

2.10 How can I make the code write STATE.OUT after every iteration?
Normally, the density and potentials are written to the file STATE.OUT only after completion of the self-consistent loop. By setting nwrite to a positive integer the file will be written during the loop every nwrite iterations.

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 XCrySDen or V_Sim. The Spacegroup utility also allows for the construction of supercells. See the Spacegroup manual 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 to .true. 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 set autormt to be .true. and the optimal radii will be determined automatically.

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 to .true. 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.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 Grace 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 OpenDX for all types of 2 and 3D plots.

4.3 Where is the Fermi level in the DOS, band structure and Fermi surface plots?
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 elk/docs/Brillouin directory). These vertices should be entered in the plot1d block. After performing the calculation you should have the files BAND.OUT and BANDLINES.OUT 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 particular atom?
Yes, use tasks=21. Band structures are then plotted to files BAND_Sss_Aaaaa.OUT 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 ngridk, or both.

4.7 How do I plot a total or partial density of states (DOS)?
Set tasks 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 and tasks 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?
The Y_lm projections of the atomic site resolved DOS are arranged in logical order in the PDOS files, namely: (l,m) = (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 projected DOS is found in blocks 5-9.

4.9 How can the partial DOS be resolved into irreducible representations (IR) like e_g and t_2g?
The partial DOS file PDOS_Sss_Aaaaa.OUT contains (l,m)-resolved projections which correspond to the IR of the site symmetry for atom aaaa 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 setting lmirep to .false.

4.10 Can I get a DOS summed over m quantum numbers?
Yes. By default, the partial density of states is resolved over (l,m) quantum numbers. If you want the partial DOS to be summed over m (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, lmaxvr, rgkmax, gmaxvr, 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 plot.

4.13 How can I obtain orbital occupation numbers after performing a LDA+U calculation?
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 see this Elk Wiki.

5.1 How is magnetism implemented in the code?
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, v_s. 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 σ⋅B_s. 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 direction.

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?
Yes.

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 more.

5.12 What are the units of the external magnetic field used by the code?
The relevant term in the Hamiltonian is
H_B = g_e/(4c) σ·B_ext,
where g_e is the electron g-factor. In terms of SI units, one unit of bfieldc is equal to 1715.255397557 Tesla.

6.1 Do forces work with spin-orbit coupling and non-collinear magnetism?
Yes.

6.2 How can I constrain the positions of particular atoms during a structural optimisation run?
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?
Yes.

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.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 reasonable range.