Example Calculation

The easiest way to get started with aiida-tbextraction is by modifying an existing example. For this reason, we will describe such an example here.

Workflow Structure

The following diagram shows the structure of the OptimizeFirstPrinciplesTightBinding workflow. In a first step, the first-principles workflow is called to create the necessary inputs for the Wannier90 calculations. Next, the WindowSearch workflow is called, which runs an optimization (using aiida-optimize) for the energy windows. The individual calculations are carried out by the RunWindow workflow, which is shown in the second part of the diagram.

First, the RunWindow workflow calculates the tight-binding model using the TightBindingCalculation workflow. This runs Wannier90, and then parses the output to the TBmodels HDF5 format. If requested, the model is also sliced and symmetrized. Then, the resulting tight-binding model is evaluated using a subclass of the ModelEvaluationBase workflow. This workflow assigns a “cost value” to the tight-binding model, which is a measure for the quality (lower is better) of the model. The cost value is what’s used by the optimization workflow to find the optimal energy windows.


Optimization workflow diagram

The workflow which runs the first-principles calculations can be different depending on which code you would like to run. Currently, only Quantum ESPRESSO is supported by default, with the QuantumEspressoFirstPrinciplesRun workflow. It first runs an SCF calculation, and then uses the QuantumEspressoReferenceBands and QuantumEspressoWannierInput workflows to calculate the reference band structure and Wannier90 input files, respectively.

Launching the Workflow

The following example runs a simple tight-binding calculation for InSb. It builds up the required inputs step by step, using a process builder. To test that your setup has worked correctly, you can try launching this example. The code is also available on GitHub. You will have to adjust some values specific to your configuration, like the queue name on the compute cluster or the names of the codes.

Once you have successfully ran this example, you can start adapting it to your use case. The Reference section can help you understand the individual input values. For production runs, you should decrease the tolerances window_tol and cost_tol (or use their default values), and use more accurate inputs for the DFT and Wannier calculations (for example, you should increase dis_num_iter).

#!/usr/bin/env runaiida
# -*- coding: utf-8 -*-

# © 2017-2019, ETH Zurich, Institut für Theoretische Physik
# Author: Dominik Gresch <greschd@gmx.ch>

import os
import pathlib

import numpy as np
from ase.io.vasp import read_vasp

from aiida import orm
from aiida.engine.launch import submit

from aiida_tbextraction.fp_run import QuantumEspressoFirstPrinciplesRun
from aiida_tbextraction.model_evaluation import BandDifferenceModelEvaluation
from aiida_tbextraction.optimize_fp_tb import OptimizeFirstPrinciplesTightBinding

# Define constants for Code and metadata
# Modify these to match your configured codes, and adapt the metadata
# as needed.
METADATA_PW = {
    'options': {
        'resources': {
            'num_machines': 1,
            'num_mpiprocs_per_machine': 4
        },
        'withmpi': True,
        'max_wallclock_seconds': 1200
    }
}
CODE_PW = orm.Code.get_from_string('pw')
METADATA_WANNIER = {
    'options': {
        'resources': {
            'num_machines': 1,
            'num_mpiprocs_per_machine': 1
        },
        'withmpi': False,
        'max_wallclock_seconds': 1200
    }
}
CODE_WANNIER = orm.Code.get_from_string('wannier90')
METADATA_PW2WANNIER = METADATA_WANNIER
CODE_PW2WANNIER = orm.Code.get_from_string('pw2wannier90')


def create_builder():
    builder = OptimizeFirstPrinciplesTightBinding.get_builder()

    # Add the input structure
    builder.structure = orm.StructureData()
    builder.structure.set_ase(read_vasp('inputs/POSCAR'))

    # Specify that QuantumEspressoFirstPrinciplesRun should be used to run the first-principles calculations
    builder.fp_run_workflow = QuantumEspressoFirstPrinciplesRun

    # Set the inputs for the QuantumEspressoFirstPrinciplesRun workflow
    common_qe_parameters = orm.Dict(
        dict=dict(
            CONTROL=dict(etot_conv_thr=1e-3),
            SYSTEM=dict(noncolin=True, lspinorb=True, nbnd=36, ecutwfc=30),
        )
    )

    pseudo_dir = pathlib.Path(__file__
                              ).parent.absolute() / 'inputs' / 'pseudos'
    pseudos = {
        'In':
        orm.UpfData(
            file=str(pseudo_dir / 'In.rel-pbe-dn-kjpaw_psl.1.0.0.UPF')
        ),
        'Sb':
        orm.UpfData(file=str(pseudo_dir / 'Sb.rel-pbe-n-kjpaw_psl.1.0.0.UPF'))
    }

    # We use the same general Quantum ESPRESSO parameters for the
    # scf, nscf, and bands calculations. The calculation type and
    # k-points will be set by the workflow.
    repeated_pw_inputs = {
        'pseudos': pseudos,
        'parameters': common_qe_parameters,
        'metadata': METADATA_PW,
        'code': CODE_PW
    }

    builder.fp_run = {
        'scf': repeated_pw_inputs,
        'bands': {
            'pw': repeated_pw_inputs
        },
        'to_wannier': {
            'nscf': repeated_pw_inputs,
            'pw2wannier': {
                'code': CODE_PW2WANNIER,
                'metadata': METADATA_PW2WANNIER
            },
            'wannier': {
                'code': CODE_WANNIER,
                'metadata': METADATA_WANNIER
            }
        }
    }

    # Setting the k-points for the reference bandstructure
    kpoints_list = []
    kvals = np.linspace(0, 0.5, 20, endpoint=False)
    kvals_rev = np.linspace(0.5, 0, 20, endpoint=False)
    for k in kvals_rev:
        kpoints_list.append((k, k, 0))  # Z to Gamma
    for k in kvals:
        kpoints_list.append((0, k, k))  # Gamma to X
    for k in kvals:
        kpoints_list.append((k, 0.5, 0.5))  # X to L
    for k in kvals_rev:
        kpoints_list.append((k, k, k))  # L to Gamma
    for k in np.linspace(0, 0.375, 21, endpoint=True):
        kpoints_list.append((k, k, 2 * k))  # Gamma to K
    builder.kpoints = orm.KpointsData()
    builder.kpoints.set_kpoints(
        kpoints_list,
        labels=[(i * 20, label)
                for i, label in enumerate(['Z', 'G', 'X', 'L', 'G', 'K'])]
    )

    # Setting the k-points mesh used to run the SCF and Wannier calculations
    builder.kpoints_mesh = orm.KpointsData()
    builder.kpoints_mesh.set_kpoints_mesh([6, 6, 6])

    builder.wannier.code = CODE_WANNIER
    builder.code_tbmodels = orm.Code.get_from_string('tbmodels')

    # Setting the workflow to evaluate the tight-binding models
    builder.model_evaluation_workflow = BandDifferenceModelEvaluation

    # Setting the additional inputs for the model evaluation workflow
    builder.model_evaluation = dict(
        code_bands_inspect=orm.Code.get_from_string('bands_inspect')
    )

    # Set the initial energy window value
    builder.initial_window = orm.List(list=[-1, 3, 10, 18])

    # Tolerance for the energy window.
    builder.window_tol = orm.Float(1.5)
    # Tolerance for the 'cost_value'.
    builder.cost_tol = orm.Float(0.3)
    # The tolerances are set higher than might be appropriate for a 'production'
    # run to make the example run more quickly.

    # Setting the parameters for Wannier90
    builder.wannier_parameters = orm.Dict(
        dict=dict(
            num_wann=14,
            num_bands=36,
            dis_num_iter=100,
            num_iter=0,
            spinors=True,
            # exclude_bands=range(1, )
        )
    )
    # Choose the Wannier90 trial orbitals
    builder.wannier_projections = orm.List(
        list=['In : s; pz; px; py', 'Sb : pz; px; py']
    )
    # Set the resource requirements for the Wannier90 run
    builder.wannier.metadata = METADATA_WANNIER

    # Set the symmetry file
    builder.symmetries = orm.SinglefileData(
        file=os.path.abspath('inputs/symmetries.hdf5')
    )

    # Pick the relevant bands from the reference calculation
    builder.slice_reference_bands = orm.List(list=list(range(12, 26)))

    return builder


if __name__ == '__main__':
    builder = create_builder()
    node = submit(builder)
    print('Submitted workflow with pk={}'.format(node.pk))

Note that the result of this example does not accurately represent the band-structure of InSb. This is due to the limitation of PBE, and could be overcome by using costlier hybrid functional calculations.

Note

Using hybrid functionals is currently not supported in the QuantumEspressoFirstPrinciplesRun workflow, due to limitations in QE itself.

Customizing the workflow

There are different ways in which you can customize the tight-binding extraction workflow to match your use case. If you already have the output from a first-principles calculation, you can directly use either the WindowSearch, RunWindow or TightBindingCalculation workflow, depending on which features you want to run. Also, you could change the workflow which is uded to run the first-principles calculations. When implementing such a workflow, you should make sure that it adheres to the interface defined in FirstPrinciplesRunBase, meaning that it should take at least the inputs defined there, and create the corresponding outputs. The easiest way to do this is by inheriting from the base class. The workflow can take additional inputs, which you can pass through the fp_run namespace in the OptimizeFirstPrinciplesTightBinding workflow inputs. Finally, you can also use a different workflow for evaluating the tight-binding model. This allows you to create a measure for the quality of the model which fits best to your specific use case.