Basic usage¶
The package is implemented as a Python extension library, which makes embedding ab initio calculations into complex workflows much easier. However, for basic uses there is a possibility to run the code in a more traditional manner.
The simplest way to run the code is to create an input config file, e.g., example.in, and pass it to the main program as follows:
> python -m gfalm.main example.in
or simply
> gfalm example.in
The format of the input file is described in this documentation.
A slightly more advanced way of using the code is to execute the following Python script:
import gfalm
from gfalm.utils import read_input_file
input_data = read_input_file('example.in')
result = gfalm.main.run(**input_data)
This way of running the code is the first step to more advanced usage of the package. Furthermore, one gets access to a result object containing various output data.
The advantage of the scripting approach is that all input parameters from the user input file, as well as default values for omitted parameters, are available in input_data dictionary and can be modified before being passed to gfalm.main.run() function.
Input file¶
Warning
The package is still in the alpha-stage, therefore the details of the file format might change in future versions. At the same time, validation of input parameters will also be improved, which implies that incorrect parameters/values will be reported to the user.
The input file has a simple config-like format. The input parameters are grouped into sections. There are five principal sections:
[Control]: various parameters controlling calculations
[Structure]: crystal structure
[Alloy]: alloy structure
[Kpoints]: k-points and Brillouin zone parameters
[Atomic]: radial-mesh and atomic-solver parameters
[Contour]: energy contour parameters
[DOS]: density of states parameters
[BSF]: parameters of the Bloch spectral function
The minimal mandatory set of parameters to start a self-consistent calculation is the following:
[Control]
job_name
[Structure]
lattice_scale
orws_radius
lattice
sites_frac
orsites_cart
[Kpoints]
mesh_size
The atom positions are defined as three float numbers followed by the label of the atomic species (in a simple case, just an element symbol) and possibly a set of site parameters. Among the important site parameters are mom specifying the initial spin splitting, lsm_mode specifying the mode of the local mangetic state (LSM here stands for the local spin moment). Obviously, these parameters are only relevant for spin polarized calculations. Possible LSM modes are described in more details in section [Structure] below.
This defines the crystal structure and the density of the \(k\)-mesh in the Brillouin zone. The linear scale of the crystal structure can be defined in two alternative ways. One is the lattice_scale parameter which is a simple multiplication factor for the matrix of translation vectors (lattice). Another is ws_radius, which defines the value of the average Wigner-Seitz radius corresponding to a sphere whose volume is equal to the average volume per atom. Both linear scales can be given either in Ånströms, which is specified by a token A following the number (e.g., 3.6 A), or in atomic units, specified by a token au (e.g., 2.7 au). When no specific token is given the default input units are used, which are Rydberg atomic units (can be controlled by parameters length_unit, energy_unit in the [Control] section).
Besides, it is important to keep in mind the following parameters:
\(l\)-cut-off (lmax in the [Structure] and [Alloy] sections);
parameters of the energy integration contour, especially the depth (depth) and the number of integration points (nz_sections).
Furthermore, a setup of an alloy compound requires section [Alloy].
Other essential options are described below in subsequent subsecions.
Control parameters¶
This section defines parameters that control execution, convergence, and method choice. The only mandatory parameter in this section is job_name. During the run of the code a log file <job_name>.log and a working directory job.<job_name> will be created. The working directory contains internal, mostly non-human-readable, data, which is mainly used for storing the information necessary to restore the state after a restart.
Parameter job_type defines the type of the calculation performed. The default is sc, which means a self-consistent DFT calculation. After such a calculation is performed, one could choose non-self-consistent modes dos or bsf to produce, respectively, the density of states and the Bloch spectral function. Two other options, check and atomic, are for debug purposes only and should never be used otherwise.
The exchange-correction (XC) type is selected by parameter xc_type. Note that this sets the XC type for the self-consistency. Once the DFT self-consistency is reached the total energy is calculated using three XC functionals: LDA, PBE, and PBEsol. This makes sense for most of the metallic systems due to the variational property of the DFT. The difference in the electronic structure and hence the charge density produced by different local XC functionals is very small for metallic systems (except for magnetic metals, see the remark of caution below), while the change in the total energy because of the charge density variation is second order in the variation and can be neglected. Thus, for example, the total energy calculated a posteriori using a PBE XC functional with the charge density obtained with LDA will be practically indistinguishable from the total energy obtained by a fully self-consistent PBE calculation.
That said, one should be careful with systems in which the choice of the XC functional affects the electronic structure itself. This happens in non-metallic systems (where the gap is definitely different in LDA and GGA) and in magnetic metallic systems, where the magnetic moment might differ (for a given volume) for an LDA or PBE self-consistent calculation. However, in our experience, the magnetic moment obtained with LDA is in a better agreement with experiment when calculated at the experimental volume. At the same time, the overall error in the total energy is larger than in PBE (or PBEsol). Therefore, we recommend to use LDA for self-consistency but then take the total energy calculated either with PBE or PBEsol XC functional. This is tantamount to performing a self-consistent PBE(sol) calculation, with the magnetic moment fixed to its LDA value.
Another important parameter is spin_polarized whose meaning should be clear by itself. Note that for a DLM calculation (see the description of [Strucutre] and [Alloy] sections on how to set it up) this option should be set to yes despite the non-spin-polarized total DOS and charge density resulting from such a calculation.
Parameter temperature sets the electronic temperature in Kelvins. Specifically, when this parameter is non-zero, the Green’s function is integrated with the Fermi function (using a special contour type Fermi that is chosen automatically for a non-zero temperature). This produces certain degree of smearing of the DOS, especially at large temperatures. The temperature parameter is also used in the LSF mode (see [Strucutre] and [Alloy] below). Usually, LSF implies finite temperature Green’s function integration.
The last parameter affecting calculation results is kappa_correction. Currently, the correction is implemented only a posteriori, i.e., it is applied after the ASA self-consistency is reached. It should generally improve the electronic structure (i.e., DOS). However, there are some unresolved issues with the total energy when the correction is enabled (the total energy is actually also improved but an unknown systematic bias result in the incorrect prediction of the equilibrium volume). For this reason, we do not recommend using the total energy with this correction. At the same time, the correction should be enabled if one wants to produce plots of the DOS or BSF.
All other parameters in the [Control] section only influence input/output and the convergence rate.
The self-consistent run performs maximum max_iter DFT iterations. The iterations usually stop earlier when the convergence criterion is achieved. The convergence criterion is that each of the three errors – total energy, chrage density, spin density – must be smaller than their respective tolerance levels, etot_tol, chd_tol, spd_tol. Within each DFT iteration a series of sub-iterations (so-called “small” or “linear” iteration) is performed. A sub-iteration solves an approximate single-site KKR problem and updates the charge density and the Fermi level. These iterations are meant to accelerate CPA calculations but they also improve the convergence in case of ordered compounds. The maximum number of these iterations is given by max_lin_iter and the convergence criterion here is that the change in the Fermi level is smaller than ef_tol.
DFT convergence also requires admixing of the previous charge density to the new one. GreenALM additionally performs linear mixing of the Kohn-Sham potential, which can be very helpful for improving the convergence of magnetic systems. The potential mixing parameter is set by pot_mix. The more important charge density mixing is controlled by a block of parameters chd_mixing. The charge density mixing is performed as follows: in the first n_init_iter iterations only linear mixing is done, with the coefficients for charge and spin densities being init_chd_mix and init_spd_mix, repsectively. Then, the main mixers are initialized. Currently, two mixers are implemented: linear (simple) and Anderson mixing. The default parameters are quite safe and work in most of the cases. If you want to change mixing parameters do it carefully, a small step at a time.
In case of spin-polarized calculations there are two modes of mixing:
Charge and spin density are mixed together. In particular, in the Adnerson mixer this implies that the difference vectors contain contributions from both densities. This mode is activated when spd_type = none, i.e., there is no separate spin density mixer. spd_params in this case is ignored (but not init_spd_mix).
There are two separate mixers for the charge and spin densities. The mode is activated when spd_type is set to anything but none. This can be useful for certain magnetic systems in which the local magnetic moments converge much slower than the charge density. One of the workarounds is to set spd_type = linear with a small mixing coefficient (0.1 or even lower).
Structure parameters¶
This section defines a crystal structure as well as parameters of the structure constant matrix (SCM), which is an indispensible part of any KKR method. A crystal structure is defined by a linear scale (lattice_scale or ws_radius), by the matrix of translation vectors, lattice, and by site positions given either by sites_frac as fractional coordinates or by sites_cart as cartesian coordinates.
The linear scale can be defined by two alternative ways. The simplest is the lattice_scale, which is just a multiplication factor for the matrix lattice. This allows straightforward definition of structures for which the lattice constant is known. For instance, the following example,
lattice_scale = 3.6 A
lattice =
0.5 0.5 0.0
0.0 0.5 0.5
0.5 0.0 0.5
defines an fcc structure with the lattice constant equal to 3.6 \(\AA\).
There are, however, calculations where it is more convenient to specify the average volume per atom rather than the lattice constant (e.g., calculate energy difference between fcc and hcp structures for a fixed volume). A linear scale corresponding to the average volume per atom, \(\overline{\Omega}\), is traditionally specified by the Wigner-Seitz (WS) radius, \(r_{WS}\), defined as \(\overline{\Omega} = 4 \pi / 3 r_{WS}^3\), i.e., this is a radius of a sphere with the volume equal to \(\overline{\Omega}\). The WS radius is specified by parameter ws_radius.
Note that one should provide only one of the linear scales. If both are given, a warning will be issued and ws_radius will be used in calculations.
The lattice traslations and the linear scale define a unit cell of the system. The atomic positions inside the unit cell are given by sites_frac or site_cart followed by a block (must be indented) in which each line specifies a site. Each site has four obligatory parameters: three coordinates (fractional or cartesian) and a species label. In the simplest case, the species label is just an element symbol, e.g.,
site_frac =
0.0 0.0 0.0 Cu
A more complex situation when a site (sublattice) is occupied by an alloy will be considered in the next section.
Furthermore, there are other site parameters that can follow the species label. First, one could provide a site-specific \(l\)-cut-off lcut. This especially makes sense when the site is occupied by an empty sphere (as a way of basis extension) which does not have its own electrons and can often be described by a lower \(l\)-cut-off.
Another important set of parameters are only relevant to spin-polarized calculations and related to the local spin moment (LSM) mode. Currently, four LSM modes are implemented:
lsm_mode = none: this corresponds to a calculation of the equilibrium state of an ordered magnet, e.g., ferromagnet or anti-ferromagnet.
lsm_mode = fxm: fixed spin moment calculations, whereby the magnitude of the local magnetic moment will be constrained to the value set by mom. Note that the desired LSM is achieved by finding a constraint local magnetic field. If the target LSM is unphysical (e.g., 4 \(\mu_B\) for Ni) the constraint magnetic field might become very large and the calculation will go awry.
lsm_mode = dlm: disordered local moment (DLM) calculation, which mean that the total site moment will be zero but the magnitude of a single up- or down-component can be non-zero. Note that this setting implies the time reversal symmetry, i.e., that there are no sites with simultaneously lsm_mode = none and mom > 0. Moreover, this option automatically generates an alloy sublattice with two components and activates thus CPA calculation. In contrast to manually setting up an alloy of up- and down-compoents, this mode will force time reversal symmetry. In practice, this means that quantities (SPO, charge density, local KS potential, etc.) are evaluated only for one component, while for another one spin-flip transformation is used.
lsm_mode = lsf: LSF model described in the Introduction is used for the local spin moment. This setting implies that parameter temperature in section [Control] is set to a non-negligible value. Another LSF parameter given by lsf_params defines an effective dimension of spin fluctuations and can be set to 1, 2, 3, or -1. Generally, for very itinerant elements, such as, e.g., Ni, one should always choose lsf_params = 3. For magnetic elements with a localized spin, such as e.g., bcc Fe, one should choose lsf_params = -1, which is a special choice corresponding to an entropy contribution of the form \(S_{mag} = \log (1 + M)\). Other options are for intermediate cases, where a particular choice of the parameter can be inferred from a series of fixed spin moment calculations, which produce the total energy as a function of the local magnetic moment.
LSM modes can be combined by giving them as a comma-separated list. The following combinations are meaningful:
lsm_mode = fxm, dlm # DLM calculation with the LSM fixed to `mom`
lsm_mode = dlm, lsf # DLM calculation with LSF for the local moment
At the same time, LSM modes fxm and lsf are mutually exclusive for obvious reasons.
For a given crystal structure a SCM is initialized in the beginning of the calculations. Apart from the crystal structure parameters, SCM is specified by the \(l\)-cut-off, lcut, at each lattice site. The default cut-off is given by parameter lmax of the [Structure] section itself. It is used for all sites for which parameter lcut is not set explicitly.
For performance reasons, the generated SCMs (both real- and reciprocal-space SCM are generated) are cached into a file. The name of this file can be set by parameter name, which is also a convenient way of naming the structure. The path for the SCM cache files can be set by parameter path. This allows one to transfer SCM cache files from one set of calculations to another.
Alloy parameters¶
This section contains definitions of pseudo-species that represent substitutional disorder. A pseudo-species is designated by a unique label (consisting of alphanumeric characters) and a list of alloy components given by atomic species along with their atomic fractions, i.e.,
FeNi =
Fe 0.7
Ni 0.3
For each component one can specify optional parameters. Most of the parameters are the same as the site parameters in the [Structure] section (except lcut which does not make sense for a single alloy component). There are, however, two additional parameters specific to CPA: al_scr, bet_scr – screening parameters corresponding to, respectively, \(\alpha_{\mu}\), \(\beta_{\mu}\) described in the CPA section of the Introduction. For instance, a more advanced version of the above example can look like this:
FeNi =
Fe 0.7 lsm_mode = dlm, lsf lsf_params = 3 al_scr = 0.75 bet_scr = 1.05
Ni 0.3 lsm_mode = dlm, lsf lsf_params = 3 al_scr = 0.71 bet_scr = 1.03
Once a pseudo-species is defined it can be used as a species in the structure definition. A reference to a pseudo-species is given in a bash-like format as ${pseudo-species-label}. For example, for the pseudo-species given above the reference is ${FeNi}.
Parameters of the Brillouin zone¶
The Brillouin zone is specified in section [Kpoints]. The only mandatory parameter in this section is mesh_size that specifies the mesh discretizing the Brillouin zone. Two other parameters are related to symmetries and should be changed only for debug purposes.
Parameters of the atomic solver¶
Section [Atomic] contains parameters related to the radial mesh and the atomic solver. The default values are suitable for almost any imaginable calculation and should be modified only in rare special cases.
Parameters of the energy contour¶
Section [Contour] specifies the energy integration contour in the complex plane. There are two types of contours: zero-temperature (circular and elliptic) and finite-temperature (fermi). A zero-temperature contour crosses the real axis at two points: close to the Fermi level, \(E_F\) and at a lower energy point referred to as the bottom of the contour, \(E_b\). The difference \(E_F - E_b\) is called the depth of the contour and is given by a parameter depth of this section.
A finite-temperature contour is needed for integration of the Green’s function multiplied by the Fermi function. The upper integration limit in this case is, formally, \(+\infty\) but because of the exponential decay of the Fermi function can be limited by an energy \(E_F + \gamma T\), where \(\gamma\) is around 10-40. The only finite-temperature contour currently implemented in GreenALM is a so-called Fermi contour and its upper limit is specified by the number of Matsubara poles, \(n_{Mat}\), corresponding to \(\gamma = 2 \pi n_{Mat}\)
Postprocessing parameters¶
Currently, two post-processsing calculations are available: the density of states (DOS) and the Bloch spectral function (BSF). Both of them are performed starting from the complete self-consistent calculation, i.e., they imply restart = forbid.
To carry out a DOS calculation one sets job_type to dos and adds a section [DOS] in the input file. The calculation is defined by three parameters in section [DOS]: energy range e_range, number of energy mesh points n_points, value of the small imaginary part im_delta. An important point is that the energies are defined with respect to the Fermi level read in from the stored state, implying that the Fermi level in the resulting energy mesh is equal to 0.0 by construction.
A BSF calculation is specified analogously by adding a section [BSF] and setting job_type to bsf. The parameters in the [BSF] section are somewhat more complicated because several modes are supported. The mode is defined by option type in the parameter block kset_params.
type = segments corresponds to a BSF calculated for a range of energies along a series of line segments in the reciprocal space. The end points of the segments should then be given as a list of \(k\)-points in parameter kpoints. Type of coordinates frac or cart, is given immediately after the equal sign. For example:
kpoints = frac
0.0 0.0 0.0
0.5 0.0 0.0
0.5 0.5 0.0
0.0 0.0 0.0
The number of points along each segment is controlled by a parameter segment_len inside parameter block kset_params.
type = sections corresponds to a BSF calculation in a 2D cross-section of the reciprocal space. Cross-sections are given by a list of cross-section specifications in a parameter sections, with the type of coordinates (frac or cart) being provided immediately after the equal sign. A cross-section specification starts with a label followed by a parameter line. Each parameter line contains origin and two vectors, k1 and k2 specifying two sides of a parallelogram. The number of \(k\)-points along each side is optionally given by a pair of integers in parameter size. All these parameters should be given in one line. If size is not given, the default value section_size from the kset_params block is used.