Subroutines

You can call a subroutine in this library as follows:

USE libtetrabz, ONLY : libtetrabz_occ

CALL libtetrabz_occ(ltetra,bvec,nb,nge,eig,ngw,wght)

Every subroutine has a name starts from libtetrabz_.

For the C program, it can be used as follows:

#include "libtetrabz.h"

libtetrabz_occ(&ltetra,bvec,&nb,nge,eig,ngw,wght)

Variables should be passed as pointers. Arrays should be declared as one dimensional arrays. Also, the communicator argument for the routine should be transformed from the C/C++’s one to the fortran’s one as follows.

comm_f = MPI_Comm_c2f(comm_c);

Total energy, charge density, occupations

\[\begin{align} \sum_{n} \int_{\rm BZ} \frac{d^3 k}{V_{\rm BZ}} \theta(\varepsilon_{\rm F} - \varepsilon_{n k}) X_{n k} \end{align}\]
CALL libtetrabz_occ(ltetra,bvec,nb,nge,eig,ngw,wght,comm)

Parameters

INTEGER,INTENT(IN) :: ltetra

Specify the type of the tetrahedron method. 1 \(\cdots\) the linear tetrahedron method. 2 \(\cdots\) the optimized tetrahedron method [1].

REAL(8),INTENT(IN) :: bvec(3,3)

Reciprocal lattice vectors (arbitrary unit). Because they are used to choose the direction of tetrahedra, only their ratio is used.

INTEGER,INTENT(IN) :: nb

The number of bands.

INTEGER,INTENT(IN) :: nge(3)

Specify the \(k\)-grid for input orbital energy.

REAL(8),INTENT(IN) :: eig(nb,nge(1),nge(2),nge(3))

The orbital energy measured from the Fermi energy ( \(\varepsilon_{\rm F} = 0\) ).

INTEGER,INTENT(IN) :: ngw(3)

Specify the \(k\)-grid for output integration weights. You can make ngw \(\neq\) nge (See Appendix).

REAL(8),INTENT(OUT) :: wght(nb,ngw(1),ngw(2),ngw(3))

The integration weights.

INTEGER,INTENT(IN),OPTIONAL :: comm

Optional argument. Communicators for MPI such as MPI_COMM_WORLD. Only for MPI / Hybrid parallelization. For C compiler without MPI, just pass NULL to omit this argment.

Fermi energy and occupations

\[\begin{align} \sum_{n} \int_{\rm BZ} \frac{d^3 k}{V_{\rm BZ}} \theta(\varepsilon_{\rm F} - \varepsilon_{n k}) X_{n k} \end{align}\]
CALL libtetrabz_fermieng(ltetra,bvec,nb,nge,eig,ngw,wght,ef,nelec,comm)

Parameters

INTEGER,INTENT(IN) :: ltetra

Specify the type of the tetrahedron method. 1 \(\cdots\) the linear tetrahedron method. 2 \(\cdots\) the optimized tetrahedron method [1].

REAL(8),INTENT(IN) :: bvec(3,3)

Reciprocal lattice vectors (arbitrary unit). Because they are used to choose the direction of tetrahedra, only their ratio is used.

INTEGER,INTENT(IN) :: nb

The number of bands.

INTEGER,INTENT(IN) :: nge(3)

Specify the \(k\)-grid for input orbital energy.

REAL(8),INTENT(IN) :: eig(nb,nge(1),nge(2),nge(3))

The orbital energy measured from the Fermi energy ( \(\varepsilon_{\rm F} = 0\) ).

INTEGER,INTENT(IN) :: ngw(3)

Specify the \(k\)-grid for output integration weights. You can make ngw \(\neq\) nge (See Appendix).

REAL(8),INTENT(OUT) :: wght(nb,ngw(1),ngw(2),ngw(3))

The integration weights.

REAL(8),INTENT(OUT) :: ef

The Fermi energy.

REAL(8),INTENT(IN) :: nelec

The number of (valence) electrons per spin.

INTEGER,INTENT(IN),OPTIONAL :: comm

Optional argument. Communicators for MPI such as MPI_COMM_WORLD. Only for MPI / Hybrid parallelization. For C compiler without MPI, just pass NULL to omit this argment.

Partial density of states

\[\begin{align} \sum_{n} \int_{\rm BZ} \frac{d^3 k}{V_{\rm BZ}} \delta(\omega - \varepsilon_{n k}) X_{n k}(\omega) \end{align}\]
CALL libtetrabz_dos(ltetra,bvec,nb,nge,eig,ngw,wght,ne,e0,comm)

Parameters

INTEGER,INTENT(IN) :: ltetra

Specify the type of the tetrahedron method. 1 \(\cdots\) the linear tetrahedron method. 2 \(\cdots\) the optimized tetrahedron method [1].

REAL(8),INTENT(IN) :: bvec(3,3)

Reciprocal lattice vectors (arbitrary unit). Because they are used to choose the direction of tetrahedra, only their ratio is used.

INTEGER,INTENT(IN) :: nb

The number of bands.

INTEGER,INTENT(IN) :: nge(3)

Specify the \(k\)-grid for input orbital energy.

REAL(8),INTENT(IN) :: eig(nb,nge(1),nge(2),nge(3))

The orbital energy measured from the Fermi energy ( \(\varepsilon_{\rm F} = 0\) ).

INTEGER,INTENT(IN) :: ngw(3)

Specify the \(k\)-grid for output integration weights. You can make ngw \(\neq\) nge (See Appendix).

REAL(8),INTENT(OUT) :: wght(ne,nb,ngw(1),ngw(2),ngw(3))

The integration weights.

INTEGER,INTENT(IN) :: ne

The number of energy where DOS is calculated.

REAL(8),INTENT(IN) :: e0(ne)

Energies where DOS is calculated.

INTEGER,INTENT(IN),OPTIONAL :: comm

Optional argument. Communicators for MPI such as MPI_COMM_WORLD. Only for MPI / Hybrid parallelization. For C compiler without MPI, just pass NULL to omit this argment.

Integrated density of states

\[\begin{align} \sum_{n} \int_{\rm BZ} \frac{d^3 k}{V_{\rm BZ}} \theta(\omega - \varepsilon_{n k}) X_{n k}(\omega) \end{align}\]
CALL libtetrabz_intdos(ltetra,bvec,nb,nge,eig,ngw,wght,ne,e0,comm)

Parameters

INTEGER,INTENT(IN) :: ltetra

Specify the type of the tetrahedron method. 1 \(\cdots\) the linear tetrahedron method. 2 \(\cdots\) the optimized tetrahedron method [1].

REAL(8),INTENT(IN) :: bvec(3,3)

Reciprocal lattice vectors (arbitrary unit). Because they are used to choose the direction of tetrahedra, only their ratio is used.

INTEGER,INTENT(IN) :: nb

The number of bands.

INTEGER,INTENT(IN) :: nge(3)

Specify the \(k\)-grid for input orbital energy.

REAL(8),INTENT(IN) :: eig(nb,nge(1),nge(2),nge(3))

The orbital energy measured from the Fermi energy ( \(\varepsilon_{\rm F} = 0\) ).

INTEGER,INTENT(IN) :: ngw(3)

Specify the \(k\)-grid for output integration weights. You can make ngw \(\neq\) nge (See Appendix).

REAL(8),INTENT(OUT) :: wght(ne,nb,ngw(1),ngw(2),ngw(3))

The integration weights.

INTEGER,INTENT(IN) :: ne

The number of energy where DOS is calculated.

REAL(8),INTENT(IN) :: e0(ne)

Energies where DOS is calculated.

INTEGER,INTENT(IN),OPTIONAL :: comm

Optional argument. Communicators for MPI such as MPI_COMM_WORLD. Only for MPI / Hybrid parallelization. For C compiler without MPI, just pass NULL to omit this argment.

Nesting function and Fr&oumlhlich parameter

\[\begin{align} \sum_{n n'} \int_{\rm BZ} \frac{d^3 k}{V_{\rm BZ}} \delta(\varepsilon_{\rm F} - \varepsilon_{n k}) \delta(\varepsilon_{\rm F} - \varepsilon'_{n' k}) X_{n n' k} \end{align}\]
CALL libtetrabz_dbldelta(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,comm)

Parameters

INTEGER,INTENT(IN) :: ltetra

Specify the type of the tetrahedron method. 1 \(\cdots\) the linear tetrahedron method. 2 \(\cdots\) the optimized tetrahedron method [1].

REAL(8),INTENT(IN) :: bvec(3,3)

Reciprocal lattice vectors (arbitrary unit). Because they are used to choose the direction of tetrahedra, only their ratio is used.

INTEGER,INTENT(IN) :: nb

The number of bands.

INTEGER,INTENT(IN) :: nge(3)

Specify the \(k\)-grid for input orbital energy.

REAL(8),INTENT(IN) :: eig1(nb,nge(1),nge(2),nge(3))

The orbital energy measured from the Fermi energy ( \(\varepsilon_{\rm F} = 0\) ). Do the same with eig2.

REAL(8),INTENT(IN) :: eig2(nb,nge(1),nge(2),nge(3))

Another orbital energy. E.g. \(\varepsilon_{k + q}\) on a shifted grid.

INTEGER,INTENT(IN) :: ngw(3)

Specify the \(k\)-grid for output integration weights. You can make ngw \(\neq\) nge (See Appendix).

REAL(8),INTENT(OUT) :: wght(nb,nb,ngw(1),ngw(2),ngw(3))

The integration weights.

INTEGER,INTENT(IN),OPTIONAL :: comm

Optional argument. Communicators for MPI such as MPI_COMM_WORLD. Only for MPI / Hybrid parallelization. For C compiler without MPI, just pass NULL to omit this argment.

A part of DFPT calculation

\[\begin{align} \sum_{n n'} \int_{\rm BZ} \frac{d^3 k}{V_{\rm BZ}} \theta(\varepsilon_{\rm F} - \varepsilon_{n k}) \theta(\varepsilon_{n k} - \varepsilon'_{n' k}) X_{n n' k} \end{align}\]
CALL libtetrabz_dblstep(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,comm)

Parameters

INTEGER,INTENT(IN) :: ltetra

Specify the type of the tetrahedron method. 1 \(\cdots\) the linear tetrahedron method. 2 \(\cdots\) the optimized tetrahedron method [1].

REAL(8),INTENT(IN) :: bvec(3,3)

Reciprocal lattice vectors (arbitrary unit). Because they are used to choose the direction of tetrahedra, only their ratio is used.

INTEGER,INTENT(IN) :: nb

The number of bands.

INTEGER,INTENT(IN) :: nge(3)

Specify the \(k\)-grid for input orbital energy.

REAL(8),INTENT(IN) :: eig1(nb,nge(1),nge(2),nge(3))

The orbital energy measured from the Fermi energy ( \(\varepsilon_{\rm F} = 0\) ). Do the same with eig2.

REAL(8),INTENT(IN) :: eig2(nb,nge(1),nge(2),nge(3))

Another orbital energy. E.g. \(\varepsilon_{k + q}\) on a shifted grid.

INTEGER,INTENT(IN) :: ngw(3)

Specify the \(k\)-grid for output integration weights. You can make ngw \(\neq\) nge (See Appendix).

REAL(8),INTENT(OUT) :: wght(nb,nb,ngw(1),ngw(2),ngw(3))

The integration weights.

INTEGER,INTENT(IN),OPTIONAL :: comm

Optional argument. Communicators for MPI such as MPI_COMM_WORLD. Only for MPI / Hybrid parallelization. For C compiler without MPI, just pass NULL to omit this argment.

Static polarization function

\[\begin{align} \sum_{n n'} \int_{\rm BZ} \frac{d^3 k}{V_{\rm BZ}} \frac{\theta(\varepsilon_{\rm F} - \varepsilon_{n k}) \theta(\varepsilon'_{n' k} - \varepsilon_{\rm F})} {\varepsilon'_{n' k} - \varepsilon_{n k}} X_{n n' k} \end{align}\]
CALL libtetrabz_polstat(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,comm)

Parameters

INTEGER,INTENT(IN) :: ltetra

Specify the type of the tetrahedron method. 1 \(\cdots\) the linear tetrahedron method. 2 \(\cdots\) the optimized tetrahedron method [1].

REAL(8),INTENT(IN) :: bvec(3,3)

Reciprocal lattice vectors (arbitrary unit). Because they are used to choose the direction of tetrahedra, only their ratio is used.

INTEGER,INTENT(IN) :: nb

The number of bands.

INTEGER,INTENT(IN) :: nge(3)

Specify the \(k\)-grid for input orbital energy.

REAL(8),INTENT(IN) :: eig1(nb,nge(1),nge(2),nge(3))

The orbital energy measured from the Fermi energy ( \(\varepsilon_{\rm F} = 0\) ). Do the same with eig2.

REAL(8),INTENT(IN) :: eig2(nb,nge(1),nge(2),nge(3))

Another orbital energy. E.g. \(\varepsilon_{k + q}\) on a shifted grid.

INTEGER,INTENT(IN) :: ngw(3)

Specify the \(k\)-grid for output integration weights. You can make ngw \(\neq\) nge (See Appendix).

REAL(8),INTENT(OUT) :: wght(nb,nb,ngw(1),ngw(2),ngw(3))

The integration weights.

INTEGER,INTENT(IN),OPTIONAL :: comm

Optional argument. Communicators for MPI such as MPI_COMM_WORLD. Only for MPI / Hybrid parallelization. For C compiler without MPI, just pass NULL to omit this argment.

Phonon linewidth

\[\begin{align} \sum_{n n'} \int_{\rm BZ} \frac{d^3 k}{V_{\rm BZ}} \theta(\varepsilon_{\rm F} - \varepsilon_{n k}) \theta(\varepsilon'_{n' k} - \varepsilon_{\rm F}) \delta(\varepsilon'_{n' k} - \varepsilon_{n k} - \omega) X_{n n' k}(\omega) \end{align}\]
CALL libtetrabz_fermigr(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,ne,e0,comm)

Parameters

INTEGER,INTENT(IN) :: ltetra

Specify the type of the tetrahedron method. 1 \(\cdots\) the linear tetrahedron method. 2 \(\cdots\) the optimized tetrahedron method [1].

REAL(8),INTENT(IN) :: bvec(3,3)

Reciprocal lattice vectors (arbitrary unit). Because they are used to choose the direction of tetrahedra, only their ratio is used.

INTEGER,INTENT(IN) :: nb

The number of bands.

INTEGER,INTENT(IN) :: nge(3)

Specify the \(k\)-grid for input orbital energy.

REAL(8),INTENT(IN) :: eig1(nb,nge(1),nge(2),nge(3))

The orbital energy measured from the Fermi energy ( \(\varepsilon_{\rm F} = 0\) ). Do the same with eig2.

REAL(8),INTENT(IN) :: eig2(nb,nge(1),nge(2),nge(3))

Another orbital energy. E.g. \(\varepsilon_{k + q}\) on a shifted grid.

INTEGER,INTENT(IN) :: ngw(3)

Specify the \(k\)-grid for output integration weights. You can make ngw \(\neq\) nge (See Appendix).

REAL(8),INTENT(OUT) :: wght(ne,nb,nb,ngw(1),ngw(2),ngw(3))

The integration weights.

INTEGER,INTENT(IN) :: ne

The number of branches of the phonon.

REAL(8),INTENT(IN) :: e0(ne)

Phonon frequencies.

INTEGER,INTENT(IN),OPTIONAL :: comm

Optional argument. Communicators for MPI such as MPI_COMM_WORLD. Only for MPI / Hybrid parallelization. For C compiler without MPI, just pass NULL to omit this argment.

Polarization function (complex frequency)

\[\begin{align} \sum_{n n'} \int_{\rm BZ} \frac{d^3 k}{V_{\rm BZ}} \frac{\theta(\varepsilon_{\rm F} - \varepsilon_{n k}) \theta(\varepsilon'_{n' k} - \varepsilon_{\rm F})} {\varepsilon'_{n' k} - \varepsilon_{n k} + i \omega} X_{n n' k}(\omega) \end{align}\]
CALL libtetrabz_polcmplx(ltetra,bvec,nb,nge,eig1,eig2,ngw,wght,ne,e0,comm)

Parameters

INTEGER,INTENT(IN) :: ltetra

Specify the type of the tetrahedron method. 1 \(\cdots\) the linear tetrahedron method. 2 \(\cdots\) the optimized tetrahedron method [1].

REAL(8),INTENT(IN) :: bvec(3,3)

Reciprocal lattice vectors (arbitrary unit). Because they are used to choose the direction of tetrahedra, only their ratio is used.

INTEGER,INTENT(IN) :: nb

The number of bands.

INTEGER,INTENT(IN) :: nge(3)

Specify the \(k\)-grid for input orbital energy.

REAL(8),INTENT(IN) :: eig1(nb,nge(1),nge(2),nge(3))

The orbital energy measured from the Fermi energy ( \(\varepsilon_{\rm F} = 0\) ). Do the same with eig2.

REAL(8),INTENT(IN) :: eig2(nb,nge(1),nge(2),nge(3))

Another orbital energy. E.g. \(\varepsilon_{k + q}\) on a shifted grid.

INTEGER,INTENT(IN) :: ngw(3)

Specify the \(k\)-grid for output integration weights. You can make ngw \(\neq\) nge (See Appendix).

COMPLEX(8),INTENT(OUT) :: wght(ne,nb,nb,ngw(1),ngw(2),ngw(3))

The integration weights.

INTEGER,INTENT(IN) :: ne

The number of imaginary frequencies where polarization functions are calculated.

COMPLEX(8),INTENT(IN) :: e0(ne)

Complex frequencies where polarization functions are calculated.

INTEGER,INTENT(IN),OPTIONAL :: comm

Optional argument. Communicators for MPI such as MPI_COMM_WORLD. Only for MPI / Hybrid parallelization. For C compiler without MPI, just pass NULL to omit this argment.