Questaal Home

Setting up the lmf ctrl file and running calculations

The behaviour of all the Questaal codes is governed by the ctrl file. This file has a common format that is understood by all the programs in the Questaal suite, although there are many keys/options that are specific to different methods. Filenames in Questaal are organised such that the root of the filename labels the kind of file and the extension is a project identifier – this makes finding files from the command line easier.

There are many tokens for the control file and they are arranged into sections: options relating to the Hamiltonian belong in the “HAM” section, structural details in “STRUC”, species data in “SPEC”, etc. The codes provide quick reference to the input tokens via the command line flag “–input”; other command line options are show by invoking the program with “–help.”

lmf --input | less 
lmf --help

Minimal ctrl file for GaN

# wurtzite GaN


 XCFUN=1  # Ceperley-Alder LDA 

 NKABC=1 1 1


 ALAT=6.02822874  # expt. lattice constant 3.19A,c/a=1.6266458
 PLAT=sqrt(3)/2 -1/2        0
      0          1          0
      0          0          1.6266458

 ATOM=Ga Z=31 R=2.0 LMX=2
 ATOM=N  Z=7  R=1.6 LMX=2

 ATOM=Ga XPOS=2/3 1/3 0
 ATOM=Ga XPOS=1/3 2/3 1/2
 ATOM=N  XPOS=2/3 1/3 0.377
 ATOM=N  XPOS=1/3 2/3 0.877


  • white space is unimportant (but subsections must be indented if on a new line)
  • tokens are, eg: STRUC_NBAS, HAM_AUTOBAS_LMTO
  • options after an ATOM relate to that species until the next ATOM is found
  • XPOS for fractions of the lattice vectors, POS for cartesian (multiples of ALAT)


  • muffin-tin radius R – guessed manually; improve later
  • basis parameters – guessed automatically using AUTOBAS
  • lmax for the basis – choice
  • resolution of the smooth grid, GMAX – still needs to be completed

Checking the structure with lmchk

lmchk is useful for printing geometry, neighbour distances, angles, etc. Also for determining good muffin-tin estimates:

lmchk --getwsr gan | less

recommends for R:

 makrm0: initial MT radii from first estat potential maximum
  site   spec            rmt       rmt-       rmt-       rold   lock
                                  <spec avg>  spec-min
    1    1:Ga           2.0153     0.0000     0.0000     2.0000
    2    1:Ga           2.0153     0.0000     0.0000     2.0000
    3    2:N            1.6679     0.0000     0.0000     1.6000
    4    2:N            1.6679     0.0000     0.0000     1.6000

use these to improve the input file

Setting up the basis using AUTOBAS

The basis consists of augmented smoothed-Hankel functions, up to ATOM_LMAX, described by RSMH and EH.

The default mechanism to setup the LMTO basis is based on fitting smoothed-Hankel basis functions to atomic eigenfunctions; these form one (kappa) set, giving RSMH and EH pairs for each upto LMAX. A second set is obtained from the first by adding a constant to EH. This approach is quite reliable and is fully automatic.

Automatic basis construction is activated by nonzero HAM_AUTOBAS_MTO and HAM_AUTOBAS_LMTO. If MTO=1,2 a one,two kappa set is generated. LMTO controls the size of the basis to be constructed (MTO=2 LMTO=5 ~ spd,spd).

The basis fitting is done by the free-atom code lmfa and is written to the file basp0.gan:

lmfa gan | tee lmfa_log

note: the main output for the codes is directed to standard out

 Species Ga:  Z=31  Qc=28  R=2.015300  Q=0
 mesh:   rmt=2.015300  rmax=49.443315  a=0.025  nr=389  nr(rmax)=517
  Pl=  4.5     4.5     4.5     4.5
  Ql=  2.0     1.0     0.0     0.0

  iter     qint         drho          vh0          rho0          vsum     beta
    1   31.000000   5.270E+03      155.0000    0.1542E+03      -62.2954   0.30
   51   31.000000   4.475E-05      302.3465    0.3447E+05     -166.2355   0.30

 sumev=-1.548744  etot=-3882.268124  eref=0.000000

 Free-atom wavefunctions:
 valence:      eval       node at      max at       c.t.p.   rho(r>rmt)
   4s      -0.67398         0.807       1.882       2.876     0.584682
   4p      -0.20078         0.869       2.519       5.145     0.830684
   4d       0.01244         1.254      33.370      49.443*    0.999995
   4f       0.01890         0.000      35.790      49.443*    1.000000

 core:        ecore       node at      max at       c.t.p.   rho(r>rmt)
   1s    -750.27749         0.000       0.032       0.065     0.000000
   2s     -92.58142         0.065       0.183       0.285     0.000000
   2p     -80.67737         0.000       0.146       0.309     0.000000
   3s     -10.83385         0.267       0.565       0.823     0.000217
   3p      -7.24866         0.240       0.567       0.979     0.000870
   3d      -1.42658         0.000       0.531       2.017     0.022565

 Optimise free-atom basis for species Ga, rmt=2.0153
 l  it    Rsm      Eh     stiffR   stiffE      Eval      Exact     Pnu    Ql
 0  24   1.916  -0.563     266.6    191.7   -0.67366  -0.67398    4.78   2.00
 1  13   1.859  -0.100     289.5   1726.4   -0.19838  -0.20078    4.63   1.00
 eigenvalue sum:  exact  -1.54874    opt basis  -1.54571    error 0.00303

 Make LMTO basis parms for species Ga to lmxb=3, rmt=2.0153  vbar=0
 l  it    Rsm       Eh        Eval      Exact     Pnu    Ql   Gmax
 0  10   1.344*  -0.251    -0.66140  -0.67398    4.78   2.00   5.5
 1  11   1.344*  -0.100*   -0.13808  -0.20078    4.63   1.00   5.9

lmfa writes to basp0.gan in order not to overwrite basp.gan. For running lmf, we need to rename the basp file. lmfa also writes free-atom densities to disk, which are overlapped to initialise the self-consistency cycle.

cat basp0.gan
mv basp0.gan basp.gan

For lmax=2 and LMTO=5, we have a (spd,spd) basis for both Ga and N.

 FREEAT:  estimate HAM_GMAX from RSMH:  GMAX=8.8

Once the basis has been constructed, it is possible to make an estimate of the GMAX necessary to describe the basis. This is printed at the end of the lmfa output. We need to edit the ctrl file with the recommended value (8.8) and to specify the k-mesh (say 10^3).

Running lmf

mpirun -np 24 lmf gan | tee lmf_log

This should converge in 9 iterations in about ~3s (parallelism is over k-points). Convergence progress can be monitored using the save file:

h ehf=-7983.012975 ehk=-7982.2720881
i ehf=-7982.9833865 ehk=-7982.2611507
i ehf=-7982.7731994 ehk=-7982.6602763
i ehf=-7982.7336772 ehk=-7982.7335436
i ehf=-7982.733418 ehk=-7982.7333624
i ehf=-7982.7333979 ehk=-7982.7333642
i ehf=-7982.7333797 ehk=-7982.7333467
i ehf=-7982.7333582 ehk=-7982.7333237
c ehf=-7982.7333627 ehk=-7982.7333263

^h – density started using overlapped free atoms ^i – iterating ^c – converged

ehf – Harris-Foulkes total energy (depends on input density and not the output density) ehk – Hohenburg-Kohn total energy (Ryd)

and the gap is 2.04087 eV. (grep gap lmf_log).

main output files:

  • atm.gan (from lmfa) – atomic charges and potentials (ASCII)
  • basp.gan (from lmfa) – basis definition (like ctrl)
  • log.gan – terse output information (ASCII)
  • mixm.gan – scf mixing history (binary)
  • moms.gan – weights of overlap matrix (binary)
  • rst.gan – restart file (binary)
  • save.gan – scf progress (ASCII)
  • wkp.gan – k-point weights (binary)

Adding local orbitals


Rerunning lmfa picks up the shallow Ga 3d:

 Search for possible local orbitals
 LO l=2: qistl=0.0226 exceeds QLOC, autogen PZ=13.91024 use PZ=13.91024
 LO l=2: e=-1.43 exceeds ELOC, autogen PZ=13.91024 use PZ=13.91024

basp0.gan reads:

 Ga RSMH= 1.344 1.344 1.344 EH= -0.1 -0.1 -0.1 RSMH2= 1.344 1.344 1.344 EH2= -0.9 -0.9 -0.9 PZ= 0 0 13.9102
 N RSMH= 1.002 0.909 1.112 EH= -0.191 -0.1 -0.1 RSMH2= 1.002 0.909 1.112 EH2= -0.991 -0.9 -0.9

copy basp0 to basp and run lmfa again (change in core charge number); delete mixm, rst.gan.

grep ^c save.gan
c ehf=-7982.7333627 ehk=-7982.7333263
c ehf=-7982.7921134 ehk=-7982.7920645

the gap becomes 1.97349 eV.

Big basis – lmax

increase lmax to 3 to setup an (spdf,spd) basis set re-run lmfa mv basp0.gan basp.gan

 Ga RSMH= 1.344 1.344 1.344 1.344 EH= -0.1 -0.1 -0.1 -0.1 RSMH2= 1.344 1.344 1.344 EH2= -0.9 -0.9 -0.9 PZ= 0 0 13.9102
 N RSMH= 1.002 0.909 1.112 1.112 EH= -0.191 -0.1 -0.1 -0.1 RSMH2= 1.002 0.909 1.112 EH2= -0.991 -0.9 -0.9

lmf now gives:

grep ^c save.gan
c ehf=-7982.7333627 ehk=-7982.7333263
c ehf=-7982.7921134 ehk=-7982.7920645
c ehf=-7982.8033203 ehk=-7982.8032586

gap 1.93483 eV.

Bigger basis – combined LMTO/APW

additional ctrl flags:


plane-wave cut-offs in Ryd.


 suham:  q-dependent PW basis with  Emin = 0 < E < 2.
         Est. min,max PW dimension = 11,23.  Use npwpad = 3 => ndham = 136
c ehf=-7982.8080754 ehk=-7982.8080078

gap 1.93285 eV.


 suham:  q-dependent PW basis with  Emin = 0 < E < 4.
         Est. min,max PW dimension = 30,51.  Use npwpad = 4 => ndham = 165
c ehf=-7982.8093846 ehk=-7982.8093305

gap 1.93235 eV.


  • usually PWEMAX=1 or 2 Ryd is effective for open structures
  • over-completeness becomes a problem (see OVEPS)
  • for large APW component, reduce R, remove MTOs.

Ctrl preprocessor directives

The Questaal ctrl is interpreted by a preprocessor which allows variable substitution, looping, arithmetic, flow control, … Preprocessor directives are preceded by the “%” symbol, which introduces constants and variables, and braces “{}” which are then expanded.

First use of this is to shadow ctrl options and make them accessible from the command line, eg, for checking the k-convergence:

% const nk=10

The variable nk can now be changed as an argument to lmf, eg:

lmf -vnk=12 gan

NKABC accepts either one number or three numbers to describe an M.P. grid – if one number, then all directions have the same number of k-divisions. If a negative number is given, this is interpreted as the desired number in the total BZ, and a set of divisions proportional to the length of the inverse cell vectors will be generated that is closest to this target. In this case, this reduces the number of divisions along the (longer) c-axis.

for i in -200 -400 -800 -2000 -4000 -8000
mpirun -np 24 lmf -vnk=$i gan |grep BZMESH


i nk=-200 ehf=-7982.8032021 ehk=-7982.8031549
c nk=-200 ehf=-7982.8032019 ehk=-7982.803155
i nk=-400 ehf=-7982.8032963 ehk=-7982.8032494
c nk=-400 ehf=-7982.8032962 ehk=-7982.8032495
i nk=-800 ehf=-7982.8033115 ehk=-7982.8032652
c nk=-800 ehf=-7982.8033099 ehk=-7982.8032638
i nk=-2e3 ehf=-7982.8033133 ehk=-7982.8032674
c nk=-2e3 ehf=-7982.8033133 ehk=-7982.8032676
i nk=-4e3 ehf=-7982.8033133 ehk=-7982.8032678
c nk=-4e3 ehf=-7982.8033131 ehk=-7982.8032678
i nk=-8e3 ehf=-7982.803313 ehk=-7982.8032679
c nk=-8e3 ehf=-7982.8033129 ehk=-7982.803268

Obviously metals require more attention; by default Questaal uses metallic occupations and tetrahedron integration: only in special circumstances would you need to change these. To display the as-preprocessed ctrl file, use: lmf –showp gan .

Automatic setup with blm

sample cif file:

cp /home/vol05/tmp15/gan_icsd_34476.cif .
cif2cell icsd_34476.cif
cif2cell icsd_34476.cif > out
cif2init out

cif2cell is a (very useful!) external tool (cif2cell at sourceforge); cif2init converts the output to an “init” file, which contains basic structural information than can be read by blm. blm generates R and sets various variables depending on the specific kind of calculation.

mv init init.gan
blm gan
mv actrl.gan ctrl.gan

(actrl is to avoid potentially clobbering ctrl)

vi ctrl.gan
cat site.gan
lmfa --usebasp gan
cat basp.gan
mpirun -np 24 lmf -vgmax=8.8 -vnkabc=10 gan |tee log_lmf

the result is different because the default basis is (spd,sp) and because XCFUN=2 (von Barth-Hedin) has been chosen.

mpirun -np 24 lmf -vgmax=8.8 -vnkabc=10 -vlxcf=1 gan |tee log_lmf2


  • check the preprocessor names/values generated by blm
  • -v is used to change the values of preprocessor variables
  • blm has lots of options – tune the setup at the point of blm invocation
  • blm –help
  • –express=0 gives more “old fashioned” input files