Post Analysis in NanoDCAL

In nanodcal, wave function of an electronic state is expressed as

\[\psi_{i}(\textbf{r})=\sum_{\nu}c_{\nu}^{i}\zeta_{\nu}(\textbf{r})\ ,\]

where the subscript \(i\) labels a complete set of quantum numbers which uniquely identify the electronic state; \(\zeta_{\nu}(\textbf{r})\) is the atomic orbital basis and the index \(\nu\) accounts for all atomic sites \(\textbf{R}_{I_{\nu}}\), multiple-zeta index \(n_{\nu}\), angular momentum quantum numbers \(l_{\nu},m_{\nu}\), and spin quantum number \(m_{s_{\mu}}\) for spin polarized system. The coefficients \(c_{\nu}^{i}\) and the corresponding energy of the state is calculated by solving the general eigen value equation

\[H_{\mu\nu}c_{\nu}^{i}=\epsilon_{i}S_{\mu\nu}c_{\nu}^{i}\ .\]

The eigenvalue \(\{\epsilon_{i}\}\) of an isolated system is discrete and forms a continuous band structure for a periodic bulk system. In general, the spectrum of an open system consists of both discrete levels and continuous bands. A state corresponding to a discrete level is normally a bound state which is localized; a state in the continuum is normally a scattering state which originally comes from an electronic channel of a device electrode and has a certain probability to be scattered into other channels of the electrodes.

This section discusses some electronic structure related post analysis of nanodcal.

Band structure

For a periodic system with primitive cell vectors (a\(_{1}\), a\(_{2}\), a\(_{3}\)), its band structure, \(\epsilon_{b}(\textbf{k})\), is calculated in nanodcal by solving the general eigen value problem

(42)\[H_{\mu\nu}(\textbf{k})c_{\nu}^{b}(\textbf{k})=\epsilon_{b}(\textbf{k})S_{\mu\nu} (\textbf{k})c_{\nu}^{b}(\textbf{k})\]

where \(b\) labels the band and

\[H_{\mu\nu}(\textbf{k})\equiv\sum_{D}e^{i\textbf{ R}_{D}\cdot\textbf{k}}H_{\mu\nu}^{D},\]
\[S_{\mu\nu}(\textbf{k})\equiv\sum_{D}e^{i\textbf{R}_{D} \cdot\textbf{k}}S_{\mu\nu}^{D}\ .\]

The index \(D=(d_{1},d_{2},d_{3})\) and \(d_{1}\), \(d_{2}\), \(d_{3}\) are integers. \(D\) labels an image cell with a displacement \(\textbf{R}_{D}\equiv\sum_{j}d_{j}\textbf{a}_{j}\), and \(H_{\mu\nu}^{D}\) (\(S_{\mu\nu}^{D}\)) is the Hamiltonian (overlapping) matrix element between an orbital \(\mu\) within the primitive cell and an orbital \(\nu\) within the image cell \(D\). In nanodcal the band structure is normally calculated along a line segment (or several such line segments) between two high symmetrical points \(P_{1}\) and \(P_{2}\) in k space,

\[\epsilon_{b}(\lambda)=\epsilon_{b}(\lambda\textbf{k} (P_{1})+(1-\lambda)\textbf{k}(P_{2}))\ ,\:\:\lambda\in[0,1].\]

Together with \(\epsilon_{b}(\textbf{k})\), nanodcal also gives out the contribution weight \(w_{\mu}^{b}(\textbf{k})\) of every orbital \(\mu\) to the band \(b\) at the point k as


where there is no summation on the \(\mu\).

Complex-band structure

In the last subsection, one calculates a band structure \(\epsilon_{b}(\textbf{k})\) by solving Eq. (42) for many different values of real numbers \(k_1,k_2,k_3\). On the other hand, the same band structure can be obtained by fixing the value of \(\epsilon_{b}\) and determining \(\textbf{k}\) that gives this \(\epsilon_b\) value. Importantly, not only real \(\textbf{k}\) can give a real \(\epsilon_b\), but also some complex \(\textbf{k}\) that can give the same real \(\epsilon_b\). This is the complex band structure which is useful for analyzing transport features in the tunneling regime.

Consider a periodic system with primitive cell vectors (a\(_{1}\), a\(_{2}\), a\(_{3}\)), the corresponding \(k\)-space cell vectors (b\(_{1}\), b\(_{2}\), b\(_{3}\)) can be defined by

\[\textbf{b}_{i}\cdot\textbf{a}_{j}=2\pi\delta_{ij}\ ,\]

and a \(k\)-space vector k can be expressed as

\[\textbf{k}=k_{j}\textbf{b}_{j}\ .\]

Extending one of the real numbers \(k_{1}\), \(k_{2}\), and \(k_{3}\) into the complex plane, say \(k_{3}\), nanodcal calculates a complex-band structure, \(k_{3}(\epsilon;k_{1},k_{2})\) which is a complex valued function of real energy \(\epsilon\) having the \(k_{1}\) and \(k_{2}\) as given parameters. The complex band structure is obtained by solving the following equation for \(k_{3}\),

\[H_{\mu\nu}(k_{3}^{b};k_{1},k_{2})c_{\nu}^{b}(k_{1},k_{2}) =\epsilon S_{\mu\nu}(k_{3}^{b};k_{1},k_{2})c_{\nu}^{b}(k_{1},k_{2})\ ,\]

where \(b\) labels the complex band and

\[H_{\mu\nu}(k_{3}^{b};k_{1},k_{2})\equiv\sum_{D}e^{i\textbf{ R}_{D}\cdot(k_{1}\textbf{b}_{1}+k_{2}\textbf{ b}_{2}+k_{3}^{b}\textbf{b}_{3})}H_{\mu\nu}^{D},\]
\[S_{\mu\nu}(k_{3}^{b};k_{1},k_{2})\equiv\sum_{D}e^{i\textbf{R}_{D} \cdot(k_{1}\textbf{b}_{1}+k_{2}\textbf{b}_{2}+k_{3}^{b} \textbf{b}_{3})}S_{\mu\nu}^{D}.\]

The \(D=(d_{1},d_{2},d_{3})\) with \(d_{1}\), \(d_{2}\) and \(d_{3}\) integers, labels an image cell with a displacement \(\textbf{R}_{D}\equiv\sum_{j}d_{j}\textbf{a}_{j}\), and \(H_{\mu\nu}^{D}\) (\(S_{\mu\nu}^{D}\)) is the Hamiltonian (overlapping) matrix element between an orbital \(\mu\) within the primitive cell and an orbital \(\nu\) within the image cell \(D\).

Scattering states

In nanodcal, the scattering states are defined for energy ranges between the left and right chemical potentials \(\mu_l\) and \(\mu_r\), respectively. To solve for these states, at a given energy \(E\), an inverse energy band structure problem [TGW01] is solved first. We then group all states as left and right propagating states depending on their group velocity. Let \(k\) denote electron momentum in the left lead and \(p\) in the right lead. For a scattering state \(\Psi^{k_n}_l\) coming from the left (\(l\)) lead, it starts as a right propagating state \(\Phi^{k_n}_l\) and gets reflected back as a left propagating state \(\phi^{k_m}_l\) with reflection amplitude \(r^{k_m,k_n}\) in the left lead, and transmitted into the right lead as a right propagating state \(\phi^{p_m}_r\) with transmission amplitude \(t^{p_m,k_n}_R\). In these quantities, the indices \(m, n\) label the sub-bands. In the numerics, the scattering states are represented as a linear combination of atomic orbital inside the device region. These allow one to write, for example, a scattering state going from left to right, as:

\[\begin{split}\Psi^{k_n}=\left \{ \begin{array}{l l} \Phi^{k_n}_l + \phi^{k_m}_l r^{k_m,k_n} & \quad \textrm{inside left lead}\\ \psi^{k_n}_d & \quad \textrm{inside device}\\ \phi^{p_m}_r t^{p_m,k_n} & \quad \textrm{inside right lead} \end{array} \right.\end{split}\]

A scattering state in the right lead can be written in a similar fashion. After the NEGF-DFT self-consistent iteration is completed, one applies the obtained device Hamiltonian to the scattering states of the above form and the reflection and transmission amplitudes can be calculated by matrix techniques. The modular square of these amplitudes give reflection and transmission coefficients [TGW01, T01].

Density of states, Density of scattering states

The total number of states calculated in nanodcal can be written as

\[N=\sum_{\mu\in\textbf{\{C\}}}\sum_{\nu,i} \frac{1}{2}[(c_{\mu}^{i*}c_{\nu}^{i}S_{\nu\mu})+ (S_{\mu\nu}c_{\nu}^{i*}c_{\mu}^{i})]\ ,\]

where \(\textbf{\{C\}}\) denotes the central cell: namely the central unit cell for closed systems or the central scattering region for open systems. The corresponding density of states (DOS) can be written as

(43)\[n(\epsilon)=\sum_{\mu\in\textbf{\{C\}}} \sum_{\nu,i}\delta(\epsilon-\epsilon_{i})\frac{1}{2} [(c_{\mu}^{i*}c_{\nu}^{i}S_{\nu\mu})+(S_{\mu\nu}c_{\nu}^{i*} c_{\mu}^{i})]\ .\]

Since the energy spectrum \(\{\epsilon_{i}\}\) in general consists of both discrete points and continuous bands, the summation over \(i\) should go over all the electronic states and must be understood as both summation over those discrete points and integration over those continuous bands.

The DOS \(n(\epsilon)\) can be decomposed in the orbital space spanned by \(\mu\) into contributions of various orbital types \(\kappa\),

\[n(\epsilon)=\sum_{\kappa} n(\epsilon;\kappa),\]


(44)\[n(\epsilon; \kappa) \equiv \sum_{\mu \in \{\kappa \} \cap \textbf{\{C\}}} \sum_{\nu,i}\delta(\epsilon-\epsilon_{i})\frac{1}{2}[c_{\mu}^{i*}c_{\nu}^{i} S_{\nu\mu}+S_{\mu\nu}c_{\nu}^{i*}c_{\mu}^{i}]\ ,\]

where an orbital type labeled by \(\kappa\) can be specified flexibly in nanodcal. Namely \(\kappa\) can be the angular momentum quantum number \(l_{\mu}\) of the orbital, or the atomic sites, or the sequential numbers of the orbital. The \(n(\epsilon;\kappa)\) is usually called projected density of states.

The \(n(\epsilon)\) and \(n(\epsilon;\kappa)\) can be further decomposed into contributions from different real space local regions \(\Omega\),

\[n(\epsilon)=\sum_{\Omega} n(\epsilon; \Omega)\ ,\]
\[n(\epsilon; \kappa)=\sum_{\Omega} n(\epsilon; \kappa,\Omega)\ .\]

Quantities \(n(\epsilon; \Omega)\) and \(n(\epsilon; \kappa,\Omega)\) are calculated as:

\[n(\epsilon; \Omega) \equiv \int_{\Omega} d\textbf{r} \sum_{\mu\in\textbf{\{C\}}}\sum_{\nu,i}\delta(\epsilon-\epsilon_{i}) \frac{1}{2} [c_{\mu}^{i*}\zeta_{\mu}(\textbf{r})c_{\nu}^{i}\zeta_{\nu} (\textbf{r})+c_{\nu}^{i*}\zeta_{\nu}(\textbf{r})c_{\mu}^{i}\zeta_{\mu} (\textbf{r})]\ ,\]
\[\begin{split}\begin{aligned} & n(\epsilon; & \kappa, \Omega)\equiv\frac{1}{2}\sum_{\mu\in\textbf{\{}\kappa\textbf{\}}\cap\textbf{\{C\}}} \int_{\Omega}d\textbf{r} \sum_{\nu}\sum_{i}\delta(\epsilon-\epsilon_{i})\nonumber \\ & & \times [c_{\mu}^{i*}\zeta_{\mu}(\textbf{r})c_{\nu}^{i}\zeta_{\nu}(\textbf{r})+ c_{\nu}^{i*}\zeta_{\nu}(\textbf{r})c_{\mu}^{i}\zeta_{\mu}(\textbf{r})],\end{aligned}\end{split}\]

where a local region \(\Omega\), with a shape of parallelepiped, can be defined flexibly in nanodcal by its position and three cell vectors. The \(n(\epsilon;\Omega)\) is usually called the local density of states, and \(n(\epsilon;\kappa,\Omega)\) may be called projected-local density of states.

Since the scattering states form a subset of the total electronic states, density of scattering states is just a portion of the (total) density of states which only counts the scattering states. The density of scattering states (DOSS) calculated in nanodcal, \(n_{sc}(\epsilon)\), or partial density of scattering states associated with a particular lead \(\alpha\), \(n_{sc}(\epsilon; \alpha)\), can be written as

\[n_{sc}(\epsilon)=\sum_{\alpha}n_{sc}(\epsilon;\alpha )\ ,\]


\[n_{sc}(\epsilon; \alpha)=\frac{1}{2}\sum_{i\in{\mathbb{S}}_{\alpha}} \sum_{\mu\in\textbf{\{C\}}} \sum_{\nu}\delta(\epsilon-\epsilon_{i})[(c_{\mu}^{i*}c_{\nu}^{i}S_{\nu\mu}) +(S_{\mu\nu}c_{\nu}^{i*}c_{\mu}^{i})],\]

where \({\mathbb{S}}_{\alpha}\) denotes the scattering states that originate from lead \(\alpha\).

Similar to the DOS, the DOSS \(n_{sc}(\epsilon;\alpha)\) can be decomposed, in the space spanned by \(\mu\), into the contributions of various orbital types \(\kappa\),

\[n_{sc}(\epsilon;\alpha;\kappa)=\frac{1}{2}\sum_{i\in{\mathbb{S}}_{\alpha}} \sum_{\mu\in\textbf{\{}\kappa\textbf{\}}\cap\textbf{\{C\}}} \sum_{\nu}\delta(\epsilon-\epsilon_{i}) [c_{\mu}^{i*}c_{\nu}^{i}S_{\nu\mu}+S_{\mu\nu}c_{\nu}^{i*}c_{\mu}^{i}]\ .\]

The quantity \(n_{sc}(\epsilon;\alpha;\kappa)\) can be further decomposed, in real space, into contributions from different local regions \(\Omega\),

\[\begin{split}\begin{aligned} & n_{sc}(\epsilon; & \alpha;\kappa,\Omega)=\frac{1}{2}\sum_{i\in{\mathbb{S}}_{\alpha}} \sum_{\mu\in\textbf{\{ }\kappa\textbf{\}}\cap\textbf{\{C\}}}\int_{\Omega}d\textbf{r}\nonumber \\ & & \sum_{\nu}\delta(\epsilon-\epsilon_{i})[c_{\mu}^{i*}\zeta_{\mu}(\textbf{r})c_{\nu}^{i} \zeta_{\nu}(\textbf{r})+c_{\nu}^{i*}\zeta_{\nu}(\textbf{r})c_{\mu}^{i}\zeta_{\mu}(\textbf{r})]\ .\end{aligned}\end{split}\]


The charge can be calculated through density matrix and real space charge density, respectively. In nanodcal, the total number of electrons calculated through the density matrix, called orbital charge \(N_o\), can be written as

\[N_o=\frac{1}{2}\sum_{\mu\in\textbf{\{C\}}}[(S\rho)_{\mu\mu}+(\rho S)_{\mu\mu}],\]

where the C means the central unit cell for closed systems, or the central scattering region for open devices. \(N_o\) can be decomposed, in orbital space spanned by \(\mu\), into contributions of various orbital types,

\[N_o=\sum_{\kappa}N_o(\kappa)\ \ ,\]


\[N_o(\kappa)\equiv\frac{1}{2}\sum_{\mu\in\textbf{\{}\kappa\textbf{\}}\cap\textbf{\{C\}}} [(S\rho)_{\mu\mu}+(\rho S)_{\mu\mu}]\ .\]

In nanodcal, the total number of electrons calculated through the real space charge density, called grid charge \(N_r\), can be written as

\[N_{r}={\displaystyle \int_{\textbf{C}}}d\textbf{r}\,\rho(\textbf{r})\ .\]

\(N_r\) can be decomposed into contributions from the grid points in the central cell,

\[N_{r} = \sum_{i} N_{r}(i)\]

where \(i\) is the real space grid label and

\[N_r(i)=\int_{i}d\textbf{r}\,\rho(\textbf{r})\ .\]

For a closed system, \(N_{o}\) and \(N_{r}\) are in principle equal to each other and their difference can be used as an indicator of the precision of the numerical computation. However, for an open system they are in principle not equal to each other by their definitions, hence the difference should not be viewed as an indicator of numerical error. For a spin polarized system, both the \(N_{o}\) and the \(N_{r}\) can be further analyzed to obtain its spin polarization vector in spin space. The direction of the spin polarization vector is along the spin polarized direction and its amplitude can be written as

\[\eta=\frac{N_{p}-N_{ap}}{N_{p}+N_{ap}}\ ,\]

where \(N_{p}\) (\(N_{ap}\)) is the number of electrons whose spin are parallel (or anti-parallel) to the spin polarized direction.

Total energy

The total energy calculated in nanodcal is written as

\[E_{{ tot}}=E_{{ bs}}+\delta E_{\delta{ H}}+\delta U_{{ xc}}+E_{{ sr}}\ ,\]


\[E_{{ bs}}\equiv\sum_{\mu\in\textbf{\{C\}}}\frac{1}{2}[(H\rho)_{\mu\mu}+(\rho H)_{\mu\mu}]\ ,\]
\[\delta E_{\delta{ H}}\equiv-\frac{1}{2}{\displaystyle \int_{\textbf{C}}}d\textbf{r}\, V_{\delta H}(\textbf{r})\left[\rho(\textbf{r})+\rho^{NA}(\textbf{r})\right]\ ,\]
\[\delta U_{xc}\equiv{\displaystyle \int_{\textbf{C}}}d\textbf{r}\,\{\epsilon_{xc}(\textbf{r}) [\rho(\textbf{r})+\rho_{c}(\textbf{r})]-V_{xc}(\textbf{r})\rho(\textbf{r})\}\ ,\]
\[E_{{ sr}}\equiv{\displaystyle \sum_{I\in\textbf{\{C\}}}}(E_{{ sr}}^{I}+\frac{1}{2}{\displaystyle \sum_{J(\neq I)}}V_{sr}^{I,J})\ ,\]
\[E_{{ sr}}^{I}=-\frac{1}{2}{\displaystyle \int\int}d\textbf{r}d\textbf{r}'\, \frac{\rho^{{ NA},I}(\textbf{r})\rho^{{ NA},I}(\textbf{r}')}{|\textbf{r}-\textbf{r}'|}\ ,\]
\[V_{sr}^{I,J}=\frac{Z_{I}Z_{J}}{|\textbf{R}_{I}-\textbf{R}_{J}|}- {\displaystyle \int\int}d\textbf{r}d\textbf{r}'\, \frac{\rho^{{ NA},I}(\textbf{r}-\textbf{R}_{I}) \rho^{{ NA},J}(\textbf{r}-\textbf{R}_{J})}{|\textbf{r}-\textbf{r}'|}\ ,\]

where the indices \(I,J\) label the atoms, \(\rho_{c}(\textbf{r})\) is the partial-core charge density [LRC82].

Transmission, Transmission channel

Considering an open system with several leads, the total probability of the electrons at energy \(\epsilon\), coming from all channels of a lead \(\alpha\), and scattered into channels of another lead \(\beta\) is described as a transmission coefficient \(T_{\alpha\beta}(\epsilon)\) , which can be calculated through Green’s function,

(45)\[T_{\alpha\beta}(\epsilon)=tr[G^{r}(\epsilon)\Gamma_{\alpha}(\epsilon)G^{a}(\epsilon) \Gamma_{\beta}(\epsilon)],\]


\[\Gamma_{\alpha}(\epsilon)\equiv i[\Sigma_{\alpha}^{r}(\epsilon)-\Sigma_{\alpha}^{a}(\epsilon)],\]

where \(\Sigma_{\alpha}^{r}(\epsilon)\) (\(\Sigma_{\alpha}^{a}(\epsilon)\)) is the retarded (advanced) self-energy of lead \(\alpha\) and \(tr[\cdots]\) stands for the trace of \([\cdots]\).

The transmission coefficient \(T_{\alpha\beta}(\epsilon)\) can also be calculated through elements of scattering matrix, where \(p\) (\(q\)) labels the channel of lead \(\alpha\) (\(\beta\)). When the system is perfect, electrons can move freely in the system. Since there is no scattering in perfect systems, \(s_{\beta q,\alpha p}=\delta_{qp}\), and the transmission coefficient \(T_{\alpha\beta}(\epsilon)\) is reduced to the number of (open) transmission channels of lead \(\alpha\). In nanodcal, the number of transmission channels is calculated by counting the number of the Bloch waves existed in the system which propagate along a given direction.

Current and Conductance

Considering an open system with several leads, the spin non-polarized electric current in lead \(\beta\) can be calculated by the Landauer formula,

\[\begin{split}\begin{aligned} I_{\beta} & = & \frac{2e^{2}}{h}\sum_{\alpha}\int d\epsilon [f_{\alpha}(\epsilon)T_{\alpha\beta}(\epsilon)-f_{\beta}(\epsilon)T_{\beta\alpha}(\epsilon)]\nonumber \\ & = & \frac{2e^{2}}{h}\sum_{\alpha}\int d\epsilon [f_{\alpha}(\epsilon)-f_{\beta}(\epsilon)]T_{\alpha\beta}(\epsilon)\ ,\end{aligned}\end{split}\]

where \(f_{\alpha}(\epsilon)\) is the distribution function of the electrons in lead \(\alpha\), and the factor of 2 comes from the spin degeneracy. The above equation can be re-written as


with the conductance \(G_{\alpha\beta}\) defined as

(46)\[G_{\alpha\beta}=\frac{2e^{2}}{h}\int d\epsilon\,\frac{f_{\alpha}(\epsilon)-f_{\beta}(\epsilon)} {V_{\beta}-V_{\alpha}}T_{\alpha\beta}(\epsilon)\ ,\]

where \(V_{\beta}\) is the bias voltage applied on lead \(\beta\).

For spin polarized transport, the spin current (spin polarized charge current) is calculated from spin resolved transmission coefficient, \(T_{\sigma}\) where \(\sigma\equiv \uparrow,\downarrow\) is the spin index. \(T_\sigma\) is still given by Eq. (45) but the matrices on the right hand side are spin resolved (the linear size of the matrix is therefore doubled). The spin-current (spin polarized charge current) for a two probe system is obtained by the following formula,

\[I_\sigma =\frac{e}{h}\int d\epsilon T_{\sigma}(\epsilon,V_b)[ f_L(\epsilon) -f_R(\epsilon)]\ ,\]

where we explicitly restored the bias voltage \(V_b\) dependence of the transmission coefficient. The integration is over the bias window between \(\mu_L\) and \(\mu_R\) which are the electrochemical potentials of the left/right leads. From the spin current, the spin resolved conductance can be obtained similar to Eq. (46).

The total charge current is given by

\[I =I_{\uparrow}+I_{\downarrow}\ .\]

AC Conductance

When an AC bias applied, the AC conductance \(\tilde{G}\) calculated in nanodcal can be written as [WWG99] (in atomic unit)

(47)\[\tilde{G}_{\alpha\beta}(\omega)=G_{\alpha\beta}^{c}(\omega)\,+ \, w_{\alpha}(\omega)G_{\beta}^{d}(\omega)\ ,\]

where the subscripts \(\alpha\), \(\beta\) denotes the leads and \(\omega\) is the AC frequency.

The first term \(G^{c}\) of Eq. (47) is the dynamic conductance due to particle current alone,

\[\begin{split}\begin{aligned} G_{\alpha\beta}^{c}(\omega) & = & -\int\frac{d\epsilon}{2\pi}\, \frac{f(\epsilon)-f(\epsilon+\omega)}{\omega} tr\left[-G^{r}(\epsilon+\omega)\Gamma_{\beta}G^{a}(\epsilon)\Gamma_{\alpha}\right. \nonumber \\ && \left. +G^{r}(\epsilon+\omega)\Gamma G^{a}(\epsilon) \Gamma_{\alpha}\delta_{\alpha\beta}-i\omega G^{r}(\epsilon+\omega) G^{a}(\epsilon)\Gamma_{\alpha}\delta_{\alpha\beta}\right].\end{aligned}\end{split}\]

The second term \(w_{\alpha}G_{\beta}^{d}\) of Eq. (47) is the contribution from the displacement current due to electrodynamics where the quantities

\[G_{\beta}^{d}(\omega)\equiv i\omega\int\frac{d\epsilon}{2\pi}\, \frac{f(\epsilon)-f(\epsilon+\omega)}{\omega} tr\left[G^{r}(\epsilon+\omega)\Gamma_{\beta}G^{a}(\epsilon)\right]\ ,\]


\[w_\alpha\ \equiv\ -\frac{\sum_\gamma G_{\alpha\gamma}^c}{\sum_\gamma G_{\gamma}^d}\ .\]

Atomic Force and Stress

The expression for the atomic force within our DFT-NEGF formalism is as follows:

\[\begin{split}{\bf F}_{K} &=& -\sum_{\mu \nu} \left[\rho_{\nu \mu} \left(\frac{\partial T_{\mu \nu}}{\partial {\bf R}_{K}} + \frac{\partial V^{NL}_{\mu \nu}}{\partial {\bf R}_{K}} \right) - \epsilon_{\nu \mu}\frac{\partial S_{\mu \nu}}{\partial {\bf R}_{K}} \right] - \sum_{\langle IJ\rangle} \frac{\partial U_{SR}^{IJ}}{\partial {\bf R}_{K}}\nonumber\\ && -\int d{\bf r}\frac{\partial V^{NA}({\bf r})}{\partial {\bf R}_{K}} \rho({\bf r})-\int d{\bf r}\frac{\partial \rho_{pc}({\bf r})}{\partial {\bf R}_{K}} V_{XC}({\bf r}) + \int d{\bf r}\frac{\partial \rho^{NA} ({\bf r})}{\partial {\bf R}_{K}} V^{\delta}_{H}({\bf r}) \nonumber\\ && - \int d{\bf r}\left[\sum_{\mu\nu}\rho_{\mu\nu}\bigg(\frac{\partial \phi_{\mu}({\bf r})}{\partial {\bf R}_{K}}\phi_{\nu}({\bf r}) + \phi_{\mu}({\bf r})\frac{\partial \phi_{\nu}({\bf r})}{\partial {\bf R}_{K}} \bigg)\right]V_{eff}({\bf r})\nonumber\\\end{split}\]

The \(T_{\mu\nu}\), etc. are kinetic matrix etc.

In a supercell calculation, the stress tensor has the form

\[\begin{split}\sigma &=& -\sum_{\mu \nu} \left[\rho_{\nu \mu} \left(\frac{\partial T_{\mu \nu}}{\partial \epsilon} + \frac{\partial V^{NL}_{\mu \nu}}{\partial \epsilon} \right) - \epsilon_{\nu \mu}\frac{\partial S_{\mu \nu}}{\partial \epsilon} \right] - \sum_{\langle IJ\rangle} \frac{\partial U_{SR}^{IJ}}{\partial \epsilon}\nonumber\\ && - \frac{\partial}{\partial\epsilon}\sum_{\langle IJ\rangle}\int d{\bf r} \rho_{I}^{NA}({\bf r})V_{J}^{NA}({\bf r}) - \int d{\bf r}\delta\rho({\bf r}) \textit{IFFT}\bigg(\frac{{\bf k k}}{k^4}\delta\rho^{k}\bigg) \nonumber\\ && - {\bf 1}\bigg[\int d{\bf r}\bigg(V^{NA}({\bf r})+\frac{1}{2}V_{\delta H} ({\bf r})\bigg)\delta\rho({\bf r}) + E_{XC}\bigg] -\int d{\bf r}\frac{\partial \rho_{pc}({\bf r})}{\partial \epsilon} V_{XC}({\bf r})\nonumber\\ && + \int d{\bf r}\frac{\partial \rho^{NA}({\bf r})}{\partial \epsilon} \bigg(V^{NA}({\bf r})+V_{\delta H}({\bf r})\bigg) - \int d{\bf r}\frac{\partial V^{NA}({\bf r})}{\partial \epsilon} \delta\rho({\bf r})\nonumber\\ && - \int d{\bf r}\left[\sum_{\mu\nu}\rho_{\mu\nu}\bigg(({\bf r}-{\bf R}_{\mu})\nabla\phi_{\mu}({\bf r})\phi_{\nu}({\bf r})+\phi_{\mu}({\bf r})({\bf r}-{\bf R}_{\nu})\nabla\phi_{\nu}({\bf r}) \bigg)\right]V_{eff}({\bf r})\nonumber\\ .. && +\int d{\bf r}\left[\sum_{\mu}2\rho_{\mu}^{occ}\phi_{\mu}({\bf r})({\bf r}-{\bf R}_{\mu})\nabla\phi_{\mu}({\bf r}) \right]\bigg(V^{NA}({\bf r})+V_{\delta H}({\bf r})\bigg)\end{split}\]

Taylor, J., Guo, H., & Wang, J. (2001). Ab initio modeling of quantum transport properties of molecular electronic devices. Physical Review B, 63(24), 245407.

  1. Datta, Electronic Transport in Mesoscopic Systems, (Cambridge University Press, New York, 1995).


Louie, S. G. S., Froyen, S., & Cohen, M. M. L. (1982). Nonlinear ionic pseudopotentials in spin-density-functional calculations. Physical Review B, 26(4), 1738–1742.


Wang, B., Wang, J., & Guo, H. (1999). Current Partition: A Nonequilibrium Green’s Function Approach. <10.1103/PhysRevLett.82.398> Physical Review Letters, 82(2), 398.


The calculation details for scattering states can be found in J. Taylor, Ph.D Thesis, McGill University (2001).