Using the extended tight-binding Hamiltonian

Now that you know how to perform an extended tight-binding calculation (xTB) on a non-periodic system, it is time to go more in-depth into the settings used for that calculation. After all, the results produced by the calculation should be accurate and reliable. Also, using the correct settings for the xTB Hamiltonian will play a crucial part during geometry optimizations, as this can improve convergence speed and results accuracy.

The best place to start is by having a look at the Hamiltonian block from the dftb_pin.hsd file generated during our previous simulation. That way, all possible settings are visible. Fore more details on each entry you can check Chapter 2.4 from the DFTB+ manual.

To make the settings block more clear, the keys can be shuffled to create more consistent blocks, which will be explained below.

Hamiltonian = xTB {
  Method = "GFN1-xTB"
  ShellResolvedSCC = Yes
  SCC = Yes
  ReadInitialCharges = No
  InitialCharges = {}
  SCCTolerance = 1.0000000000000001E-005
  ConvergentSCCOnly = Yes
  SpinPolarisation = {}
  ElectricField = {}
  Solver = RelativelyRobust {}
  Charge = 0.0000000000000000
  MaxSCCIterations = 100
  Dispersion = {}
  Solvation = {}
  Electrostatics = GammaFunctional {}
  Differentiation = FiniteDiff {
    Delta = 1.2207031250000000E-004
  }
  ForceEvaluation = "Traditional"
  Mixer = Broyden {
    MixingParameter = 0.20000000000000001
    InverseJacobiWeight = 1.0000000000000000E-002
    MinimalWeight = 1.0000000000000000
    MaximalWeight = 100000.00000000000
    WeightFactor = 1.0000000000000000E-002
  }
  Filling = Fermi {
    Temperature = 9.5004460357391702E-004
  }
}

Hamiltonian = xTB {
  Method = "GFN1-xTB"
  SCC = Yes
  SCCTolerance = 1.0000000000000001E-005
  MaxSCCIterations = 100
  ShellResolvedSCC = Yes
  ConvergentSCCOnly = Yes
  Solver = RelativelyRobust {}
  Mixer = Broyden {
    MixingParameter = 0.20000000000000001
    InverseJacobiWeight = 1.0000000000000000E-002
    MinimalWeight = 1.0000000000000000
    MaximalWeight = 100000.00000000000
    WeightFactor = 1.0000000000000000E-002
  }
  Filling = Fermi {
    Temperature = 9.5004460357391702E-004
  }
  Charge = 0.0000000000000000
  InitialCharges = {}
  ReadInitialCharges = No
  SpinPolarisation = {}
  ElectricField = {}
  Solvation = {}
  Dispersion = {}
  Electrostatics = GammaFunctional {}
  ForceEvaluation = "Traditional"
  Differentiation = FiniteDiff {
    Delta = 1.2207031250000000E-004
  }
}

Extended tight-binding models

Method = "GFN1-xTB"

The first tag, Method, is used to specify the extended tight-binding method used for the calculations. Currently, DFTB+ supports three models:

• GFN1-xTB
• IPEA-xTB
• GFN2-xTB

Much like DFT functionals, which one is best depends on the system studied and on the chemical species involved. GFN1-xTB is the first extended tight-binding model proposed by Grimme et al. and covers parameters for all atomic species with Z ≤ 86. IPEA-xTB is a reparameterization of GFN1-xTB aimed at predicting better ionization potentials. To do so, it expands the basis set of some elements, while keeping most identical to GFN1-xTB. GFN2-xTB is the newest model available, featuring better descriptions for halogens and hydrogen, while also using D4 dispersion corrections as opposed to the D3 corrections used in GFN1-xTB and IPEA-xTB. The advantages of using GFN2-xTB over GFN1-xTB depend on the system, but it is generally observed that GFN2-xTB has more convergence problems for bad geometries compared to its predecessor. The table below contains information about the basis sets used for each atomic species:

The elements that have an asterisk use incorrect basis sets in some of the three models. These are known bugs that were part of the original parameterizations and thus require significant changes to the models to correct. As such, they will remain unchanged.

Self-convergent charges (SCC)

The most important block in the Hemiltonian block is the SCC part, which directly impacts the performance and accuracy of the calculations.

SCC = Yes
SCCTolerance = 1.0000000000000001E-005
MaxSCCIterations = 100
ShellResolvedSCC = Yes
ConvergentSCCOnly = Yes
Solver = RelativelyRobust {}
Mixer = Broyden {
  MixingParameter = 0.20000000000000001
  InverseJacobiWeight = 1.0000000000000000E-002
  MinimalWeight = 1.0000000000000000
  MaximalWeight = 100000.00000000000
  WeightFactor = 1.0000000000000000E-002
}

The first tag, SCC, enables the self-convergent charge procedure during the calculations. Without it, the calculation would use the charge distribution of neutral atoms, which for most systems are not representative. This setting is enabled by default and should only be disabled if you know the neutral atoms provide a decent description. If enabled, the SCC procedure has two stop conditions:

• When the change in charges from cycle to cycle is less than the SCCTolerance;
• When MaxSCCIterations cycles have been performed.

Next, the tag ShellResolvedSCC determines if all distinct Hubbard U values for the different atomic angular momenta shells are used. However, this does not seem to have an impact on xTB calculations. The tag ConvergentSCCOnly indicates that forces and other properties are only computed after the SCC procedure converges, not after each cycle.

The last two tags in this section, Solver and Mixer, affect the charge calculation procedure. Solver represents the eigensolver used by DFTB+, which is set to RelativelyRobust by default. Alternatives to this are DivideAndConqer (which requires twice the memory of the other solvers) and QR (which is based on QR decomposition). Depending on how the code is compiled, more solvers can be available, such as MAGAM (GPU solver), ELPA (efficient and GPU), OMM, PEXSI, or NTPoly. These have to be added manually at compile-time and have their own advantages and disadvantages, which should be assessed before usage.

On the other hand, Mixer determines that way past solutions from the SCC procedure are "mixed" to create a new solution. You may be tempted to think that during the self-convergence procedure, the charges (or electron density for SCF) computed for one cycle are reused as the starting point of the next cycle. While this can be accomplished by using the Simple mixer and setting the MixingParameter to 1.0, most quantum calculation software use complex methods to determine the best starting point for the next cycle. Both DFTB+ and VASP use Broyden mixing by default. Alternatives to this, which require less memory, are Anderson and DIIS. I never used them myself, but it can be a good idea to try them as well to see if it improves convergence and performance.

Electron filling and smearing

The Filing block specifies the partial occupation method and electron smearing used in the calculations. Currently, two methods are available: Fermi and MethfesselPaxton. Both methods take as argument Temperature, which controls the smearing in the system much like SIGMA does in VASP, with the difference that DFTB+ uses Hartree as default units. By default, DFTB+ uses a smearing of 300 K, which can be considered low enough for some applications. Nevertheless, according to Basiuk et al., if HOMO-LUMO gaps or other electronic properties are of interest, the smearing can be lowered to 30 K for more accurate results.

Filling = Fermi {
  Temperature = 9.5004460357391702E-004
}

Filling = MethfesselPaxton {
  Order = 1
  Temperature = 7.3498586E-003 
}

Which of the two fractional occupation methods produces the best results depends on the system studied. DFTB+ defaults to Fermi filling with 300 K smearing, while VASP defaults to 1st order MethfesselPaxton filling with 0.2 eV (2321 K) smearing. Nevertheless, the latter can be inaccurate for insulators due to unphysical partial occupancies.

System charges

This part of the input controls the charges inside the system. The Charge property sets the system charge, which is important when ionic systems are (e.g. acetate, chloride, etc.). For neutral systems, it should be left unchanged. For ionic system, this should be paired with the following tag, InitialCharges, which assigns the correct charges to the respective atoms.

Charge = 0.0000000000000000
InitialCharges = {}
ReadInitialCharges = No

The following tag, ReadInitialCharges, allows for previously computed results to be read, which provides the basis of the restart system in DFTB+. As mentioned previously, simulations create a charges.bin file that stores information about the charges computed during the simulation. Loading these charges at start-time would result in the simulation continuing from those charges, which is effectively a restart.

Spin-polarized calculations

While this will be explore in depth in another tutorial, the SpinPolarisation tag enables spin-polarised calculations. By default no value is assigned to it, forcing DFTB+ to ignore spin polarisation. By changing it to Colinear or NonColinear, the user can introduce polarisation in the system, which might be necessary to accurately model some species (e.g. manganese).

SpinPolarisation = {}

Implicit effects

As implicit effects, one can introduce an electric field or a solvent in the simulation, altering the behavior of the system. For zwitterionic systems the introduction of solvation is recommended as it improves the convergence of the SCC cycle. I am not currently interested in simulations with either properties, so unfortunately I cannot provide more information than what is already stated in the DFTB+ manual.

ElectricField = {}
Solvation = {}

Hamiltonian modifications

This part of the settings provides the ability to alter the Hamiltonian used in extended tight-binding, but in practice has little utility. While not specifically stated anywhere, the default dispersion schemes used in GFN1-xTB, IPEA-xTB, and GFN2-xTB are already part of the base parameterization, and thus does not require manual additions using the Dispersion tag. Meanwhile, the Electrostatics tag should only be changed from GammaFunctional to Poisson during transport calculations.

Dispersion = {}
Electrostatics = GammaFunctional {}

Force calculation

The last part covers how the forces are being calculated during the simulation. ForceEvaluation can be performed in three ways: traditional, dynamics, and dynamicsT0. The first one is the default version used to calculate forces in DFTB, while the other two are meant to be used for non-converged charges. In each case, the differentiation can be performed using either finite difference (FiniteDiff), or using Richardson extrapolation, which is the most accurate method to determine forces.

ForceEvaluation = "Traditional"
Differentiation = FiniteDiff {
  Delta = 1.2207031250000000E-004
}