1 2! Copyright (C) 2002-2011 J. K. Dewhurst, S. Sharma and C. Ambrosch-Draxl. 3! This file is distributed under the terms of the GNU General Public License. 4! See the file COPYING for license details. 5 6! main routine for the Elk code 7program elk 8use modmain 9use modmpi 10use modomp 11use modvars 12use modramdisk 13implicit none 14! local variables 15logical exist 16integer itask 17! initialise MPI execution environment 18call mpi_init(ierror) 19! duplicate mpi_comm_world 20call mpi_comm_dup(mpi_comm_world,mpicom,ierror) 21! determine the number of MPI processes 22call mpi_comm_size(mpicom,np_mpi,ierror) 23! determine the local MPI process number 24call mpi_comm_rank(mpicom,lp_mpi,ierror) 25! determine if the local process is the master 26if (lp_mpi.eq.0) then 27 mp_mpi=.true. 28 write(*,*) 29 write(*,'("Elk code version ",I1.1,".",I1.1,".",I2.2," started")') version 30else 31 mp_mpi=.false. 32end if 33! read input files 34call readinput 35! initialise OpenMP variables 36call omp_init 37! initialise the MKL library 38call mkl_init 39! initialise the OpenBLAS library 40call oblas_init 41! initialise the BLIS library 42call blis_init 43if (mp_mpi) then 44 write(*,*) 45 write(*,'("Number of MPI processes : ",I6)') np_mpi 46 write(*,'("Number of OpenMP threads per MPI process : ",I4)') maxthd 47 write(*,'("Total number of threads (MPI x OpenMP) : ",I6)') np_mpi*maxthd 48 write(*,'("Maximum OpenMP nesting level : ",I4)') maxlvl 49 write(*,'("Number of OpenMP threads at first nesting level : ",I4)') maxthd1 50 write(*,'("Number of MKL threads : ",I4)') maxthdmkl 51end if 52! delete the VARIABLES.OUT file 53call delvars 54! write version number to VARIABLES.OUT 55call writevars('version',nv=3,iva=version) 56! check if Elk is already running in this directory 57if (mp_mpi) then 58 inquire(file='RUNNING',exist=exist) 59 if (exist) then 60 write(*,*) 61 write(*,'("Info(elk): several copies of Elk may be running in this path")') 62 write(*,'("(this could be intentional, or result from a previous crash,")') 63 write(*,'(" or arise from an incorrect MPI compilation)")') 64 else 65 open(50,file='RUNNING') 66 close(50) 67 end if 68end if 69! initialise the RAM disk if required 70if (ramdisk) then 71 call initrd 72 if (mp_mpi) then 73 write(*,*) 74 write(*,'("RAM disk enabled")') 75 end if 76end if 77! perform the tasks 78do itask=1,ntasks 79 task=tasks(itask) 80 if (mp_mpi) then 81 write(*,*) 82 write(*,'("Info(elk): current task : ",I6)') task 83 end if 84! synchronise MPI processes 85 call mpi_barrier(mpicom,ierror) 86! check if task can be run with MPI 87 if (lp_mpi.gt.0) then 88 if (any(task.eq.[0,1,2,3,5,15,16,28,29,61,62,63,110,120,135,136,162,170, & 89 180,185,188,200,201,205,208,240,241,270,271,300,320,330,331,350,351,371, & 90 372,373,440,460,461,462,463,471,478,600,610,620,630,640,700,701,800, & 91 801])) then 92 continue 93 else 94 write(*,'("Info(elk): MPI process ",I6," idle for task ",I6)') lp_mpi,task 95 cycle 96 end if 97 end if 98! write task to VARIABLES.OUT 99 call writevars('task',iv=task) 100 select case(task) 101 case(0,1) 102 call gndstate 103 case(2,3) 104 call geomopt 105 case(5) 106 call hartfock 107 case(10) 108 call writedos 109 case(14) 110 call writesf 111 case(15,16) 112 call writelsj 113 case(20,21,22,23) 114 call bandstr 115 case(25) 116 call effmass 117 case(28,29) 118 call mae 119 case(31,32,33) 120 call rhoplot 121 case(41,42,43) 122 call potplot 123 case(51,52,53) 124 call elfplot 125 case(61,62,63,162) 126 call wfplot 127 case(65) 128 call wfcrplot 129 case(68) 130 call rdstatus 131 case(71,72,73,81,82,83,141,142,143,151,152,153) 132 call vecplot 133 case(91,92,93) 134 call dbxcplot 135 case(100,101) 136 call fermisurf 137 case(102) 138 call fermisurfbxsf 139 case(105) 140 call nesting 141 case(110) 142 call mossbauer 143 case(115) 144 call writeefg 145 case(120) 146 call writepmat 147 case(121) 148 call dielectric 149 case(122) 150 call moke 151 case(125) 152 call nonlinopt 153 case(130) 154 call writeexpmat 155 case(135) 156 call writewfpw 157 case(140) 158 call elnes 159 case(150) 160 call writeevsp 161 case(160) 162 call torque 163 case(170) 164 call writeemd 165 case(171,172,173) 166 call emdplot 167 case(180) 168 call writeepsinv 169 case(185) 170 call writehmlbse 171 case(186) 172 call writeevbse 173 case(187) 174 call dielectric_bse 175 case(190) 176 call geomplot 177 case(195) 178 call sfacrho 179 case(196) 180 call sfacmag 181 case(200,201,202) 182 call phononsc 183 case(205) 184 call phonon 185 case(208,209) 186 call bornechg 187 case(210) 188 call phdos 189 case(220) 190 call phdisp 191 case(230) 192 call writephn 193 case(240,241) 194 call ephcouple 195 case(245) 196 call phlwidth 197 case(250) 198 call alpha2f 199 case(260) 200 call eliashberg 201 case(270,271) 202 call gndsteph 203 case(280) 204 call ephdos 205 case(285) 206 call aceplot 207 case(300) 208 call rdmft 209 case(320) 210 call tddftlr 211 case(330,331) 212 call tddftsplr 213 case(341,342,343) 214 call wxcplot 215 case(350,351,352) 216 call spiralsc 217 case(371,372,373) 218 call jprplot 219 case(400) 220 call writetmdu 221 case(430) 222 call writestrain 223 case(440) 224 call writestress 225 case(450) 226 call genafieldt 227 case(460,461,462,463) 228 call tddft 229 case(471) 230 call rhosplot 231 case(478) 232 call bornectd 233 case(480,481) 234 call dielectric_tdrt 235 case(500) 236 call testcheck 237 case(550) 238 call writew90 239 case(600) 240 call gwsefm 241 case(610) 242 call gwspecf 243 case(620) 244 call gwbandstr 245 case(630) 246 call gwscrho 247 case(640) 248 call gwdmat 249 case(700,701) 250 call gndstulr 251 case(731,732,733) 252 call rhouplot 253 case(741,742,743) 254 call potuplot 255 case(771,772,773) 256 call maguplot 257 case(800,801) 258!*********** 259 260 case default 261 write(*,*) 262 write(*,'("Error(elk): task not defined : ",I8)') task 263 write(*,*) 264 stop 265 end select 266! reset the OpenMP thread variables 267 call omp_reset 268! close all opened files 269 call closefiles 270end do 271if (mp_mpi) then 272 open(50,file='RUNNING') 273 close(50,status='DELETE') 274 write(*,*) 275 write(*,'("Elk code stopped")') 276end if 277! terminate MPI execution environment 278call mpi_finalize(ierror) 279end program 280 281!BOI 282! !TITLE: {\huge{\sc The Elk Code Manual}}\\ \Large{\sc Version 7.2.42}\\ \vskip 20pt \includegraphics[height=1cm]{elk_silhouette.pdf} 283! !AUTHORS: {\sc J. K. Dewhurst, S. Sharma} \\ {\sc L. Nordstr\"{o}m, F. Cricchio, O. Gr\aa n\"{a}s} \\ {\sc E. K. U. Gross} 284! !AFFILIATION: 285! !INTRODUCTION: Introduction 286! Welcome to the Elk Code! Elk is an all-electron full-potential linearised 287! augmented-plane-wave (FP-LAPW) code for determining the properties of 288! crystalline solids. It was developed originally at the 289! Karl-Franzens-Universit\"{a}t Graz as part of the EXCITING EU Research and 290! Training Network project\footnote{EXCITING code developed under the Research 291! and Training Network EXCITING funded by the EU, contract No. 292! HPRN-CT-2002-00317}. The guiding philosophy during the implementation of the 293! code was to keep it as simple as possible for both users and developers 294! without compromising on its capabilities. All the routines are released 295! under either the GNU General Public License (GPL) or the GNU Lesser General 296! Public License (LGPL) in the hope that they may inspire other scientists to 297! implement new developments in the field of density functional theory and 298! beyond. 299! 300! \section{Acknowledgments} 301! Lots of people contributed to the Elk code with ideas, checking and testing, 302! writing code or documentation and general encouragement. They include 303! Claudia Ambrosch-Draxl, Clas Persson, Fredrik Bultmark, Christian Brouder, 304! Rickard Armiento, Andrew Chizmeshya, Per Anderson, Igor Nekrasov, Sushil 305! Auluck, Frank Wagner, Fateh Kalarasse, J\"{u}rgen Spitaler, Stefano 306! Pittalis, Nektarios Lathiotakis, Tobias Burnus, Stephan Sagmeister, 307! Christian Meisenbichler, S\'{e}bastien Leb\`{e}gue, Yigang Zhang, Fritz 308! K\"{o}rmann, Alexey Baranov, Anton Kozhevnikov, Shigeru Suehara, Frank 309! Essenberger, Antonio Sanna, Tyrel McQueen, Tim Baldsiefen, Marty Blaber, 310! Anton Filanovich, Torbj\"{o}rn Bj\"{o}rkman, Martin Stankovski, Jerzy 311! Goraus, Markus Meinert, Daniel Rohr, Vladimir Nazarov, Kevin Krieger, Pink 312! Floyd, Arkardy Davydov, Florian Eich, Aldo Romero Castro, Koichi Kitahara, 313! James Glasbrenner, Konrad Bussmann, Igor Mazin, Matthieu Verstraete, David 314! Ernsting, Stephen Dugdale, Peter Elliott, Marcin Dulak, Jos\'{e} A. Flores 315! Livas, Stefaan Cottenier, Yasushi Shinohara, Michael Fechner, Yaroslav 316! Kvashnin, Tristan M\"uller, Arsenii Gerasimov, Manh Duc Le, Jon Lafuente 317! Bartolom\'{e}, Ren\'{e} Wirnata, Jagdish Kumar, Andrew Shyichuk, Nisha 318! Singh, Pietro Bonfa, Ronald Cohen and Alyn James. Special mention of David 319! Singh's very useful book on the LAPW method\footnote{D. J. Singh, 320! {\it Planewaves, Pseudopotentials and the LAPW Method} (Kluwer Academic 321! Publishers, Boston, 1994).} must also be made. Finally we would like to 322! acknowledge the generous support of Karl-Franzens-Universit\"{a}t Graz, the 323! EU Marie-Curie Research Training Networks initiative, the Max Born Institute 324! and the Max Planck Society. 325! 326! \vspace{24pt} 327! Kay Dewhurst\newline 328! Sangeeta Sharma\newline 329! Lars Nordstr\"{o}m\newline 330! Francesco Cricchio\newline 331! Oscar Gr\aa n\"{a}s\newline 332! Hardy Gross 333! 334! \vspace{12pt} 335! Berlin, Halle, Jerusalem and Uppsala, June 2021 336! \newpage 337! 338! \section{Units} 339! Unless explicitly stated otherwise, Elk uses atomic units. In this system 340! $\hbar=1$, the electron mass $m=1$, the Bohr radius $a_0=1$ and the electron 341! charge $e=1$ (note that the electron charge is positive, so that the atomic 342! numbers $Z$ are negative). Thus the atomic unit of length is 343! 0.529177210903(80) \AA, and the atomic unit of energy is the Hartree which 344! equals 27.211386245988(53) eV. The unit of the external magnetic fields is 345! defined such that one unit of magnetic field in {\tt elk.in} equals 346! 1715.255541 Tesla. 347! 348! \section{Compiling and running Elk} 349! \subsection{Compiling the code} 350! Unpack the code from the archive file. Run the command 351! \begin{verbatim} 352! setup 353! \end{verbatim} 354! in the {\tt elk} directory and select the appropriate system and compiler. 355! We highly recommend that you edit the file {\tt make.inc} and tune the 356! compiler options for your computer system. In particular, use of 357! machine-optimised BLAS/LAPACK libraries can result in significant increase 358! in performance, but make sure they are of version $3.x$. Following this, run 359! \begin{verbatim} 360! make 361! \end{verbatim} 362! This will hopefully compile the entire code and all the libraries into one 363! executable, {\tt elk}, located in the {\tt elk/src} directory. It will also 364! compile two useful auxiliary programs, namely {\tt spacegroup} for producing 365! crystal geometries from spacegroup data and {\tt eos} for fitting equations 366! of state to energy-volume data. If you want to compile everything all over 367! again, then run {\tt make clean} from the {\tt elk} directory, followed by 368! {\tt make}. 369! \subsubsection{Parallelism in Elk} 370! Three forms of parallelism are implemented in Elk, and all can be used in 371! combination with each other, with efficiency depending on the particular 372! task, crystal structure and computer system. You may need to contact your 373! system administrator for assistance with running Elk in parallel. 374! \begin{enumerate} 375! \item 376! OpenMP works for symmetric multiprocessors, i.e. computers that have many 377! cores with the same unified memory accessible to each. It is enabled by 378! setting the appropriate command-line options (e.g. {\tt -qopenmp} for the 379! Intel compiler) before compiling, and also at runtime by the environment 380! variable 381! \begin{verbatim} 382! export OMP_NUM_THREADS=n 383! \end{verbatim} 384! where n is the number of cores available on a particular node. The same can 385! be accomplished in {\tt elk.in} with 386! \begin{verbatim} 387! maxthd 388! n 389! \end{verbatim} 390! In addition, some vendor-supplied BLAS/LAPACK libraries use OpenMP 391! internally. The maximum number of threads used for LAPACK operations by 392! Intel's MKL can be set with 393! \begin{verbatim} 394! maxthdmkl 395! n 396! \end{verbatim} 397! \item 398! The message passing interface (MPI) is particularly suitable for running 399! Elk across multiple nodes of a cluster, with scaling to hundreds of 400! processors possible. To enable MPI, comment out the lines indicated in 401! {\tt elk/make.inc}. Then run {\tt make clean} followed by {\tt make}. If 402! $y$ is the number of nodes and $x$ is the number of cores per node, then at 403! runtime envoke 404! \begin{verbatim} 405! mpirun -np z ./elk 406! \end{verbatim} 407! where $z=x y$ is the total number of cores available on the machine. 408! Highest efficiency is obtained by using hybrid parallelism with OpenMP on 409! each node and MPI across nodes. This can be done by compiling the code 410! using the MPI Fortran compiler in combination with the OpenMP command-line 411! option. At runtime set {\tt export OMP\_NUM\_THREADS=x} and start the MPI 412! run with {\em one process per node} as follows 413! \begin{verbatim} 414! mpirun -pernode -np y ./elk 415! \end{verbatim} 416! The number of MPI processes is reported in the file {\tt INFO.OUT} which 417! serves as a check that MPI is running correctly. Note that version 2 of the 418! MPI libraries is required to run Elk. 419! \item 420! Phonon calculations use a simple form of parallelism by just examining the 421! run directory for dynamical matrix files. These files are of the form 422! \begin{verbatim} 423! DYN_Qqqqq_qqqq_qqqq_Sss_Aaa_Pp.OUT 424! \end{verbatim} 425! and contain a single row of a particular dynamical matrix. Elk simply finds 426! which {\tt DYN} files do not exist, chooses one and runs it. This way many 427! independent runs of Elk can be started in the same directory on a networked 428! file system (NFS), and will run until all the dynamical matrices files are 429! completed. Should a particular run crash, then delete the associated empty 430! {\tt DYN} file and rerun Elk. 431! \end{enumerate} 432! 433! \subsection{Memory requirements} 434! Elk is a memory-bound code and runs best on processors with large caches and 435! a large number of memory channels per core. Some tasks in Elk require a 436! considerable amount of memory which can exceed the physical memory of the 437! computer. In such cases, the number of threads at the first nesting level 438! can be reduced with (for example) 439! \begin{verbatim} 440! maxthd1 441! -4 442! \end{verbatim} 443! which restricts the number of threads at the first nesting level to 444! maxthd/4. Deeper nesting levels, which generally require less memory, will 445! still utilise the full compliment of available threads. 446! \subsubsection{Stack space} 447! The latest versions of Elk use stack space aggressively. This is because 448! accessing variables is faster on the stack than on the heap. This can, 449! however, result in the code crashing as threads run out of their stack 450! space. To avoid this, increase the stack size for each OpenMP thread with 451! (for example) 452! \begin{verbatim} 453! export OMP_STACKSIZE=64M 454! \end{verbatim} 455! before running the code. 456! 457! \subsection{Linking with the Libxc functional library} 458! Libxc is the ETSF library of exchange-correlation functionals. Elk can use 459! the complete set of LDA and GGA functionals available in Libxc as well as 460! the potential-only metaGGA's. In order to enable this, first download and 461! compile Libxc version 5. This should have produced the files {\tt libxc.a} 462! and {\tt libxcf90.a}. Copy these files to the {\tt elk/src} directory and 463! then uncomment the lines indicated for Libxc in the file {\tt elk/make.inc}. 464! Once this is done, run {\tt make clean} followed by {\tt make}. To select a 465! particular functional of Libxc, use the block 466! \begin{verbatim} 467! xctype 468! 100 nx nc 469! \end{verbatim} 470! where {\tt nx} and {\tt nc} are, respectively, the numbers of the exchange 471! and correlation functionals in the Libxc library. See the file 472! {\tt elk/src/libxcf90.f90} for a list of the functionals and their 473! associated numbers. 474! 475! \subsection{Running the code} 476! As a rule, all input files for the code are in lower case and end with the 477! extension {\tt .in}. All output files are uppercase and have the extension 478! {\tt .OUT}. For most cases, the user will only need to modify the file 479! {\tt elk.in}. In this file input parameters are arranged in blocks. 480! Each block consists of a block name on one line and the block variables on 481! subsequent lines. Almost all blocks are optional: the code uses reasonable 482! default values in cases where they are absent. Blocks can appear in any 483! order, if a block is repeated then the second instance is used. Comment 484! lines can be included in the input file and begin with the {\tt !} 485! character. 486! 487! \subsubsection{Species files} 488! The only other input files are those describing the atomic species which go 489! into the crystal. These files are found in the {\tt species} directory and 490! are named with the element symbol and the extension {\tt .in}, for example 491! {\tt Sb.in}. They contain parameters like the atomic charge, mass, 492! muffin-tin radius, occupied atomic states and the type of linearisation 493! required. Here as an example is the copper species file {\tt Cu.in}: 494! \begin{verbatim} 495! 'Cu' : spsymb 496! 'copper' : spname 497! -29.0000 : spzn 498! 115837.2716 : spmass 499! 0.371391E-06 2.0000 34.8965 500 : rminsp, rmt, rmaxsp, nrmt 500! 10 : nstsp 501! 1 0 1 2.00000 T : nsp, lsp, ksp, occsp, spcore 502! 2 0 1 2.00000 T 503! 2 1 1 2.00000 T 504! 2 1 2 4.00000 T 505! 3 0 1 2.00000 T 506! 3 1 1 2.00000 F 507! 3 1 2 4.00000 F 508! 3 2 2 4.00000 F 509! 3 2 3 6.00000 F 510! 4 0 1 1.00000 F 511! 1 : apword 512! 0.1500 0 F : apwe0, apwdm, apwve 513! 1 : nlx 514! 2 2 : lx, apword 515! 0.1500 0 T : apwe0, apwdm, apwve 516! 0.1500 1 T 517! 4 : nlorb 518! 0 2 : lorbl, lorbord 519! 0.1500 0 F : lorbe0, lorbdm, lorbve 520! 0.1500 1 F 521! 1 2 522! 0.1500 0 F 523! 0.1500 1 F 524! 2 2 525! 0.1500 0 F 526! 0.1500 1 F 527! 1 3 528! 0.1500 0 F 529! 0.1500 1 F 530! -2.8652 0 T 531! \end{verbatim} 532! The input parameters are defined as follows: 533! \vskip 6pt 534! {\tt spsymb} \\ 535! The symbol of the element. 536! \vskip 6pt 537! {\tt spname} \\ 538! The name of the element. 539! \vskip 6pt 540! {\tt spzn} \\ 541! Nuclear charge: should be negative since the electron charge is taken to be 542! postive in the code; it can also be fractional for purposes of doping. 543! \vskip 6pt 544! {\tt spmass} \\ 545! Nuclear mass in atomic units. 546! \vskip 6pt 547! {\tt rminsp}, {\tt rmt}, {\tt rmaxsp}, {\tt nrmt} \\ 548! Respectively, the minimum radius on logarithmic radial mesh; muffin-tin 549! radius; effective infinity for atomic radial mesh; and number of radial mesh 550! points to muffin-tin radius. 551! \vskip 6pt 552! {\tt nstsp} \\ 553! Number of atomic states. 554! \vskip 6pt 555! {\tt nsp}, {\tt lsp}, {\tt ksp}, {\tt occsp}, {\tt spcore} \\ 556! Respectively, the principal quantum number of the radial Dirac equation; 557! quantum number $l$; quantum number $k$ ($l$ or $l+1$); occupancy of atomic 558! state (can be fractional); {\tt .T.} if state is in the core and therefore 559! treated with the Dirac equation in the spherical part of the muffin-tin 560! Kohn-Sham potential. 561! \vskip 6pt 562! {\tt apword} \\ 563! Default APW function order, i.e. the number of radial functions and 564! therefore the order of the radial derivative matching at the muffin-tin 565! surface. 566! \vskip 6pt 567! {\tt apwe0}, {\tt apwdm}, {\tt apwve} \\ 568! Respectively, the default APW linearisation energy; the order of the energy 569! derivative of the APW radial function $\partial^m u(r)/\partial E^m$; and 570! {\tt .T.} if the linearisation energy is allowed to vary. 571! \vskip 6pt 572! {\tt nlx} \\ 573! The number of exceptions to the default APW configuration. These should be 574! listed on subsequent lines for particular angular momenta. In this example, 575! the fixed energy APW with angular momentum $d$ ({\tt lx} $=2$) is replaced 576! with a LAPW, which has variable linearisation energy. 577! \vskip 6pt 578! {\tt nlorb} \\ 579! Number of local-orbitals. 580! \vskip 6pt 581! {\tt lorbl}, {\tt lorbord} \\ 582! Respectively, the angular momentum $l$ of the local-orbital; and the order 583! of the radial derivative which goes to zero at the muffin-tin surface. 584! \vskip 6pt 585! {\tt lorbe0}, {\tt lorbdm}, {\tt lorbve} \\ 586! Respectively, the default local-orbital linearisation energy; the order of 587! the energy derivative of the local-orbital radial function; and {\tt .T.} if 588! the linearisation energy is allowed to vary. 589! 590! \subsubsection{Examples} 591! The best way to learn to use Elk is to run the examples included with the 592! package. These can be found in the {\tt examples} directory and use many of 593! the code's capabilities. The following section which describes all the input 594! parameters will be of invaluable assistance. 595! 596! \section{Input blocks} 597! This section lists all the input blocks available. It is arranged with the 598! name of the block followed by a table which lists each parameter name, what 599! the parameter does, its type and default value. A horizontal line in the 600! table indicates a new line in {\tt elk.in}. Below the table is a brief 601! overview of the block's function. 602! 603! \block{atoms}{ 604! {\tt nspecies} & number of species & integer & 0 \\ 605! \hline 606! {\tt spfname(i)} & species filename for species $i$ & string & - \\ 607! \hline 608! {\tt natoms(i)} & number of atoms for species $i$ & integer & - \\ 609! \hline 610! {\tt atposl(j,i)} & atomic position in lattice coordinates for atom $j$ 611! & real(3) & - \\ 612! {\tt bfcmt(j,i)} & muffin-tin external magnetic field in Cartesian 613! coordinates for atom $j$ & real(3) & -} 614! Defines the atomic species as well as their positions in the unit cell and 615! the external magnetic field applied throughout the muffin-tin. These fields 616! are used to break spin symmetry and should be considered infinitesimal as 617! they do not contribute directly to the total energy. Collinear calculations 618! are more efficient if the field is applied in the $z$-direction. One could, 619! for example, set up an antiferromagnetic crystal by pointing the field on 620! one atom in the positive $z$-direction and in the opposite direction on 621! another atom. If {\tt molecule} is {\tt .true.} then the atomic positions 622! are assumed to be in Cartesian coordinates. See also {\tt sppath}, 623! {\tt bfieldc} and {\tt molecule}. 624! 625! \block{autokpt}{ 626! {\tt autokpt} & {\tt .true.} if the $k$-point set is to be determined 627! automatically & logical & {\tt .false.}} 628! See {\tt radkpt} for details. 629! 630! \block{autolinengy}{ 631! {\tt autolinengy} & {\tt .true.} if the fixed linearisation energies are 632! to be determined automatically & logical & {\tt .false.}} 633! See {\tt dlefe} for details. 634! 635! \block{autoswidth}{ 636! {\tt autoswidth} & {\tt .true.} if the smearing parameter {\tt swidth} 637! should be determined automatically & logical & {\tt .false.}} 638! Calculates the smearing width from the $k$-point density, $V_{\rm BZ}/n_k$; 639! the valence band width, $W$; and an effective mass parameter, $m^{*}$; 640! according to 641! $$ \sigma=\frac{\sqrt{2W}}{m^{*}}\left(\frac{3}{4\pi} 642! \frac{V_{\rm BZ}}{n_k}\right)^{1/3}. $$ 643! The variable {\tt mstar} then replaces {\tt swidth} as the control parameter 644! of the smearing width. A large value of $m^{*}$ gives a narrower smearing 645! function. Since {\tt swidth} is adjusted according to the fineness of the 646! ${\bf k}$-mesh, the smearing parameter can then be eliminated. It is not 647! recommended that {\tt autoswidth} be used in conjunction with the 648! Fermi-Dirac smearing function, since the electronic temperature will then be 649! a function of the $k$-point mesh. See T. Bj\"orkman and O. Gr\aa n\"as, 650! {\it Int. J. Quant. Chem.} DOI: 10.1002/qua.22476 (2010) for details. See 651! also {\tt stype} and {\tt swidth}. 652! 653! \block{avec}{ 654! {\tt avec(1)} & first lattice vector & real(3) & $(1.0,0.0,0.0)$ \\ 655! \hline 656! {\tt avec(2)} & second lattice vector & real(3) & $(0.0,1.0,0.0)$ \\ 657! \hline 658! {\tt avec(3)} & third lattice vector & real(3) & $(0.0,0.0,1.0)$} 659! Lattice vectors of the crystal in atomic units (Bohr). 660! 661! \block{beta0}{ 662! {\tt beta0} & adaptive mixing parameter & real & $0.05$} 663! This determines how much of the potential from the previous self-consistent 664! loop is mixed with the potential from the current loop. It should be made 665! smaller if the calculation is unstable. See {\tt betamax} and also the 666! routine {\tt mixadapt}. 667! 668! \block{betamax}{ 669! {\tt betamax} & maximum adaptive mixing parameter & real & $0.5$} 670! Maximum allowed mixing parameter used in routine {\tt mixadapt}. 671! 672! \block{bfieldc}{ 673! {\tt bfieldc} & global external magnetic field in Cartesian coordinates & 674! real(3) & $(0.0,0.0,0.0)$} 675! This is a constant magnetic field applied throughout the entire unit cell 676! and enters the second-variational Hamiltonian as 677! $$ \frac{g_e}{4c}\,\vec{\sigma}\cdot{\bf B}_{\rm ext}, $$ 678! where $g_e$ is the electron $g$-factor. This field is normally used to break 679! spin symmetry for spin-polarised calculations and considered to be 680! infinitesimal with no direct contribution to the total energy. In cases 681! where the magnetic field is finite (for example when computing magnetic 682! response) the external ${\bf B}$-field energy reported in {\tt INFO.OUT} 683! should be added to the total by hand. This field is applied throughout the 684! entire unit cell. To apply magnetic fields in particular muffin-tins use the 685! {\tt bfcmt} vectors in the {\tt atoms} block. Collinear calculations are 686! more efficient if the field is applied in the $z$-direction. 687! 688! \block{broydpm}{ 689! {\tt broydpm} & Broyden mixing parameters $\alpha$ and $w_0$ & real & 690! $(0.4,0.15)$} 691! See {\tt mixtype} and {\tt mixsdb}. 692! 693! \block{c\_tb09}{ 694! {\tt c\_tb09} & Tran-Blaha constant $c$ & real & -} 695! Sets the constant $c$ in the Tran-Blaha '09 functional. Normally this is 696! calculated from the density, but there may be situations where this needs to 697! be adjusted by hand. See {\it Phys. Rev. Lett.} {\bf 102}, 226401 (2009). 698! 699! \block{chgexs}{ 700! {\tt chgexs} & excess electronic charge & real & $0.0$} 701! This controls the amount of charge in the unit cell beyond that required to 702! maintain neutrality. It can be set positive or negative depending on whether 703! electron or hole doping is required. 704! 705! \block{cmagz}{ 706! {\tt cmagz} & .true. if $z$-axis collinear magnetism is to be enforced & 707! logical & {\tt .false.}} 708! This variable can be set to .true. in cases where the magnetism is 709! predominantly collinear in the $z$-direction, for example a ferromagnet with 710! spin-orbit coupling. This will make the calculation considerably faster at 711! the slight expense of precision. 712! 713! \block{deltaem}{ 714! {\tt deltaem} & the size of the ${\bf k}$-vector displacement used when 715! calculating numerical derivatives for the effective mass tensor & real & 716! $0.025$} 717! See {\tt ndspem} and {\tt vklem}. 718! 719! \block{deltaph}{ 720! {\tt deltaph} & size of the atomic displacement used for calculating 721! dynamical matrices & real & $0.01$} 722! Phonon calculations are performed by constructing a supercell corresponding 723! to a particular ${\bf q}$-vector and making a small periodic displacement of 724! the atoms. The magnitude of this displacement is given by {\tt deltaph}. 725! This should not be made too large, as anharmonic terms could then become 726! significant, neither should it be too small as this can introduce numerical 727! error. 728! 729! \block{deltast}{ 730! {\tt deltast} & size of the change in lattice vectors used for calculating 731! the stress tensor & real & $0.005$} 732! The stress tensor is computed by changing the lattice vector matrix $A$ by 733! $$ A\rightarrow (1+\delta t\,e_i)A, $$ 734! where $dt$ is an infinitesimal equal in practice to {\tt deltast} and $e_i$ 735! is the $i^{\rm th}$ strain tensor. Numerical finite differences are used to 736! compute the stress tensor as the derivative of the total energy $dE_i/dt$. 737! 738! \block{dft+u}{ 739! {\tt dftu} & type of DFT+$U$ calculation & integer & 0 \\ 740! {\tt inpdftu} & type of input for DFT+U calculation & integer & 1 \\ 741! \hline 742! {\tt is} & species number & integer & - \\ 743! {\tt l} & angular momentum value & integer & -1 \\ 744! {\tt u} & the desired $U$ value & real & $0.0$ \\ 745! {\tt j} & the desired $J$ value & real & $0.0$} 746! This block contains the parameters required for an DFT+$U$ calculation, with 747! the list of parameters for each species terminated with a blank line. The 748! type of double counting required is set with the parameter {\tt dftu}. 749! Currently implemented are: 750! \vskip 6pt 751! \begin{tabularx}{\textwidth}[h]{lX} 752! 0 & No DFT+$U$ calculation \\ 753! 1 & Fully localised limit (FLL) \\ 754! 2 & Around mean field (AFM) \\ 755! 3 & An interpolation between FLL and AFM \\ 756! \end{tabularx} 757! \vskip 6pt 758! The type of input parameters is set with the parameter {\tt inpdftu}. 759! The current possibilities are: 760! \vskip 6pt 761! \begin{tabularx}{\textwidth}[h]{lX} 762! 1 & U and J \\ 763! 2 & Slater parameters \\ 764! 3 & Racah parameters \\ 765! 4 & Yukawa screening length \\ 766! 5 & U and determination of corresponding Yukawa screening length 767! \end{tabularx} 768! \vskip 6pt 769! See (amongst others) {\it Phys. Rev. B} {\bf 67}, 153106 (2003), 770! {\it Phys. Rev. B} {\bf 52}, R5467 (1995), {\it Phys. Rev. B} {\bf 60}, 771! 10763 (1999), and {\it Phys. Rev. B} {\bf 80}, 035121 (2009). 772! 773! \block{dlefe}{ 774! {\tt dlefe} & difference between the fixed linearisation energy and the 775! Fermi energy & real & $-0.1$} 776! When {\tt autolinengy} is {\tt .true.} then the fixed linearisation energies 777! are set to the Fermi energy plus {\tt dlefe}. 778! 779! \block{dncgga}{ 780! {\tt dncgga} & small constant used to stabilise non-collinear GGA & 781! real & $1\times 10^{-8}$} 782! This small constant, $d$, is required in order to remove the infinite 783! gradients obtained when using `Kubler's trick' in conjunction with GGA and 784! non-collinear magnetism. It is applied by calculating the up and down 785! densities as 786! $$ \rho^{\uparrow}({\bf r})=\rho({\bf r})+\widetilde{m}({\bf r}) 787! \qquad \rho^{\downarrow}({\bf r})=\rho({\bf r})-\widetilde{m}({\bf r}), $$ 788! where $\widetilde{m}({\bf r})=\sqrt{{\bf m}^2({\bf r})+d}$, 789! and should be taken as the smallest value for which the exchange-correlation 790! magnetic field ${\bf B}_{\rm xc}$ is smooth. 791! 792! \block{dosmsum}{ 793! {\tt dosmsum} & {\tt .true.} if the partial DOS is to be summed over $m$ & 794! logical & {\tt .false.}} 795! By default, the partial density of states is resolved over $(l,m)$ quantum 796! numbers. If {\tt dosmsum} is set to {\tt .true.} then the partial DOS is 797! summed over $m$, and thus depends only on $l$. 798! 799! \block{dosssum}{ 800! {\tt dosssum} & {\tt .true.} if the partial DOS is to be summed over spin & 801! logical & {\tt .false.}} 802! By default, the partial density of states for spin-polarised systems is spin 803! resolved. 804! 805! \block{dtimes}{ 806! {\tt dtimes} & time step used in time evolution run & real & $0.1$} 807! See also {\tt tstime}. 808! 809! \block{epsband}{ 810! {\tt epsband} & convergence tolerance for determining band energies & real & 811! $1\times 10^{-12}$} 812! APW and local-orbital linearisation energies are determined from the band 813! energies. This is done by first searching upwards in energy until the radial 814! wavefunction at the muffin-tin radius is zero. This is the energy at the top 815! of the band, denoted $E_{\rm t}$. A downward search is now performed from 816! $E_{\rm t}$ until the slope of the radial wavefunction at the muffin-tin 817! radius is zero. This energy, $E_{\rm b}$, is at the bottom of the band. The 818! band energy is taken as $(E_{\rm t}+E_{\rm b})/2$. If either $E_{\rm t}$ or 819! $E_{\rm b}$ is not found, then the band energy is set to the default value. 820! 821! \block{epschg}{ 822! {\tt epschg} & maximum allowed error in the calculated total charge beyond 823! which a warning message will be issued & real & $1\times 10^{-3}$} 824! 825! \block{epsengy}{ 826! {\tt epsengy} & convergence criterion for the total energy & real & 827! $1\times 10^{-4}$} 828! See {\tt epspot}. 829! 830! \block{epsforce}{ 831! {\tt epsforce} & convergence tolerance for the forces during a geometry 832! optimisation run & real & $2\times 10^{-3}$} 833! If the mean absolute value of the atomic forces is less than {\tt epsforce} 834! then the geometry optimisation run is ended. See also {\tt tasks} and 835! {\tt latvopt}. 836! 837! \block{epslat}{ 838! {\tt epslat } & vectors with lengths less than this are considered zero & 839! real & $10^{-6}$} 840! Sets the tolerance for determining if a vector or its components are zero. 841! This is to account for any numerical error in real or reciprocal space 842! vectors. 843! 844! \block{epsocc}{ 845! {\tt epsocc} & smallest occupancy for which a state will contribute to the 846! density & real & $1\times 10^{-8}$} 847! 848! \block{epspot}{ 849! {\tt epspot} & convergence criterion for the Kohn-Sham potential and field & 850! real & $1\times 10^{-6}$} 851! If the RMS change in the Kohn-Sham potential and magnetic field is smaller 852! than {\tt epspot} and the absolute change in the total energy is less than 853! {\tt epsengy}, then the self-consistent loop is considered converged 854! and exited. For geometry optimisation runs this results in the forces being 855! calculated, the atomic positions updated and the loop restarted. See also 856! {\tt epsengy} and {\tt maxscl}. 857! 858! \block{epsstress}{ 859! {\tt epsstress} & convergence tolerance for the stress tensor during a 860! geometry optimisation run with lattice vector relaxation & real & 861! $5\times 10^{-4}$} 862! See also {\tt epsforce} and {\tt latvopt}. 863! 864! \block{emaxelnes}{ 865! {\tt emaxelnes} & maximum allowed initial-state eigenvalue for ELNES 866! calculations & real & $-1.2$} 867! 868! \block{emaxrf}{ 869! {\tt emaxrf} & energy cut-off used when calculating Kohn-Sham response 870! functions & real & $10^6$} 871! A typical Kohn-Sham response function is of the form 872! \begin{align*} 873! \chi_s({\bf r},{\bf r}',\omega) 874! \equiv\frac{\delta\rho({\bf r},\omega)}{\delta v_s({\bf r}',\omega)} 875! =\frac{1}{N_k}\sum_{i{\bf k},j{\bf k}'}(f_{i{\bf k}}-f_{j{\bf k}'}) 876! \frac{\langle i{\bf k}|\hat{\rho}({\bf r})|j{\bf k}'\rangle 877! \langle j{\bf k}'|\hat{\rho}({\bf r}')|i{\bf k}\rangle} 878! {w+(\varepsilon_{i{\bf k}}-\varepsilon_{j{\bf k}'})+i\eta}, 879! \end{align*} 880! where $\hat{\rho}$ is the density operator; $N_k$ is the number of 881! $k$-points; $\varepsilon_{i{\bf k}}$ and $f_{i{\bf k}}$ are the eigenvalues 882! and occupation numbers, respectively. The variable {\tt emaxrf} is an energy 883! window which limits the summation over states in the formula above so that 884! $|\varepsilon_{i{\bf k}}-\varepsilon_{\rm Fermi}|<{\tt emaxrf}$. Reducing 885! this can result in a faster calculation at the expense of accuracy. 886! 887! \block{fracinr}{ 888! {\tt fracinr} & fraction of the muffin-tin radius up to which {\tt lmaxi} 889! is used as the angular momentum cut-off & real & $0.01$} 890! If {\tt fracinr} is negative then the fraction is determined from 891! $f=\sqrt{({\tt lmaxi}+1)^2/({\tt lmaxo}+1)^2}$ in order to 892! maintain a minimum density of points throughout the muffin-tin. See 893! {\tt lmaxi} and {\tt lmaxo}. 894! 895! \block{fsmtype}{ 896! {\tt fsmtype} & 0 for no fixed spin moment (FSM), 1 for total FSM, 2 for 897! local muffin-tin FSM, and 3 for both total and local FSM & integer & 0} 898! Set to 1, 2 or 3 for fixed spin moment calculations. To fix only the 899! direction and not the magnitude set to $-1$, $-2$ or $-3$. See also 900! {\tt momfix}, {\tt mommtfix}, {\tt taufsm} and {\tt spinpol}. 901! 902! \block{ftmtype}{ 903! {\tt ftmtype} & 1 to enable a fixed tensor moment (FTM) calculation, 904! 0 otherwise & integer & 0} 905! If {\tt ftmtype} is $-1$ then the symmetry corresponding to the tensor 906! moment is broken but no FTM calculation is performed. See {\it Phys. Rev. 907! Lett.} {\bf 103}, 107202 (2009) and also {\tt tmomfix}. 908! 909! \block{fxclrc}{ 910! {\tt fxclrc} & parameters for the dynamical long-range contribution (LRC) to 911! the TDDFT exchange-correlation kernel & real(2) & $(0.0,0.0)$} 912! These are the parameters $\alpha$ and $\beta$ for the kernel proposed in 913! {\it Phys. Rev. B} {\bf 72}, 125203 (2005), namely 914! $$ f_{xc}({\bf G},{\bf G}',{\bf q},\omega)=-\frac{\alpha+\beta\omega^2}{q^2} 915! \delta_{{\bf G},{\bf G}'}\delta_{{\bf G},{\bf 0}}. $$ 916! 917! \block{fxctype}{ 918! {\tt fxctype} & integer defining the type of exchange-correlation kernel 919! $f_{\rm xc}$ & integer & $-1$} 920! The acceptable values are: 921! \vskip 6pt 922! \begin{tabularx}{\textwidth}[h]{lX} 923! $-1$ & $f_{\rm xc}$ defined by {\tt xctype} \\ 924! 0,1 & RPA ($f_{\rm xc}=0$) \\ 925! 200 & Long-range contribution (LRC) kernel, S. Botti {\it et al.}, 926! {\it Phys. Rev. B} {\bf 72}, 125203 (2005); see {\tt fxclrc} \\ 927! 210 & `Bootstrap' kernel, S. Sharma, J. K. Dewhurst, A. Sanna and 928! E. K. U. Gross, {\it Phys. Rev. Lett.} {\bf 107}, 186401 (2011) \\ 929! 211 & Single iteration bootstrap 930! \end{tabularx} 931! 932! \block{gmaxrf}{ 933! {\tt gmaxrf} & maximum length of $|{\bf G}|$ for computing response 934! functions & real & $3.0$} 935! 936! \block{gmaxvr}{ 937! {\tt gmaxvr} & maximum length of $|{\bf G}|$ for expanding the interstitial 938! density and potential & real & $12.0$} 939! This variable has a lower bound which is enforced by the code as follows: 940! $$ {\rm gmaxvr}\rightarrow\max\,({\rm gmaxvr},2\times{\rm gkmax} 941! +{\rm epslat}) $$ 942! See {\tt rgkmax}. 943! 944! \block{hdbse}{ 945! {\tt hdbse} & {\tt .true.} if the direct term is to be included in the BSE 946! Hamiltonian & logical & {\tt .true.}} 947! 948! \block{highq}{ 949! {\tt highq} & {\tt .true.} if a high quality parameter set should be used & 950! logical & {\tt .false.}} 951! Setting this to {\tt .true.} results in some default parameters being 952! changed to ensure good convergence in most situations. These changes can be 953! overruled by subsequent blocks in the input file. See also {\tt vhighq}. 954! 955! \block{hmaxvr}{ 956! {\tt hmaxvr} & maximum length of ${\bf H}$-vectors & real & $6.0$} 957! The ${\bf H}$-vectors are used for calculating X-ray and magnetic structure 958! factors. They are also used in linear response phonon calculations for 959! expanding the density and potential in plane waves. See also {\tt gmaxvr}, 960! {\tt vhmat}, {\tt reduceh}, {\tt wsfac} and {\tt hkmax}. 961! 962! \block{hxbse}{ 963! {\tt hxbse} & {\tt .true.} if the exchange term is to be included in the BSE 964! Hamiltonian & {\tt .true.}} 965! 966! \block{hybrid}{ 967! {\tt hybrid} & {\tt .true} if a hybrid functional is to be used when running 968! a Hartree-Fock calculation & logical & {\tt .false}} 969! See also {\tt hybridc} and {\tt xctype}. 970! 971! \block{hybridc}{ 972! {\tt hybridc} & hybrid functional mixing coefficient & real & $1.0$} 973! 974! \block{intraband}{ 975! {\tt intraband} & {\tt .true.} if the intraband (Drude-like) contribution is 976! to be added to the dieletric tensor & logical & {\tt .false.}} 977! 978! \block{isgkmax}{ 979! {\tt isgkmax} & species for which the muffin-tin radius will be used for 980! calculating {\tt gkmax} & integer & $-1$} 981! The APW cut-off is determined from ${\tt gkmax}={\tt rgkmax}/R$. The 982! variable {\tt isgkmax} determines which muffin-tin radius is to be used for 983! $R$. These are the options: 984! \vskip 6pt 985! \begin{tabularx}{\textwidth}[h]{lX} 986! -4 & Use the largest radius \\ 987! -3 & Use the smallest radius \\ 988! -2 & Use the fixed value $R=2.0$ \\ 989! -1 & Use the average of the muffin-tin radii \\ 990! $n\ge 1$ & Use the radius of species $n$ 991! \end{tabularx} 992! 993! \block{kstlist}{ 994! {\tt kstlist(i)} & $i$th $k$-point and state pair & integer(2) & $(1,1)$} 995! This is a user-defined list of $k$-point and state index pairs which are 996! those used for plotting wavefunctions and writing ${\bf L}$, ${\bf S}$ and 997! ${\bf J}$ expectation values. Only the first pair is used by the 998! aforementioned tasks. The list should be terminated by a blank line. 999! 1000! \block{latvopt}{ 1001! {\tt latvopt} & type of lattice vector optimisation to be performed during 1002! structural relaxation & integer & 0} 1003! Optimisation of the lattice vectors will be performed with ${\tt task}=2,3$ 1004! when ${\tt latvopt}\ne 0$. When ${\tt latvopt}=1$ the lattice vector 1005! optimisation will be constrained only by symmetry. Optimisation over all 1006! symmetry-preserving strains except isotropic scaling is performed when 1007! ${\tt latvopt}=2$. If ${\tt latvopt}<0$ then the optimisation will be over 1008! strain number $|{\tt latvopt}|$. The list of symmetric strain tensors can be 1009! produced with ${\tt task}=430$. By default (${\tt latvopt}=0$) no lattice 1010! vector optimisation is performed during structural relaxation. See also 1011! {\tt tau0latv} and {\tt atpopt}. 1012! 1013! \block{lmaxapw}{ 1014! {\tt lmaxapw} & angular momentum cut-off for the APW functions & integer & 1015! $8$} 1016! 1017! \block{lmaxdos}{ 1018! {\tt lmaxdos} & angular momentum cut-off for the partial DOS plot & 1019! integer & $3$} 1020! 1021! \block{lmaxi}{ 1022! {\tt lmaxi} & angular momentum cut-off for the muffin-tin density and 1023! potential on the inner part of the muffin-tin & integer & 2} 1024! Close to the nucleus, the density and potential is almost spherical and 1025! therefore the spherical harmonic expansion can be truncated a low angular 1026! momentum. See also {\tt fracinr}. 1027! 1028! \block{lmaxo}{ 1029! {\tt lmaxo} & angular momentum cut-off for the muffin-tin density and 1030! potential & integer & 6} 1031! 1032! \block{lmirep}{ 1033! {\tt lmirep} & {\tt .true.} if the $Y_{lm}$ basis is to be transformed 1034! into the basis of irreducible representations of the site symmetries for 1035! DOS plotting & logical & {\tt .true.}} 1036! When lmirep is set to .true., the spherical harmonic basis is transformed 1037! into one in which the site symmetries are block diagonal. Band characters 1038! determined from the density matrix expressed in this basis correspond to 1039! irreducible representations, and allow the partial DOS to be resolved into 1040! physically relevant contributions, for example $e_g$ and $t_{2g}$. 1041! 1042! \block{lorbcnd}{ 1043! {\tt lorbcnd} & {\tt .true.} if conduction state local-orbitals are to be 1044! automatically added to the basis & logical & {\tt .false.}} 1045! Adding these higher energy local-orbitals can improve calculations which 1046! rely on accurate unoccupied states, such as the response function. See also 1047! {\tt lorbordc}. 1048! 1049! \block{lorbordc}{ 1050! {\tt lorbordc} & the order of the conduction state local-orbitals & 1051! integer & 2} 1052! See {\tt lorbcnd}. 1053! 1054! \block{lradstp}{ 1055! {\tt lradstp} & radial step length for determining coarse radial mesh & 1056! integer & 4} 1057! Some muffin-tin functions (such as the density) are calculated on a coarse 1058! radial mesh and then interpolated onto a fine mesh. This is done for the 1059! sake of efficiency. {\tt lradstp} defines the step size in going from the 1060! fine to the coarse radial mesh. If it is too large, loss of precision may 1061! occur. 1062! 1063! \block{maxitoep}{ 1064! {\tt maxitoep} & maximum number of iterations when solving the exact 1065! exchange integral equations & integer & 300} 1066! See {\tt tau0oep}. 1067! 1068! \block{maxscl}{ 1069! {\tt maxscl } & maximum number of self-consistent loops allowed & integer & 1070! 200} 1071! This determines after how many loops the self-consistent cycle will 1072! terminate if the convergence criterion is not met. If {\tt maxscl} is $1$ 1073! then the density and potential file, {\tt STATE.OUT}, will {\bf not} be 1074! written to disk at the end of the loop. See {\tt epspot}. 1075! 1076! \block{mixtype}{ 1077! {\tt mixtype } & type of mixing required for the potential & integer & 1} 1078! Currently implemented are: 1079! \vskip 6pt 1080! \begin{tabularx}{\textwidth}[h]{lX} 1081! 0 & Linear mixing \\ 1082! 1 & Adaptive linear mixing \\ 1083! 3 & Broyden mixing, {\it J. Phys. A: Math. Gen.} {\bf 17}, L317 (1984) 1084! \end{tabularx} 1085! 1086! \block{mixsdb}{ 1087! {\tt mixsdb} & subspace dimension for Broyden mixing & integer & 5} 1088! This is the number of mixing vectors which define the subspace in which the 1089! Hessian matrix is calculated. See {\tt mixtype} and {\tt broydpm}. 1090! 1091! \block{molecule}{ 1092! {\tt molecule} & {\tt .true.} if the system is an isolated molecule & 1093! logical & {\tt .false.}} 1094! If {\tt molecule} is {\tt .true.}, then the atomic positions, ${\bf a}$, 1095! given in the {\tt atoms} block are assumed to be in Cartesian coordinates. 1096! 1097! \block{momfix}{ 1098! {\tt momfix} & the desired total moment for a FSM calculation & 1099! real(3) & $(0.0,0.0,0.0)$} 1100! Note that all three components must be specified (even for collinear 1101! calculations). See {\tt fsmtype}, {\tt taufsm} and {\tt spinpol}. 1102! 1103! \block{mommtfix}{ 1104! {\tt is} & species number & integer & 0 \\ 1105! {\tt ia} & atom number & integer & 0 \\ 1106! {\tt mommtfix} & the desired muffin-tin moment for a FSM calculation & 1107! real(3) & $(0.0,0.0,0.0)$} 1108! The local muffin-tin moments are specified for a subset of atoms, with the 1109! list terminated with a blank line. Note that all three components must be 1110! specified (even for collinear calculations). See {\tt fsmtype}, {\tt taufsm} 1111! and {\tt spinpol}. 1112! 1113! \block{msmooth}{ 1114! {\tt msmooth} & amount of smoothing to be applied to the 1115! exchange-correlation potentials and magnetic field & integer & 0} 1116! Smoothing operations can be applied to the exchange-correlation potentials 1117! $v_{xc}$, $w_{xc}$ and the magnetic field ${\bf B}_{xc}$ in order to improve 1118! convergence. In the muffin-tin, this smoothing takes the form of $m$ 1119! successive three-point running averages applied to the radial component. In 1120! the interstitial region, the potential is first Fourier transformed to 1121! $G$-space, then a low-pass filter of the form $\exp[-2m(G/G_{\rm max})^8]$ 1122! is applied and the function is transformed back to real-space. 1123! 1124! \block{mstar}{ 1125! {\tt mstar} & value of the effective mass parameter used for adaptive 1126! determination of {\tt swidth} & real & $10.0$} 1127! See {\tt autoswidth}. 1128! 1129! \block{mustar}{ 1130! {\tt mustar} & Coulomb pseudopotential, $\mu^*$, used in the 1131! McMillan-Allen-Dynes equation & real & $0.15$} 1132! This is used when calculating the superconducting critical temperature with 1133! the formula {\it Phys. Rev. B 12, 905 (1975)} 1134! $$ T_c=\frac{\omega_{\rm log}}{1.2 k_B}\exp\left[\frac{-1.04(1+\lambda)} 1135! {\lambda-\mu^*(1+0.62\lambda)}\right], $$ 1136! where $\omega_{\rm log}$ is the logarithmic average frequency and $\lambda$ 1137! is the electron-phonon coupling constant. 1138! 1139! \block{ncbse}{ 1140! {\tt ncbse} & number of conduction states to be used for BSE calculations & 1141! integer & 3} 1142! See also {\tt nvbse}. 1143! 1144! \block{ndspem}{ 1145! {\tt ndspem} & the number of {\bf k}-vector displacements in each direction 1146! around {\tt vklem} when computing the numerical derivatives for the 1147! effective mass tensor & integer & 1} 1148! See {\tt deltaem} and {\tt vklem}. 1149! 1150! \block{nempty}{ 1151! {\tt nempty} & the number of empty states per atom and spin & real & $4.0$ } 1152! Defines the number of eigenstates beyond that required for charge 1153! neutrality. When running metals it is not known {\it a priori} how many 1154! states will be below the Fermi energy for each $k$-point. Setting 1155! {\tt nempty} greater than zero allows the additional states to act as a 1156! buffer in such cases. Furthermore, magnetic calculations use the 1157! first-variational eigenstates as a basis for setting up the 1158! second-variational Hamiltonian, and thus {\tt nempty} will determine the 1159! size of this basis set. Convergence with respect to this quantity should be 1160! checked. 1161! 1162! \block{ngridk}{ 1163! {\tt ngridk } & the $k$-point mesh sizes & integer(3) & $(1,1,1)$} 1164! The ${\bf k}$-vectors are generated using 1165! $$ {\bf k}=(\frac{i_1+v_1}{n_1},\frac{i_2+v_2}{n_2},\frac{i_3+v_3}{n_3}), $$ 1166! where $i_j$ runs from 0 to $n_j-1$ and $0\le v_j<1$ for $j=1,2,3$. The 1167! vector ${\bf v}$ is given by the variable {\tt vkloff}. See also 1168! {\tt reducek}. 1169! 1170! \block{ngridq}{ 1171! {\tt ngridq } & the phonon $q$-point mesh sizes & integer(3) & $(1,1,1)$} 1172! Same as {\tt ngridk}, except that this mesh is for the phonon $q$-points 1173! and other tasks. See also {\tt reduceq}. 1174! 1175! \block{nosource}{ 1176! {\tt nosource} & when set to {\tt .true.}, source fields are projected out 1177! of the exchange-correlation magnetic field & logical & {\tt .false.}} 1178! Experimental feature. 1179! 1180! \block{notes}{ 1181! {\tt notes(i)} & the $i$th line of the notes & string & -} 1182! This block allows users to add their own notes to the file {\tt INFO.OUT}. 1183! The block should be terminated with a blank line, and no line should exceed 1184! 80 characters. 1185! 1186! \block{npmae}{ 1187! {\tt npmae } & number or distribution of directions for MAE calculations & 1188! integer & $-1$} 1189! Automatic determination of the magnetic anisotropy energy (MAE) requires 1190! that the total energy is determined for a set of directions of the total 1191! magnetic moment. This variable controls the number or distribution of these 1192! directions. The convention is: 1193! \vskip 6pt 1194! \begin{tabularx}{\textwidth}[h]{lX} 1195! $-4,-3,-2,-1$ & Cardinal directions given by the primitive translation 1196! vectors $n_1{\bf A}_1+n_2{\bf A}_2+n_3{\bf A}_3$, where 1197! $1\le n_i\le|{\tt npmae}|$ \\ 1198! 2 & Cartesian $x$ and $z$ directions \\ 1199! 3 & Cartesian $x$, $y$ and $z$ directions \\ 1200! $4,5,\ldots$ & Even distribution of {\tt npmae} directions 1201! \end{tabularx} 1202! 1203! \block{ntemp}{ 1204! {\tt ntemp} & number of temperature steps & integer & 40} 1205! This is the number of temperature steps to be used in the Eliashberg gap 1206! and thermodynamic properties calculations. 1207! 1208! \block{nvbse}{ 1209! {\tt nvbse} & number of valence states to be used for BSE calculations & 1210! integer & 2} 1211! See also {\tt ncbse}. 1212! 1213! \block{nwrite}{ 1214! {\tt nwrite} & number of self-consistent loops after which {\tt STATE.OUT} 1215! is to be written & integer & 0} 1216! Normally, the density and potentials are written to the file {\tt STATE.OUT} 1217! only after completion of the self-consistent loop. By setting {\tt nwrite} 1218! to a positive integer the file will instead be written every {\tt nwrite} 1219! loops. 1220! 1221! \block{nxoapwlo}{ 1222! {\tt nxoapwlo} & extra order of radial functions to be added to the existing 1223! APW and local-orbital set & integer & 0} 1224! Setting this variable will result in the APWs and local-orbitals for all 1225! species becoming higher order with corresponding increase in derivative 1226! matching at the muffin-tin surface. For example, setting {\tt nxoapwlo}=1 1227! turns all APWs into LAPWs. 1228! 1229! \block{optcomp}{ 1230! {\tt optcomp} & the components of the first- or second-order optical tensor 1231! to be calculated & integer(3) & $(1,1,1)$} 1232! This selects which components of the optical tensor you would like to plot. 1233! Only the first two are used for the first-order tensor. Several components 1234! can be listed one after the other with a blank line terminating the list. 1235! 1236! \block{phwrite}{ 1237! {\tt nphwrt} & number of $q$-points for which phonon modes are to be found & 1238! integer & 1 \\ 1239! \hline 1240! {\tt vqlwrt(i)} & the $i$th $q$-point in lattice coordinates & real(3) & 1241! $(0.0,0.0,0.0)$} 1242! This is used in conjunction with {\tt task}=230. The code will write the 1243! phonon frequencies and eigenvectors to the file {\tt PHONON.OUT} for all the 1244! $q$-points in the list. The $q$-points can be anywhere in the Brillouin zone 1245! and do not have to lie on the mesh defined by {\tt ngridq}. Obviously, all 1246! the dynamical matrices have to be computed first using {\tt task}=200. 1247! 1248! \block{plot1d}{ 1249! {\tt nvp1d} & number of vertices & integer & 2 \\ 1250! {\tt npp1d} & number of plotting points & integer & 200 \\ 1251! \hline 1252! {\tt vvlp1d(i)} & lattice coordinates for vertex $i$ & real(3) & 1253! $(0.0,0.0,0.0)\rightarrow(1.0,1.0,1.0)$} 1254! Defines the path in either real or reciprocal space along which the 1D plot 1255! is to be produced. The user should provide {\tt nvp1d} vertices in lattice 1256! coordinates. 1257! 1258! \block{plot2d}{ 1259! {\tt vclp2d(0)} & zeroth corner (origin) & real(3) & $(0.0,0.0,0.0)$ \\ 1260! \hline 1261! {\tt vclp2d(1)} & first corner & real(3) & $(1.0,0.0,0.0)$ \\ 1262! \hline 1263! {\tt vclp2d(2)} & second corner & real(3) & $(0.0,1.0,0.0)$ \\ 1264! \hline 1265! {\tt np2d} & number of plotting points in both directions & integer(2) & 1266! $(40,40)$} 1267! Defines the corners of a parallelogram and the grid size used for producing 1268! 2D plots. 1269! 1270! \block{plot3d}{ 1271! {\tt vclp3d(0)} & zeroth corner (origin) & real(3) & $(0.0,0.0,0.0)$ \\ 1272! \hline 1273! {\tt vclp3d(1)} & first corner & real(3) & $(1.0,0.0,0.0)$ \\ 1274! \hline 1275! {\tt vclp3d(2)} & second corner & real(3) & $(0.0,1.0,0.0)$ \\ 1276! \hline 1277! {\tt vclp3d(3)} & third corner & real(3) & $(0.0,0.0,1.0)$ \\ 1278! \hline 1279! {\tt np3d} & number of plotting points each direction & integer(3) & 1280! $(20,20,20)$} 1281! Defines the corners of a box and the grid size used for producing 3D plots. 1282! 1283! \block{primcell}{ 1284! {\tt primcell} & {\tt .true.} if the primitive unit cell should be found 1285! & logical & {\tt .false.}} 1286! Allows the primitive unit cell to be determined automatically from the 1287! conventional cell. This is done by searching for lattice vectors among all 1288! those which connect atomic sites, and using the three shortest which produce 1289! a unit cell with non-zero volume. 1290! 1291! \block{pulse}{ 1292! {\tt n} & number of pulses & integer & - \\ 1293! \hline 1294! {\tt a0(i)} & polarisation vector (including amplitude) & real(3) & - \\ 1295! {\tt w(i)} & frequency & real & - \\ 1296! {\tt phi(i)} & phase in degrees & real & - \\ 1297! {\tt rc(i)} & chirp rate & real & - \\ 1298! {\tt t0(i)} & peak time & real & - \\ 1299! {\tt d(i)} & full-width at half-maximum & real & -} 1300! Parameters used to generate a time-dependent vector potential ${\bf A}(t)$ 1301! representing a laser pulse. The total vector potential is the sum of 1302! individual pulses and is given by the formula 1303! $$ {\bf A}(t)=\sum_{i=1}^n {\bf A}_0^i\exp 1304! \left[-(t-t_0^i)^2/2\sigma_i^2\right] 1305! \sin\left[w_i(t-t_0^i)+\phi_i+r_{\rm c}^i t^2/2\right], $$ 1306! where $\sigma=d/2\sqrt{2\ln 2}$. See also {\tt ramp}. 1307! 1308! \block{radkpt}{ 1309! {\tt radkpt } & radius of sphere used to determine $k$-point density & 1310! real & $40.0$} 1311! Used for the automatic determination of the $k$-point mesh. If {\tt autokpt} 1312! is set to {\tt .true.} then the mesh sizes will be determined by 1313! $n_i=R_k|{\bf B}_i|+1$, where ${\bf B}_i$ are the primitive reciprocal 1314! lattice vectors. 1315! 1316! \block{ramp}{ 1317! {\tt n} & number of ramps & integer & - \\ 1318! \hline 1319! {\tt a0(i)} & polarisation vector (including amplitude) & real(3) & - \\ 1320! {\tt t0(i)} & ramp start time & real & - \\ 1321! {\tt c1(i)} & linear coefficient of ${\bf A}(t)$ & real & - \\ 1322! {\tt c2(i)} & quadratic coefficient & real & -} 1323! Parameters used to generate a time-dependent vector potential ${\bf A}(t)$ 1324! representing a constant or linearly increasing electric field 1325! ${\bf E}(t)=-\partial{\bf A}(t)/\partial t$. The vector potential is given 1326! by 1327! $$ {\bf A}(t)=\sum_{i=1}^n {\bf A}_0^i 1328! \left[c_1(t-t_0)+c_2(t-t_0)^2\right]\Theta(t-t_0). $$ 1329! 1330! \block{readadu}{ 1331! {\tt readadu} & set to {\tt .true.} if the interpolation constant for 1332! DFT+$U$ should be read from file rather than calculated & logical & 1333! {\tt .false.}} 1334! When {\tt dftu}=3, the DFT+$U$ energy and potential are interpolated 1335! between FLL and AFM. The interpolation constant, $\alpha$, is normally 1336! calculated from the density matrix, but can also be read in from the file 1337! {\tt ALPHADU.OUT}. This allows the user to fix $\alpha$, but is also 1338! necessary when calculating forces, since the contribution of the potential 1339! of the variation of $\alpha$ with respect to the density matrix is not 1340! computed. See {\tt dft+u}. 1341! 1342! \block{reducebf}{ 1343! {\tt reducebf} & reduction factor for the external magnetic fields & real & 1344! $1.0$} 1345! After each self-consistent loop, the external magnetic fields are multiplied 1346! with {\tt reducebf}. This allows for a large external magnetic field at the 1347! start of the self-consistent loop to break spin symmetry, while at the end 1348! of the loop the field will be effectively zero, i.e. infinitesimal. See 1349! {\tt bfieldc} and {\tt atoms}. 1350! 1351! \block{reduceh}{ 1352! {\tt reduceh} & set to {\tt .true.} if the reciprocal ${\bf H}$-vectors 1353! should be reduced by the symmorphic crystal symmetries & logical & .true.} 1354! See {\tt hmaxvr} and {\tt vmat}. 1355! 1356! \block{reducek}{ 1357! {\tt reducek} & type of reduction of the $k$-point set & integer & 1} 1358! Types of reduction are defined by the symmetry group used: 1359! \vskip 6pt 1360! \begin{tabularx}{\textwidth}[h]{lX} 1361! 0 & no reduction \\ 1362! 1 & reduce with full crystal symmetry group (including non-symmorphic 1363! symmetries) \\ 1364! 2 & reduce with symmorphic symmetries only 1365! \end{tabularx} 1366! \vskip 6pt 1367! See also {\tt ngridk} and {\tt vkloff}. 1368! 1369! \block{reduceq}{ 1370! {\tt reduceq} & type of reduction of the $q$-point set & integer & 1} 1371! See {\tt reducek} and {\tt ngridq}. 1372! 1373! \block{rgkmax}{ 1374! {\tt rgkmax} & $R^{\rm MT}_{\rm min}\times\max\{|{\bf G}+{\bf k}|\}$ & 1375! real & $7.0$} 1376! This sets the maximum length for the ${\bf G}+{\bf k}$ vectors, defined as 1377! {\tt rgkmax} divided by the average muffin-tin radius. See {\tt isgkmax}. 1378! 1379! \block{rotavec}{ 1380! {\tt axang} & axis-angle representation of lattice vector rotation & 1381! real(4) & $(0.0,0.0,0.0,0.0)$} 1382! This determines the rotation matrix which is applied to the lattice vectors 1383! prior to any calculation. The first three components specify the axis and 1384! the last component is the angle in degrees. The `right-hand rule' convention 1385! is followed. 1386! 1387! \block{scale}{ 1388! {\tt scale } & lattice vector scaling factor & real & $1.0$} 1389! Scaling factor for all three lattice vectors. Applied in conjunction with 1390! {\tt scale1}, {\tt scale2} and {\tt scale3}. 1391! 1392! \block{scale1/2/3}{ 1393! {\tt scale1/2/3 } & separate scaling factors for each lattice vector & 1394! real & $1.0$} 1395! 1396! \block{scissor}{ 1397! {\tt scissor} & the scissor correction & real & $0.0$} 1398! This is the scissor shift applied to states above the Fermi energy 1399! {\it Phys. Rev. B} {\bf 43}, 4187 (1991). Affects optics calculations only. 1400! 1401! \block{scrpath}{ 1402! {\tt scrpath} & scratch space path & string & null} 1403! This is the scratch space path where the eigenvector files {\tt EVALFV.OUT} 1404! and {\tt EVALSV.OUT} will be written. If the run directory is accessed via a 1405! network then {\tt scrpath} can be set to a directory on the local disk, for 1406! example {\tt /tmp/}. Note that the forward slash {\tt /} at the end of the 1407! path must be included. 1408! 1409! \block{socscf}{ 1410! {\tt socscf} & scaling factor for the spin-orbit coupling term in the 1411! Hamiltonian & real & $1.0$} 1412! This can be used to enhance the effect of spin-orbit coupling in order to 1413! accurately determine the magnetic anisotropy energy (MAE). 1414! 1415! \block{spincore}{ 1416! {\tt spincore} & set to {\tt .true.} if the core should be spin-polarised 1417! & logical & {\tt .false.}} 1418! 1419! \block{spinorb}{ 1420! {\tt spinorb} & set to {\tt .true.} if a spin-orbit coupling is required 1421! & logical & {\tt .false.}} 1422! If {\tt spinorb} is {\tt .true.}, then a $\boldsymbol\sigma\cdot{\bf L}$ 1423! term is added to the second-variational Hamiltonian. See {\tt spinpol}. 1424! 1425! \block{spinpol}{ 1426! {\tt spinpol} & set to {\tt .true.} if a spin-polarised calculation is 1427! required & logical & {\tt .false.}} 1428! If {\tt spinpol} is {\tt .true.}, then the spin-polarised Hamiltonian is 1429! solved as a second-variational step using two-component spinors in the 1430! Kohn-Sham magnetic field. The first variational scalar wavefunctions are 1431! used as a basis for setting this Hamiltonian. 1432! 1433! \block{spinsprl}{ 1434! {\tt spinsprl} & set to {\tt .true.} if a spin-spiral calculation is 1435! required & logical & {\tt .false.}} 1436! Experimental feature for the calculation of spin-spiral states. See 1437! {\tt vqlss} for details. 1438! 1439! \block{sppath}{ 1440! {\tt sppath} & path where the species files can be found & string & null} 1441! Note that the forward slash {\tt /} at the end of the path must be included. 1442! 1443! \block{ssdph}{ 1444! {\tt ssdph} & set to {\tt .true.} if a complex de-phasing factor is to be 1445! used in spin-spiral calculations & logical & {\tt .true.}} 1446! If this is {\tt .true.} then spin-spiral wavefunctions in each muffin-tin at 1447! position ${\bf r}_{\alpha}$ are de-phased by the matrix 1448! $$ \begin{pmatrix} e^{-i{\bf q}\cdot{\bf r}_{\alpha}/2} & 0 \\ 1449! 0 & e^{i{\bf q}\cdot{\bf r}_{\alpha}/2} \end{pmatrix}. $$ 1450! In simple situations, this has the advantage of producing magnon dynamical 1451! matrices which are already in diagonal form. This option should be used with 1452! care, and a full understanding of the spin-spiral configuration is required. 1453! See {\tt spinsprl}. 1454! 1455! \block{stype}{ 1456! {\tt stype} & integer defining the type of smearing to be used & integer & 1457! $3$} 1458! A smooth approximation to the Dirac delta function is needed to compute the 1459! occupancies of the Kohn-Sham states. The variable {\tt swidth} determines 1460! the width of the approximate delta function. Currently implemented are 1461! \vskip 6pt 1462! \begin{tabularx}{\textwidth}[h]{lX} 1463! 0 & Gaussian \\ 1464! 1 & Methfessel-Paxton order 1, Phys. Rev. B {\bf 40}, 3616 (1989) \\ 1465! 2 & Methfessel-Paxton order 2 \\ 1466! 3 & Fermi-Dirac 1467! \end{tabularx} 1468! \vskip 6pt 1469! See also {\tt autoswidth}, {\tt swidth} and {\tt tempk}. 1470! 1471! \block{swidth}{ 1472! {\tt swidth} & width of the smooth approximation to the Dirac delta 1473! function & real & $0.001$} 1474! See {\tt stype} for details and the variable {\tt tempk}. 1475! 1476! \newpage 1477! \block{tasks}{ 1478! {\tt task(i) } & the $i$th task & integer & $-1$} 1479! A list of tasks for the code to perform sequentially. The list should be 1480! terminated with a blank line. Each task has an associated integer as 1481! follows: 1482! \vskip 6pt 1483! \begin{tabularx}{\textwidth}[h]{lX} 1484! -1 & Write out the version number of the code. \\ 1485! 0 & Ground state run starting from the atomic densities. \\ 1486! 1 & Resumption of ground-state run using density in {\tt STATE.OUT}. \\ 1487! 2 & Geometry optimisation run starting from the atomic densities, with 1488! atomic positions written to {\tt GEOMETRY.OUT}. \\ 1489! 3 & Resumption of geometry optimisation run using density in {\tt STATE.OUT} 1490! but with positions from {\tt elk.in}. \\ 1491! 5 & Ground state Hartree-Fock run. \\ 1492! 10 & Total, partial and interstitial density of states (DOS). \\ 1493! 14 & Plots the smooth Dirac delta and Heaviside step functions used by the 1494! code to calculate occupancies. \\ 1495! 15 & Output ${\bf L}$, ${\bf S}$ and ${\bf J}$ total expectation values. \\ 1496! 16 & Output ${\bf L}$, ${\bf S}$ and ${\bf J}$ expectation values for each 1497! $k$-point and state in {\tt kstlist}. \\ 1498! 20 & Band structure plot. \\ 1499! 21 & Band structure plot which includes total and angular momentum 1500! characters for every atom. \\ 1501! 22 & Band structure plot which includes $(l,m)$ character for every atom. \\ 1502! 23 & Band structure plot which includes spin character for every atom. \\ 1503! 25 & Compute the effective mass tensor at the $k$-point given by 1504! {\tt vklem}. \\ 1505! 31, 32, 33 & 1/2/3D charge density plot. \\ 1506! 41, 42, 43 & 1/2/3D exchange-correlation and Coulomb potential plots. \\ 1507! 51, 52, 53 & 1/2/3D electron localisation function (ELF) plot. \\ 1508! 61, 62, 63 & 1/2/3D wavefunction plot: 1509! $\left|\Psi_{i{\bf k}}({\bf r})\right|^2$. \\ 1510! 65 & Write the core wavefunctions to file for plotting. \\ 1511! 71, 72, 73 & 1/2/3D plot of magnetisation vector field, 1512! ${\bf m}({\bf r})$. \\ 1513! 81, 82, 83 & 1/2/3D plot of exchange-correlation magnetic vector field, 1514! ${\bf B}_{\rm xc}({\bf r})$. \\ 1515! 91, 92, 93 & 1/2/3D plot of $\nabla\cdot{\bf B}_{\rm xc}({\bf r})$. \\ 1516! 100 & 3D Fermi surface plot using the scalar product 1517! $p({\bf k})=\Pi_i(\epsilon_{i{\bf k}}-\epsilon_{\rm F})$. \\ 1518! 101 & 3D Fermi surface plot using separate bands (minus the Fermi 1519! energy). \\ 1520! 102 & 3D Fermi surface which can be plotted with XCrysDen. \\ 1521! 105 & 3D nesting function plot. \\ 1522! 110 & Calculation of M\"{o}ssbauer contact charge densities and magnetic 1523! fields at the nuclear sites. \\ 1524! 115 & Calculation of the electric field gradient (EFG) at the nuclear 1525! sites. \\ 1526! 120 & Output of the momentum matrix elements 1527! $\langle\Psi_{i{\bf k}}|-i\nabla|\Psi_{j{\bf k}}\rangle$. \\ 1528! 121 & Linear optical dielectric response tensor calculated within the random 1529! phase approximation (RPA) and in the $q\rightarrow 0$ limit, with no 1530! microscopic contributions. \\ 1531! 122 & Magneto optical Kerr effect (MOKE) angle. \\ 1532! 125 & Non-linear optical second harmonic generation. 1533! \end{tabularx} 1534! 1535! \begin{tabularx}{\textwidth}[h]{lX} 1536! 130 & Output matrix elements of the type 1537! $\langle\Psi_{i{\bf k+q}}|\exp[i({\bf G+q})\cdot{\bf r}]| 1538! \Psi_{j{\bf k}}\rangle$. \\ 1539! 135 & Output all wavefunctions expanded in the plane wave basis up to a 1540! cut-off defined by {\tt rgkmax}. \\ 1541! 140 & Energy loss near edge structure (ELNES). \\ 1542! 141, 142, 143 & 1/2/3D plot of the electric field 1543! ${\bf E}({\bf r})\equiv\nabla V_{\rm C}({\bf r})$. \\ 1544! 151, 152, 153 & 1/2/3D plot of 1545! ${\bf m}({\bf r})\times{\bf B}_{\rm xc}({\bf r})$. \\ 1546! 162 & Scanning-tunneling microscopy (STM) image. \\ 1547! 180 & Generate the RPA inverse dielectric function with local contributions 1548! and write it to file. \\ 1549! 185 & Write the Bethe-Salpeter equation (BSE) Hamiltonian to file. \\ 1550! 186 & Diagonalise the BSE Hamiltonian and write the eigenvectors and 1551! eigenvalues to file. \\ 1552! 187 & Output the BSE dielectric response function. \\ 1553! 190 & Write the atomic geometry to file for plotting with XCrySDen and 1554! V\_Sim. \\ 1555! 195 & Calculation of X-ray density structure factors. \\ 1556! 196 & Calculation of magnetic structure factors. \\ 1557! 200 & Calculation of phonon dynamical matrices on a $q$-point set defined by 1558! {\tt ngridq} using the supercell method. \\ 1559! 202 & Phonon dry run: just produce a set of empty DYN files. \\ 1560! 205 & Calculation of phonon dynamical matrices using density functional 1561! perturbation theory (DFPT). \\ 1562! 210 & Phonon density of states. \\ 1563! 220 & Phonon dispersion plot. \\ 1564! 230 & Phonon frequencies and eigenvectors for an arbitrary $q$-point. \\ 1565! 240, 241 & Generate the ${\bf q}$-dependent phonon linewidths and 1566! electron-phonon coupling constants and write them to file. \\ 1567! 245 & Phonon linewidths plot. \\ 1568! 250 & Eliashberg function $\alpha^2F(\omega)$, electron-phonon coupling 1569! constant $\lambda$, and the McMillan-Allen-Dynes critical temperature 1570! $T_c$. \\ 1571! 300 & Reduced density matrix functional theory (RDMFT) calculation. \\ 1572! 320 & Time-dependent density functional theory (TDDFT) calculation of the 1573! dielectric response function including microscopic contributions. \\ 1574! 330, 331 & TDDFT calculation of the spin-polarised response function for 1575! arbitrary ${\bf q}$-vectors. Task 331 writes the entire response function 1576! $\overleftrightarrow{\chi}({\bf G},{\bf G}',q,\omega)$ to file. \\ 1577! 400 & Calculation of tensor moments and corresponding DFT+$U$ Hartree-Fock 1578! energy contributions. \\ 1579! 450 & Generates a laser pulse in the form of a time-dependent vector 1580! potential ${\bf A}(t)$ and writes it to AFIELDT.OUT. \\ 1581! 460 & Time evolution run using TDDFT under the influence of ${\bf A}(t)$. 1582! \end{tabularx} 1583! 1584! \block{tau0atp}{ 1585! {\tt tau0atp} & the step size to be used for atomic position optimisation & 1586! real & $0.25$} 1587! The position of atom $\alpha$ is updated on step $m$ of a geometry 1588! optimisation run using 1589! $$ {\bf r}_{\alpha}^{m+1}={\bf r}_{\alpha}^m+\tau_{\alpha}^m 1590! \left({\bf F}_{\alpha}^m+{\bf F}_{\alpha}^{m-1}\right), $$ 1591! where $\tau_{\alpha}$ is set to {\tt tau0atp} for $m=0$, and incremented by 1592! the same amount if the atom is moving in the same direction between steps. 1593! If the direction changes then $\tau_{\alpha}$ is reset to {\tt tau0atp}. 1594! 1595! \block{tau0latv}{ 1596! {\tt tau0latv} & the step size to be used for lattice vector optimisation & 1597! real & $0.25$} 1598! This parameter is used for lattice vector optimisation in a procedure 1599! identical to that for atomic position optimisation. See {\tt tau0atp} and 1600! {\tt latvopt}. 1601! 1602! \block{tauoep}{ 1603! {\tt tauoep} & step length and starting fraction for the OEP iterative 1604! solver & real(2) & $(0.05, 0.95)$} 1605! The optimised effective potential is determined using an interative method 1606! [Phys. Rev. Lett. 98, 196405 (2007)]. This variable sets the step length 1607! described in the article. The second number defines the fraction of the 1608! existing potential to be used as the starting point of the iterative 1609! procedure. These parameters are also used for inverting the Kohn-Sham 1610! equations as required by the self-consistent density $GW$ method. See 1611! {\tt maxitoep}. 1612! 1613! \block{taufsm}{ 1614! {\tt taufsm} & the step size to be used when finding the effective magnetic 1615! field in fixed spin moment calculations & real & $0.01$} 1616! An effective magnetic field, ${\bf B}_{\rm FSM}$, is required for fixing the 1617! spin moment to a given value, ${\bf M}_{\rm FSM}$. This is found by adding a 1618! vector to the field which is proportional to the difference between the 1619! moment calculated in the $i$th self-consistent loop and the required moment: 1620! $$ {\bf B}_{\rm FSM}^{i+1}={\bf B}_{\rm FSM}^i+\lambda\left({\bf M}^i 1621! -{\bf M}_{\rm FSM}\right), $$ 1622! where $\lambda$ is proportional to {\tt taufsm}. See also {\tt fsmtype}, 1623! {\tt momfix} and {\tt spinpol}. 1624! 1625! \block{tempk}{ 1626! {\tt tempk} & temperature $T$ of the electronic system in kelvin & real & -} 1627! Assigning a value to this variable sets {\tt stype} to 3 (Fermi-Dirac) and 1628! the smearing width to $k_{\rm B}T$. 1629! 1630! \block{tforce}{ 1631! {\tt tforce} & set to {\tt .true.} if the force should be calculated at the 1632! end of the self-consistent cycle & logical & {\tt .false.}} 1633! This variable is automatically set to {\tt .true.} when performing geometry 1634! optimisation. 1635! 1636! \block{tefvit}{ 1637! {\tt tefvit} & set to {\tt .true.} if the first-variational eigenvalue 1638! equation should be solved iteratively & logical & {\tt .false.}} 1639! 1640! \block{tefvr}{ 1641! {\tt tefvr} & set to {\tt .true.} if a real symmetric eigenvalue solver 1642! should be used for crystals which have inversion symmetry & logical & 1643! {\tt .true.}} 1644! For crystals with inversion symmetry, the first-variational Hamiltonian and 1645! overlap matrices can be made real by using appropriate transformations. In 1646! this case, a real symmetric (instead of complex Hermitian) eigenvalue solver 1647! can be used. This makes the calculation about three times faster. 1648! 1649! \block{tmomfix}{ 1650! {\tt ntmfix} & number of tensor moments (TM) to be fixed & integer & 0 \\ 1651! \hline 1652! {\tt is(i)} & species number for entry $i$ & integer & - \\ 1653! {\tt ia(i)} & atom number & integer & - \\ 1654! {\tt (l, n)(i)} & $l$ and $n$ indices of TM & integer & - \\ 1655! \hline 1656! {\tt (k, p, x, y)(i)} or & & & \\ 1657! {\tt (k, p, r, t)(i)} & indices for the 2-index or 3-index TM, 1658! respectively & integer & - \\ 1659! \hline 1660! {\tt z(i)} & complex TM value & complex & - \\ 1661! \hline 1662! {\tt p(i)} & parity of spatial rotation & integer & - \\ 1663! {\tt aspl(i)} & Euler angles of spatial rotation & real(3) & - \\ 1664! {\tt aspn(i)} & Euler angles of spin rotation & real(3) & - } 1665! This block sets up the fixed tensor moment (FTM). There should be as many 1666! TM entries as {\tt ntmfix}. See {\it Phys. Rev. Lett.} {\bf 103}, 107202 1667! (2009) for the tensor moment indexing convention. This is a highly 1668! complicated feature of the code, and should only be attempted with a full 1669! understanding of tensor moments. 1670! 1671! \block{tmwrite}{ 1672! {\tt tmwrite} & set to {\tt .true.} if the tensor moments and the 1673! corresponding decomposition of DFT+$U$ energy should be calculated 1674! at every loop of the self-consistent cycle & logical & {\tt .false.}} 1675! This variable is useful to check the convergence of the tensor moments in 1676! DFT+$U$ caculations. Alternatively, with {\tt task} equal to 400, one can 1677! calculate the tensor moments and corresponding DFT+$U$ energy contributions 1678! from a given density matrix and set of Slater parameters at the end of the 1679! self-consistent cycle. 1680! 1681! \block{tsediag}{ 1682! {\tt tsediag} & set to {\tt .true.} if the self-energy matrix should be 1683! treated as diagonal & logical & {\tt .true.}} 1684! When this variable is {\tt .true.}, the self-energy used in a $GW$ 1685! calculation $\Sigma_{ij}({\bf k},\omega)$ is taken to be diagonal in the 1686! Kohn-Sham state indices $i$ and $j$. When {\tt tsediag} is {\tt .false.}, 1687! the entire matrix is used. See also {\tt twdiag}. 1688! 1689! \block{tshift}{ 1690! {\tt tshift} & set to {\tt .true.} if the crystal can be shifted so that the 1691! atom closest to the origin is exactly at the origin & 1692! logical & {\tt .true.}} 1693! 1694! \block{tstime}{ 1695! {\tt tstime} & total simulation time of time evolution run & real & 1696! $1000.0$} 1697! See also {\tt dtimes}. 1698! 1699! \block{twdiag}{ 1700! {\tt twdiag} & set to {\tt .true.} if the screened interaction matrix should 1701! be treated as diagonal & logical & {\tt .false.}} 1702! When this variable is {\tt .true.}, the screened interaction used in a $GW$ 1703! calculation $W({\bf G},{\bf G}',{\bf q},\omega)$ is taken to be diagonal in 1704! the plane wave indices ${\bf G}$ and ${\bf G}'$. See also {\tt tsediag}. 1705! 1706! \block{vhmat}{ 1707! {\tt vhmat(1)} & matrix row 1 & real(3) & $(1.0,0.0,0.0)$ \\ 1708! \hline 1709! {\tt vhmat(2)} & matrix row 2 & real(3) & $(0.0,1.0,0.0)$ \\ 1710! \hline 1711! {\tt vhmat(3)} & matrix row 3 & real(3) & $(0.0,0.0,1.0)$} 1712! This is the transformation matrix $M$ applied to every vector $\bf H$ in the 1713! structure factor output files {\tt SFACRHO.OUT} and {\tt SFACMAG.OUT}. It is 1714! stored in the usual row-column setting and applied directly as 1715! ${\bf H}'=M{\bf H}$ to every vector but {\em only} when writing the output 1716! files. See also {\tt hmaxvr} and {\tt reduceh}. 1717! 1718! \block{vhighq}{ 1719! {\tt vhighq} & {\tt .true.} if a very high quality parameter set should be 1720! used & logical & {\tt .false.}} 1721! Setting this to {\tt .true.} results in some default parameters being 1722! changed to ensure excellent convergence in most situations. See also 1723! {\tt highq}. 1724! 1725! \block{vklem}{ 1726! {\tt vklem} & the $k$-point in lattice coordinates at which to compute the 1727! effective mass tensors & real(3) & $(0.0,0.0,0.0)$} 1728! See {\tt deltaem} and {\tt ndspem}. 1729! 1730! \block{vkloff}{ 1731! {\tt vkloff } & the $k$-point offset vector in lattice coordinates & 1732! real(3) & $(0.0,0.0,0.0)$} 1733! See {\tt ngridk}. 1734! 1735! \block{vqlss}{ 1736! {\tt vqlss} & the ${\bf q}$-vector of the spin-spiral state in lattice 1737! coordinates & real(3) & $(0.0,0.0,0.0)$} 1738! Spin-spirals arise from spinor states assumed to be of the form 1739! $$ \Psi^{\bf q}_{\bf k}({\bf r})= 1740! \left( \begin{array}{c} 1741! U^{{\bf q}\uparrow}_{\bf k}({\bf r})e^{i({\bf k+q/2})\cdot{\bf r}} \\ 1742! U^{{\bf q}\downarrow}_{\bf k}({\bf r})e^{i({\bf k-q/2})\cdot{\bf r}} \\ 1743! \end{array} \right). $$ 1744! These are determined using a second-variational approach, and give rise to a 1745! magnetisation density of the form 1746! $$ {\bf m}^{\bf q}({\bf r})=(m_x({\bf r})\cos({\bf q \cdot r}), 1747! m_y({\bf r})\sin({\bf q \cdot r}),m_z({\bf r})), $$ 1748! where $m_x$, $m_y$ and $m_z$ are lattice periodic. See also {\tt spinsprl}. 1749! 1750! \block{wmaxgw}{ 1751! {\tt wmaxgw} & maximum Matsubara frequency for $GW$ calculations & real & 1752! $-5.0$} 1753! This defines the cut-off of the Matsubara frequencies on the imaginary 1754! axis for calculating the $GW$ self-energy and solving the Dyson equation. 1755! If this number is negative then the cut-off is taken to be 1756! $|{\tt wmaxgw}|\times\Delta\epsilon$, where $\Delta\epsilon$ is the 1757! difference between the largest and smallest Kohn-Sham valence eigenvalues. 1758! 1759! \block{wplot}{ 1760! {\tt nwplot} & number of frequency/energy points in the DOS or optics plot & 1761! integer & $500$ \\ 1762! {\tt ngrkf} & fine $k$-point grid size used for integrating functions in the 1763! Brillouin zone & integer & $100$ \\ 1764! {\tt nswplot} & level of smoothing applied to DOS/optics output & integer & 1765! $1$ \\ 1766! \hline 1767! {\tt wplot} & frequency/energy window for the DOS or optics plot & real(2) & 1768! $(-0.5,0.5)$} 1769! DOS and optics plots require integrals of the kind 1770! $$ g(\omega_i)=\frac{\Omega}{(2\pi)^3}\int_{\rm BZ} f({\bf k}) 1771! \delta(\omega_i-e({\bf k}))d{\bf k}. $$ 1772! These are calculated by first interpolating the functions $e({\bf k})$ and 1773! $f({\bf k})$ with the trilinear method on a much finer mesh whose size is 1774! determined by {\tt ngrkf}. Then the $\omega$-dependent histogram of the 1775! integrand is accumulated over the fine mesh. If the output function is noisy 1776! then either {\tt ngrkf} should be increased or {\tt nwplot} decreased. 1777! Alternatively, the output function can be artificially smoothed up to a 1778! level given by {\tt nswplot}. This is the number of successive 3-point 1779! averages to be applied to the function $g$. 1780! 1781! \block{wsfac}{ 1782! {\tt wsfac} & energy window to be used when calculating density or magnetic 1783! structure factors & real(2) & $(-10^6,10^6)$} 1784! Only those states with eigenvalues within this window will contribute to the 1785! density or magnetisation. See also {\tt hmaxvr} and {\tt vhmat}. 1786! 1787! \block{xctype}{ 1788! {\tt xctype} & integers defining the type of exchange-correlation functional 1789! to be used & integer(3) & $(3,0,0)$} 1790! Normally only the first value is used to define the functional type. The 1791! other value may be used for external libraries. Currently implemented are: 1792! \vskip 6pt 1793! \begin{tabularx}{\textwidth}[h]{lX} 1794! $-n$ & Exact-exchange optimised effective potential (EXX-OEP) method with 1795! correlation energy and potential given by functional number $n$ \\ 1796! 1 & No exchange-correlation funtional ($E_{\rm xc}\equiv 0$) \\ 1797! 2 & LDA, Perdew-Zunger/Ceperley-Alder, {\it Phys. Rev. B} {\bf 23}, 5048 1798! (1981) \\ 1799! 3 & LSDA, Perdew-Wang/Ceperley-Alder, {\it Phys. Rev. B} {\bf 45}, 13244 1800! (1992) \\ 1801! 4 & LDA, X-alpha approximation, J. C. Slater, {\it Phys. Rev.} {\bf 81}, 385 1802! (1951) \\ 1803! 5 & LSDA, von Barth-Hedin, {\it J. Phys. C} {\bf 5}, 1629 (1972) \\ 1804! 20 & GGA, Perdew-Burke-Ernzerhof, {\it Phys. Rev. Lett.} {\bf 77}, 3865 1805! (1996) \\ 1806! 21 & GGA, Revised PBE, Zhang-Yang, {\it Phys. Rev. Lett.} {\bf 80}, 890 1807! (1998) \\ 1808! 22 & GGA, PBEsol, Phys. Rev. Lett. 100, 136406 (2008) \\ 1809! 26 & GGA, Wu-Cohen exchange (WC06) with PBE correlation, {\it Phys. Rev. B} 1810! {\bf 73}, 235116 (2006) \\ 1811! 30 & GGA, Armiento-Mattsson (AM05) spin-unpolarised functional, 1812! {\it Phys. Rev. B} {\bf 72}, 085108 (2005) \\ 1813! 100 & Libxc functionals; the second and third values of {\tt xctype} define 1814! the exchange and correlation functionals in the Libxc library, 1815! respectively \\ 1816! \end{tabularx} 1817! 1818! \section{Contributing to Elk} 1819! Please bear in mind when writing code for the Elk project that it should be 1820! an exercise in physics and not software engineering. All code should 1821! therefore be kept as simple and concise as possible, and above all it should 1822! be easy for anyone to locate and follow the Fortran representation of the 1823! original mathematics. We would also appreciate the following conventions 1824! being adhered to: 1825! \begin{itemize} 1826! \item Strict Fortran 2003 should be used. Features which are marked as 1827! obsolescent in Fortran 2003 should be avoided. These include assigned 1828! format specifiers, labeled do-loops, computed goto statements and statement 1829! functions. 1830! \item Modules should be used in place of common blocks for declaring 1831! global variables. Use the existing modules to declare new global variables. 1832! \item Any code should be written in lower-case free form style, starting 1833! from column one. Try and keep the length of each line to fewer than 80 1834! characters using the \& character for line continuation. 1835! \item Every function or subroutine, no matter how small, should be in its 1836! own file named {\tt routine.f90}, where {\tt routine} is the function or 1837! subroutine name. It is recommended that the routines are named so as to 1838! make their purpose apparent from the name alone. 1839! \item Use of {\tt implicit none} is mandatory. Remember also to define the 1840! {\tt intent} of any passed arguments. 1841! \item Local allocatable arrays must be deallocated on exit of the routine to 1842! prevent memory leakage. Use of automatic arrays should be limited to arrays 1843! of small size. 1844! \item Every function or subroutine must be documented with the Protex source 1845! code documentation system. This should include a short \LaTeX\ description 1846! of the algorithms and methods involved. Equations which need to be 1847! referenced should be labeled with {\tt routine\_1}, {\tt routine\_2}, etc. 1848! The authorship of each new piece of code or modification should be 1849! indicated in the {\tt REVISION HISTORY} part of the header. See the Protex 1850! documentation for details. 1851! \item Ensure as much as possible that a routine will terminate the program 1852! when given improper input instead of continuing with erroneous results. 1853! Specifically, functions should have a well-defined domain for which they 1854! return accurate results. Input outside that domain should result in an 1855! error message and termination. 1856! \item Report errors prior to termination with a short description, for 1857! example: 1858! \begin{verbatim} 1859! write(*,*) 1860! write(*,'("Error(readinput): natoms <= 0 : ",I8)') natoms(is) 1861! write(*,'(" for species ",I4)') is 1862! write(*,*) 1863! stop 1864! \end{verbatim} 1865! \item Wherever possible, real numbers outputted as ASCII data should be 1866! formatted with the {\tt G18.10} specifier. 1867! \item Avoid redundant or repeated code: check to see if the routine you need 1868! already exists, before writing a new one. 1869! \item All reading in of ASCII data should be done in the subroutine 1870! {\tt readinput}. For binary data, separate routines for reading and writing 1871! should be used (for example {\tt writestate} and {\tt readstate}). 1872! \item Input filenames should be in lowercase and have the extension 1873! {\tt .in} . All output filenames should be in uppercase with the extension 1874! {\tt .OUT} . 1875! \item All internal units should be atomic. Input and output units should be 1876! atomic by default and clearly stated otherwise. Rydbergs should not be used 1877! under any circumstances. 1878! \end{itemize} 1879! \subsection{Licensing} 1880! Routines which constitute the main part of the code are released under the 1881! GNU General Public License (GPL). Library routines are released under the 1882! less restrictive GNU Lesser General Public License (LGPL). Both licenses 1883! are contained in the file {\tt COPYING}. Any contribution to the code must 1884! be licensed at the authors' discretion under either the GPL or LGPL. 1885! Author(s) of the code retain the copyrights. Copyright and (L)GPL 1886! information must be included at the beginning of every file, and no code 1887! will be accepted without this. 1888! 1889!EOI 1890 1891