Ab initio nonequilibrium quantum transport and forces with the real-space projector augmented wave method

We present an efficient implementation of a nonequilibrium Green's function method for self-consistent calculations of electron transport and forces in nanostructured materials. The electronic structure is described at the level of density functional theory using the projector augmented wave method to describe the ionic cores and an atomic orbital basis set for the valence electrons. External bias and gate voltages are treated in a self-consistent manner and the Poisson equation with appropriate boundary conditions is solved in real space. Contour integration of the Green's function and parallelization over $k$ points and real space makes the code highly efficient and applicable to systems containing several hundreds of atoms. The method is applied to a number of different systems, demonstrating the effects of bias and gate voltages, multiterminal setups, nonequilibrium forces, and spin transport.

Figure 1

A scattering region including the nanostructure of interest (e.g., a molecule) and part of the electrode atoms is sandwiched between two semi-infinite electrodes. Periodic boundary conditions are used in the $x,y$ directions and open boundary conditions in the $z$ direction. The electron potential in the electrodes is periodic and can be obtained from a ground-state DFT calculation employing periodic boundary conditions in all directions. The Hartree potential at the scattering region boundary, which is used as boundary condition for the Poisson equation, is also obtained from the electrode calculation. The whole system can be subject to an external bias or gate voltage, and the electronic structure of the scattering region is calculated self-consistently in the presence of such external fields.

Figure 2

The contour used for the Green's function integral in the complex energy plane. The coordinates of the indicated points are as follows: A $(E_{\min },0)$, where $E_{\min }$ is less than all the eigenenergies of the system, which is usually taken as ${\mu }_L-100eV$ since we only calculate the valence electrons states; B $({\mu }_L+m{k}_{B}T,\Delta )$, where $m$ satisfies $e^{-m}\approx 0$ (a typical value of $m$ is 10), $k_B$ is the Boltzmann constant, $T$ is the electron temperature, $\Delta$ is normally between ${\epsilon }_n$ and ${\epsilon }_{n+1}$, where ${\epsilon }_i=(2i-1)\pi$ (with $i$ a positive integer) is a pole of the Fermi-Dirac distribution function, thus the singulars below $\Delta$ should be counted when summing up the residues; C $({\mu }_L-m{k}_{B}T,\Delta )$; D $({\mu }_R+m{k}_{B}T,\Delta )$; E $({\mu }_L-m{k}_{B}T,\eta )$, $\eta$ is infinitesimal to avoid the inversion divergence; F $({\mu }_R+m{k}_{B}T,\eta )$.

Figure 3

Schematic of the matrix blocks. (a) The Hamiltonian or overlap matrix: The blue and green blocks represent the on-site and coupling submatrices for the different quenching layers, respectively. (b) The Green's function matrix when evaluating the density, $\text{DOS,}$ or transmission: The red and orange blocks represent the submatrices needed to calculate the density matrix or $\text{DOS}$, and pink blocks are for the transmission coefficient. (c) The red blocks in the retarded Green's function matrix are necessary when calculating the Keldysh Green's function.

Figure 4

Upper panel: Capacitor model with the electrodes made of bcc sodium with the voltage drop along the (100) direction. The vacuum distance between the two electrodes is $8\text{Å}$, and the size of the unit cell in the transverse directions is $12.7\text{Å}\times 12.7\text{Å}$. The rectangle represents the scattering region. Lower panel: The nonequilibrium part of the electron density (solid blue) and the induced effective potential (dashed red) under a bias voltage of 5 V, with the zero-bias values as the reference. The calculated values are averaged over the transverse plane.

Figure 5

The bias voltage effect on a nitrogen molecule between two gold electrodes. Upper panel: The atomic structure of the system. The arrows represent the directions of the atomic forces generated by the bias voltage. Middle panel: $\text{PDOS}$(partial density of states) of the nitrogen molecule at 0- and 1-V bias voltage. The Fermi level is located at zero and the red dashed lines show the location of the bias window in the 1-V case. Lower panel: The magnitude of the nonequilibrium atomic forces as a function of bias voltage.

Figure 6

Effect of electrostatic gating of a benzene-1,4-dithiolate molecule between two gold fcc(111) electrodes. Upper panel: The atomic structure of the system. Middle panel: The applied external potential. Lower panel: The effective additional potential after self-consistency.

Figure 7

Gate tuning the conductance of a molecular transistor. Upper panel: The Au-anthraquinone-Au structure. Middle panel: The gate effect on the potential. Lower panel: The transmission coefficient at 0 and $-2$ V gate voltages.

Figure 8

Spin transport in a zigzag graphene nanoribbon (ZGNR) bridge connecting to two semi-infinite graphene flakes. Upper panel: The atomic structure of the graphene/ZGNR/graphene system. Middle panel: The $\text{PDOS}$ of a C atom at the ribbon's zigzag edge under zero bias; the red and blue solid lines represent the spin-up and spin-down $\text{PDOS}$ of the C atom at the center of the zigzag edge (marked with X), the red and blue dashed lines represent the spin-up and spin-down $\text{PDOS}$ of the C atom next to the previous one at the zigzag edge (marked with Y), and the green dashed line is the Fermi level. Lower panel: The $\text{PDOS}$ of a C atom at the ribbon's zigzag edge under bias voltage $V=2.4\text{V}$.

Figure 9

Calculated current (blue squares) and edge magnetic moment per C atom (red circles) as a function of bias voltage for the ZGNR structure shown in Fig. 8.

Figure 10

Left: The C${}_{60}$-6-terminal structure and real-space $\text{DOS}$ at the Fermi level. Right: The averaged effective potential projected onto a plane cutting through the ${\text{C}}_{60}$ molecule.

Figure 11

The transmission function evaluated at the Fermi energy between the six leads. The element $(i,j)$ refers to the transmission from lead $i$ to lead $j$; the diagonal element is the reflection coefficient.