API Reference

This section provides detailed API documentation for all gauNEGF modules.

Core Modules

NEGF Base Class

Self-consistent field (SCF) implementation for Non-Equilibrium Green’s Function calculations.

This module provides the base NEGF class for performing self-consistent DFT+NEGF calculations using Gaussian quantum chemistry package. It implements the energy-independent self-energy approach developed by Damle et al., which provides an efficient approximation for molecular transport calculations. The module handles:

• Integration with Gaussian for DFT calculations

• SCF convergence with Pulay mixing [2]

• Contact self-energy calculations

• Voltage bias and electric field effects

• Spin-polarized calculations

The implementation follows the standard NEGF formalism where the density matrix is calculated self-consistently with the Fock matrix from Gaussian DFT calculations. Convergence is accelerated using the direct inversion in the iterative subspace (DIIS) method developed by Pulay [2]. The core NEGF-DFT implementation is based on the ANT.Gaussian approach developed by Palacios et al. [3], which pioneered the integration of NEGF with Gaussian-based DFT calculations.

References

[1] Damle, P., Ghosh, A. W., & Datta, S. (2002). First-principles analysis of molecular

conduction using quantum chemistry software. Chemical Physics, 281(2-3), 171-187. DOI: 10.1016/S0301-0104(02)00496-2

[2] Pulay, P. (1980). Convergence acceleration of iterative sequences. The case of SCF

iteration. Chemical Physics Letters, 73(2), 393-398. DOI: 10.1016/0009-2614(80)80396-4

[3] Palacios, J. J., Pérez-Jiménez, A. J., Louis, E., & Vergés, J. A. (2002).

Fullerene-based molecular nanobridges: A first-principles study. Physical Review B, 66(3), 035322. DOI: 10.1103/PhysRevB.66.035322

class gauNEGF.scf.NEGF(fn, basis='chkbasis', func='hf', spin='r', fullSCF=True, route=None, section=None, nPulay=4)

Bases: object

Non-Equilibrium Green’s Function calculator integrated with Gaussian DFT.

This class implements the energy-independent NEGF approach developed by Damle et al. [1] for efficient molecular transport calculations. It manages the self-consistent field calculation between NEGF transport calculations and DFT electronic structure calculations using Gaussian.

The energy-independent approximation assumes constant self-energies, which significantly reduces computational cost while maintaining accuracy for many molecular systems.

The class handles: - Interaction with Gaussian - Management of density and Fock matrices - Pulay mixing for convergence [2] - Constant (energy-independent) contact self-energies - Voltage bias effects

For energy-dependent calculations, see the NEGFE subclass.

Parameters:

• fn (str) – Base filename for Gaussian input/output files (without extension)

• basis (str, optional) – Gaussian basis set name (default: ‘chkbasis’)

• func (str, optional) – DFT functional to use (default: ‘b3pw91’)

• spin ({'r', 'u', 'ro', 'g'}, optional) – Spin configuration: - ‘r’: restricted - ‘u’: unrestricted - ‘ro’: restricted open-shell - ‘g’: generalized open-shell (non-collinear) (default: ‘r’)

• fullSCF (bool, optional) – Whether to run full SCF or use Harris approximation (default: True)

• route (str, optional) – Additional Gaussian route commands (default: ‘’)

• nPulay (int, optional) – Number of previous iterations to use in Pulay mixing (default: 4)

F

Fock matrix in eV

Type:

ndarray

P

Density matrix

Type:

ndarray

S

Overlap matrix

Type:

ndarray

fermi

Fermi energy in eV

Type:

float

nelec

Number of electrons

Type:

float

References

[1] Damle, P., Ghosh, A. W., & Datta, S. (2002). First-principles analysis of molecular

conduction using quantum chemistry software. Chemical Physics, 281(2-3), 171-187. DOI: 10.1016/S0301-0104(02)00496-2

[2] Pulay, P. (1980). DOI: 10.1016/0009-2614(80)80396-4

FockToP()

Calculate density matrix from Fock matrix using energy-independent approach.

This method implements the energy-independent density matrix calculation from Damle et al. (2002). By assuming constant self-energies, the density matrix can be calculated analytically without energy integration, significantly reducing computational cost.

The method: 1. Transforms Fock and Gamma matrices to orthogonal basis 2. Diagonalizes the transformed Fock matrix 3. Updates Fermi energy if needed 4. Calculates density matrix analytically 5. Transforms back to non-orthogonal basis

Returns:

(eigenvalues, occupations) sorted by energy

Return type:

tuple

References

[1] Damle, P., Ghosh, A. W., & Datta, S. (2002). First-principles analysis of molecular

conduction using quantum chemistry software. Chemical Physics, 281(2-3), 171-187. DOI: 10.1016/S0301-0104(02)00496-2

PMix(damping, Pulay=False)

Mix old and new density matrices using damping or Pulay DIIS method [2].

The Pulay mixing method (also known as DIIS - Direct Inversion in the Iterative Subspace) uses information from previous iterations to predict the optimal density matrix. This method is particularly effective for systems with challenging convergence behavior, and closely follows ANT.Gaussian approaches.

Parameters:

• damping (float) – Mixing parameter between 0 and 1

• Pulay (bool, optional) – Whether to use Pulay mixing (default: False)

Returns:

(RMSDP, MaxDP) - RMS and maximum density matrix differences

Return type:

tuple

Notes

The Pulay DIIS method [2] minimizes the error in the iterative subspace spanned by previous density matrices. This often provides faster and more stable convergence compared to simple damping, especially for systems with strong electron correlation or near degeneracies [3].

References

[2]

Pulay, P. (1980). DOI: 10.1016/0009-2614(80)80396-4

[3]

Palacios, J. J., et al. (2002). DOI: 10.1103/PhysRevB.66.035322

PToFock()

Calculate new Fock matrix from current density matrix using Gaussian.

Returns:

Energy difference from previous iteration

Return type:

float

SCF(conv=1e-05, damping=0.02, maxcycles=100, checkpoint=True, pulay=True)

Run self-consistent field calculation until convergence.

The SCF cycle alternates between Fock matrix construction and density matrix updates until convergence is reached. Convergence acceleration is achieved through either simple damping or Pulay DIIS mixing [2]. The Pulay method is applied every nPulay iterations (where nPulay is set in __init__) [3].

References

[3]

Palacios, J. J., et al. (2002). DOI: 10.1103/PhysRevB.66.035322

Parameters:

• conv (float, optional) – Convergence criterion for energy and density (default: 1e-5)

• damping (float, optional) – Mixing parameter between 0 and 1 (default: 0.02)

• maxcycles (int, optional) – Maximum number of SCF cycles (default: 100)

• checkpoint (bool, optional) – Save density matrix at each iteration and load if job interrupted (default: True)

• pulay (bool, optional) – Whether to use Pulay DIIS mixing [2] (default: True)

Returns:

(count, PP, TotalE) - cycle number, number of electrons, and DFT energy at each cycle

Return type:

tuple

Notes

Convergence is determined by three criteria: 1. Energy change (dE) 2. RMS density matrix difference (RMSDP) 3. Maximum density matrix difference (MaxDP)

All three must be below the convergence threshold.

The Pulay DIIS method [2] is applied every nPulay iterations when enabled, which often provides faster and more stable convergence compared to simple damping, especially for challenging systems.

References

[2]

Pulay, P. (1980). DOI: 10.1016/0009-2614(80)80396-4

__init__(fn, basis='chkbasis', func='hf', spin='r', fullSCF=True, route=None, section=None, nPulay=4)

Initialize NEGF calculator and run initial DFT calculation.

Sets up the calculator with specified parameters and runs an initial DFT calculation using Gaussian to obtain the starting Fock and overlap matrices.

Parameters are documented in class docstring.

getHOMOLUMO()

Calculate HOMO and LUMO energies.

Diagonalizes the Fock matrix in orthogonalized basis to get orbital energies, then identifies HOMO and LUMO based on electron occupation and spin configuration.

Returns:

Array of [HOMO, LUMO] energies in eV

Return type:

ndarray

getSigma(E=0)

runDFT(fullSCF=True)

Run DFT calculation using Gaussian.

Performs either a full SCF calculation or generates an initial Harris guess using Gaussian. Updates the Fock matrix and orbital indices after completion.

Parameters:

fullSCF (bool, optional) – If True, runs full SCF to convergence. If False, uses Harris guess only. (default: True)

Notes

• Attempts to load from checkpoint file first

• Falls back to full calculation if checkpoint fails

• Updates self.F and self.locs with new Fock matrix

saveMAT(matfile='out.mat')

Save calculation results to MATLAB format file.

Parameters:

matfile (str, optional) – Output filename (default: “out.mat”)

Returns:

Fock matrix in orthogonalized basis

Return type:

ndarray

setContacts(lContact=None, rContact=None)

Set contact atoms and get corresponding orbital indices.

Identifies the orbital indices corresponding to the specified contact atoms. If no contacts specified, uses all orbitals.

Parameters:

• lContact (int or array-like, optional) – Atom number(s) for left contact. None means all atoms. (default: None)

• rContact (int or array-like, optional) – Atom number(s) for right contact. None means all atoms. (default: None)

Returns:

(left_orbital_indices, right_orbital_indices)

Return type:

tuple of ndarrays

setDen(P_)

Set the density matrix and update dependent quantities.

Updates the density matrix, stores it in Gaussian format, recalculates the number of electrons, and updates the Fock matrix.

Parameters:

P (ndarray) – New density matrix

setFock(F_)

Set the Fock matrix, converting from eV to atomic units.

Parameters:

F (ndarray) – Fock matrix in Hartree units

setSigma(lContact=None, rContact=None, sig=-0 - 0.1j, sig2=None)

Set self-energies for left and right contacts.

Configures the contact self-energies, handling various input formats and spin configurations. Self-energies can be scalar, vector, or matrix, with automatic handling of spin degrees of freedom.

Parameters:

• lContact (array-like) – Atom numbers for left contact, all atoms if None (default: None)

• rContact (array-like) – Atom numbers for right contact, all atoms if None (default: None)

• sig (scalar or array-like, optional) – Self-energy for left contact. Can be: - scalar: same value for all orbitals - vector: one value per orbital - matrix: full self-energy matrix (default: -0.1j)

• sig2 (scalar, array-like, or None, optional) – Self-energy for right contact. If None, uses sig. Same format options as sig. (default: None)

Returns:

Left and right contact orbital indices

Return type:

tuple

Notes

• Handles spin configurations (‘r’, ‘u’, ‘ro’, ‘g’)

• Automatically expands scalar/vector inputs to full matrices

• Verifies matrix dimensions match Fock matrix

• Updates self.sigma1, self.sigma2, self.sigma12

Raises:

Exception – If matrix dimensions don’t match or invalid input format

setVoltage(qV, fermi=nan, Emin=None, Eminf=None)

Set voltage bias and Fermi energy, updating electric field.

Applies a voltage bias between contacts and updates the chemical potentials and electric field. Can optionally set the Fermi energy and integration limits.

Parameters:

• qV (float) – Voltage bias in eV

• fermi (float, optional) – Fermi energy in eV. If not provided, will be calculated or use existing value (default: np.nan)

• Emin (float, optional) – Minimum energy for integration in eV (default: None)

• Eminf (float, optional) – Lower bound energy in eV (default: None)

Notes

• Requires contacts to be set first

• Updates chemical potentials as fermi ± qV/2

• Calculates and applies electric field between contacts

• If fermi not provided, uses (HOMO+LUMO)/2 for initial guess

updateN()

Update the total number of electrons from the density matrix.

Calculates the number of electrons by tracing the product of the density and overlap matrices. For restricted calculations, multiplies by 2 to account for spin degeneracy.

Returns:

Total number of electrons

Return type:

float

writeChk()

Write current state to Gaussian checkpoint file.

Energy-Dependent NEGF

Energy-dependent extensions to the NEGF class for quantum transport calculations.

This module extends the base NEGF class to handle energy-dependent self-energies, temperature effects, and advanced Fermi energy search methods. It provides support for Bethe lattice and 1D chain contacts [1] with proper energy integration.

References

[1]

Jacob, D. & Palacios, J. J. Chem. Phys. 134, 044118 (2011)

class gauNEGF.scfE.NEGFE(fn, basis='chkbasis', func='hf', spin='r', fullSCF=True, route=None, section=None, nPulay=4)

Bases: NEGF

Extended NEGF class with energy-dependent self-energies and temperature effects.

This class extends the base NEGF implementation to handle: - Energy-dependent contact self-energies - Finite temperature effects - Advanced Fermi energy search methods - Integration methods for the density matrix

Inherits all attributes and methods from NEGF class.

FockToP()

Calculate density matrix using energy-dependent integration.

Performs complex contour integration for the equilibrium part and real axis integration for the non-equilibrium part. Updates Fermi energy using specified method if required.

Returns:

(energies, occupations) - Sorted eigenvalues and occupations

Return type:

tuple

PToFock()

Calculate new Fock matrix and update surfG object.

Returns:

Energy difference from previous iteration

Return type:

float

getSigma(E)

Get contact self-energies at specified energy.

Parameters:

E (float) – Energy in eV

Returns:

(left_sigma, right_sigma) - Contact self-energies

Return type:

tuple

integralCheck(cycles=10, damp=0.02, pauseFermi=False)

Check and optimize integration parameters.

Parameters:

• tol (float, optional) – Tolerance for integration (default: 1e-4)

• cycles (int, optional) – Number of SCF cycles to run (default: 10)

• damp (float, optional) – Damping parameter (default: 0.02)

• pauseFermi (bool, optional) – Whether to pause Fermi updates (default: False)

setContact1D(contactList, tauList=None, stauList=None, alphas=None, aOverlaps=None, betas=None, bOverlaps=None, neList=None, eta=1e-09, T=0)

Set energy-dependent 1D chain contacts.

Parameters:

• contactList (list) – List of atom indices for contacts

• tauList (list or None, optional) – Coupling matrices or atom indices (default: None)

• stauList (list or None, optional) – Overlap matrices (default: None)

• alphas (array_like or None, optional) – On-site energies (default: None)

• aOverlaps (array_like or None, optional) – On-site overlaps (default: None)

• betas (array_like or None, optional) – Hopping energies (default: None)

• bOverlaps (array_like or None, optional) – Hopping overlaps (default: None)

• neList (list or None, optional) – Number of electrons per unit cell (default: None)

• eta (float, optional) – Broadening parameter in eV (default: 1e-9)

• T (float, optional) – Temperature in Kelvin (default: 0)

Returns:

Left and right contact orbital indices

Return type:

tuple

setContactBethe(contactList, latFile='Au', eta=1e-09, T=0)

Set energy-dependent Bethe lattice contacts.

Parameters:

• contactList (list) – List of atom indices for contacts

• latFile (str, optional) – Lattice parameter file (default: ‘Au’)

• eta (float, optional) – Broadening parameter in eV (default: 1e-9)

• T (float, optional) – Temperature in Kelvin (default: 0)

Returns:

Left and right contact orbital indices

Return type:

tuple

setIntegralLimits(N1=None, N2=None, Nnegf=None, tol=0.0001, Emin=None)

Set integration parameters for density calculation.

Parameters:

• N1 (int) – Number of points for complex contour

• N2 (int) – Number of points for real axis

• Emin (float or None, optional) – Minimum energy for integration (default: None)

setSigma(lContact=None, rContact=None, sig=-0 - 0.1j, sig2=None, T=0)

Set constant self-energy contacts with temperature.

Parameters:

• lContact (list) – Atom numbers for left contact, all atoms if None (default: None)

• rContact (array-like) – Atom numbers for right contact, all atoms if None (default: None)

• sig (complex, optional) – Left contact self-energy (default: -0.1j)

• sig2 (complex or None, optional) – Right contact self-energy (default: None)

• T (float, optional) – Temperature in Kelvin (default: 0)

Returns:

Left and right contact orbital indices

Return type:

tuple

setVoltage(qV, fermi=nan, Emin=None, Eminf=None, fermiMethod='muller')

Set voltage bias and Fermi search method.

Parameters:

• qV (float) – Applied voltage in eV

• fermi (float, optional) – Fermi energy in eV (default: np.nan)

• Emin (float, optional) – Minimum energy for integration (default: None)

• Eminf (float, optional) – Minimum energy for Fermi search (default: None)

• fermiMethod (str, optional) – Method for Fermi search: ‘muller’, ‘secant’, or ‘predict’ (default: ‘muller’)

Density Module

Density matrix calculation methods for quantum transport simulations.

This module provides functions for calculating density matrices in quantum transport calculations using various integration methods:

• Analytical integration for energy-independent self-energies

• Complex contour integration for equilibrium calculations

• Real-axis integration for non-equilibrium calculations

• Parallel processing support for large systems

Notes

The module supports both serial and parallel computation modes, automatically selecting the most efficient method based on system size and available resources. Temperature effects are included through Fermi-Dirac statistics.

gauNEGF.density.DOSg(F, S, g, E)

Calculate density of states at given energy.

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• E (float) – Energy in eV

Returns:

Density of states at energy E

Return type:

float

gauNEGF.density.Gr(F, S, g, E)

Calculate retarded Green’s function at given energy.

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• E (float) – Energy in eV

Returns:

Retarded Green’s function G(E) = [ES - F - Σ(E)]^(-1)

Return type:

ndarrayh

gauNEGF.density.bisectFermi(V, Vc, D, Gam, Nexp, conv=0.001, Eminf=-1000000.0)

Find Fermi energy using bisection method.

Uses bisection to find the Fermi energy that gives the expected number of electrons. The search is performed using the analytical density matrix calculation.

Parameters:

• V (ndarray) – Eigenvectors of Fock matrix in orthogonalized basis

• Vc (ndarray) – Inverse conjugate transpose of V

• D (ndarray) – Eigenvalues of Fock matrix

• Gam (ndarray) – Broadening matrix Γ = i[Σ - Σ†] in orthogonalized basis

• Nexp (float) – Expected number of electrons

• conv (float, optional) – Convergence criterion for electron number (default: 1e-3)

• Eminf (float, optional) – Lower bound for integration in eV (default: -1e6)

Returns:

Fermi energy in eV that gives the expected electron count

Return type:

float

Notes

The search is bounded by the minimum and maximum eigenvalues. Warns if maximum iterations reached without convergence.

gauNEGF.density.calcEmin(F, S, g, tol=1e-05, maxN=1000)

gauNEGF.density.calcFermi(g, ne, Emin, Emax, fermiGuess=0, N1=100, N2=50, Eminf=-1000000.0, T=0, tol=0.0001, maxcycles=20, nOrbs=0)

Calculate Fermi energy using bisection method.

Parameters:

• g (surfG object) – Surface Green’s function calculator

• ne (float) – Target number of electrons

• Emin (float) – Lower bound for complex contour in eV

• Emax (float) – Upper bound for search in eV

• fermiGuess (float, optional) – Initial guess for Fermi energy in eV (default: 0)

• N1 (int, optional) – Number of complex contour points (default: 100)

• N2 (int, optional) – Number of real axis points (default: 50)

• Eminf (float, optional) – Lower bound for real axis integration in eV (default: -1e6)

• tol (float, optional) – Convergence tolerance (default: 1e-4)

• maxcycles (int, optional) – Maximum number of iterations (default: 20)

• nOrbs (int, optional) – Number of orbitals to consider, 0 for all (default: 0)

Returns:

(fermi, Emin, N1, N2) - Optimized parameters: - fermi: Calculated Fermi energy in eV - Emin: Lower bound for complex contour - N1: Number of complex contour points - N2: Number of real axis points

Return type:

tuple

gauNEGF.density.calcFermiBisect(g, ne, Emin, Ef, N, tol=0.0001, conv=0.001, maxcycles=10, T=0)

Calculate Fermi energy of system using bisection

gauNEGF.density.calcFermiMuller(g, ne, Emin, Ef, N, tol=0.0001, conv=0.001, maxcycles=10, T=0)

Calculate Fermi energy using Muller’s method, starting with 3 initial points

gauNEGF.density.calcFermiSecant(g, ne, Emin, Ef, N, tol=0.0001, conv=0.001, maxcycles=10, T=0)

Calculate Fermi energy using Secant method, updating dE at each step

gauNEGF.density.density(V, Vc, D, Gam, Emin, mu)

Calculate density matrix using analytical integration for energy-independent self-energies.

Implements the analytical integration method from Eq. 27 in PRB 65, 165401 (2002). The method assumes energy-independent self-energies and uses the spectral function representation of the density matrix.

Parameters:

• V (ndarray) – Eigenvectors of Fock matrix in orthogonalized basis

• Vc (ndarray) – Inverse conjugate transpose of V

• D (ndarray) – Eigenvalues of Fock matrix

• Gam (ndarray) – Broadening matrix Γ = i[Σ - Σ†] in orthogonalized basis

• Emin (float) – Lower bound for integration in eV

• mu (float) – Chemical potential in eV

Returns:

Density matrix in orthogonalized basis

Return type:

ndarray

Notes

The integration is performed analytically using the residue theorem. The result includes contributions from poles below the Fermi energy.

gauNEGF.density.densityComplex(F, S, g, Emin, mu, tol=0.001, T=300, debug=False)

Calculate equilibrium density matrix using complex contour integration.

Performs numerical integration along a complex contour that encloses the poles of the Fermi function. More efficient than real-axis integration for equilibrium calculations.

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• Emin (float) – Lower bound for integration in eV

• mu (float) – Chemical potential in eV

• N (int, optional) – Number of integration points (default: 100)

• T (float, optional) – Temperature in Kelvin (default: 300)

• parallel (bool, optional) – Whether to use parallel processing (default: False)

• numWorkers (int, optional) – Number of worker processes for parallel mode

• showText (bool, optional) – Whether to print progress messages (default: True)

• method ({'ant', 'legendre', 'chebyshev'}, optional) – Integration method to use (default: ‘ant’)

• use_adaptive (bool, optional) – If True, use dyadic adaptive integration on the contour integral. Requires len(w) = 2^n + 1. Defaults to False.

• adaptive_tol (float, optional) – Convergence tolerance for adaptive integration (default: 1e-3)

Returns:

Equilibrium density matrix

Return type:

ndarray

Notes

The ‘ant’ method uses a modified Gauss-Chebyshev quadrature optimized for transport calculations, matching the ANT.Gaussian implementation.

gauNEGF.density.densityComplexN(F, S, g, Emin, mu, N=100, T=300, parallel=False, numWorkers=None, showText=True, method='ant')

Calculate equilibrium density matrix using complex contour integration.

Performs numerical integration along a complex contour that encloses the poles of the Fermi function. More efficient than real-axis integration for equilibrium calculations.

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• Emin (float) – Lower bound for integration in eV

• mu (float) – Chemical potential in eV

• N (int, optional) – Number of integration points (default: 100)

• T (float, optional) – Temperature in Kelvin (default: 300)

• parallel (bool, optional) – Whether to use parallel processing (default: False)

• numWorkers (int, optional) – Number of worker processes for parallel mode

• showText (bool, optional) – Whether to print progress messages (default: True)

• method ({'ant', 'legendre', 'chebyshev'}, optional) – Integration method to use (default: ‘ant’)

• use_adaptive (bool, optional) – If True, use dyadic adaptive integration on the contour integral. Requires len(w) = 2^n + 1. Defaults to False.

• adaptive_tol (float, optional) – Convergence tolerance for adaptive integration (default: 1e-3)

Returns:

Equilibrium density matrix

Return type:

ndarray

Notes

The ‘ant’ method uses a modified Gauss-Chebyshev quadrature optimized for transport calculations, matching the ANT.Gaussian implementation.

gauNEGF.density.densityGrid(F, S, g, mu1, mu2, ind=None, tol=0.001, T=300, debug=False)

Calculate non-equilibrium density matrix using real-axis integration.

Performs numerical integration for the non-equilibrium part of the density matrix when a bias voltage is applied. Uses ANT modified Gauss-Chebyshev quadrature.

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• mu1 (float) – Left contact chemical potential in eV

• mu2 (float) – Right contact chemical potential in eV

• ind (int or None, optional) – Contact index (None for total) (default: None)

• N (int, optional) – Number of integration points (default: 100)

• T (float, optional) – Temperature in Kelvin (default: 300)

• parallel (bool, optional) – Whether to use parallel processing (default: False)

• numWorkers (int, optional) – Number of worker processes for parallel mode

• showText (bool, optional) – Whether to print progress messages (default: True)

Returns:

Non-equilibrium contribution to density matrix

Return type:

ndarray

gauNEGF.density.densityGridN(F, S, g, mu1, mu2, ind=None, N=100, T=300, parallel=False, numWorkers=None, showText=True)

Calculate non-equilibrium density matrix using real-axis integration.

Performs numerical integration for the non-equilibrium part of the density matrix when a bias voltage is applied. Uses Gauss-Legendre quadrature.

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• mu1 (float) – Left contact chemical potential in eV

• mu2 (float) – Right contact chemical potential in eV

• ind (int or None, optional) – Contact index (None for total) (default: None)

• N (int, optional) – Number of integration points (default: 100)

• T (float, optional) – Temperature in Kelvin (default: 300)

• parallel (bool, optional) – Whether to use parallel processing (default: False)

• numWorkers (int, optional) – Number of worker processes for parallel mode

• showText (bool, optional) – Whether to print progress messages (default: True)

Returns:

Non-equilibrium contribution to density matrix

Return type:

ndarray

gauNEGF.density.densityGridTrap(F, S, g, mu1, mu2, ind=None, N=100, T=300)

gauNEGF.density.densityReal(F, S, g, Emin, mu, tol=0.001, T=0, maxN=1000, debug=False)

Calculate equilibrium density matrix using adaptive real-axis integration.

Wrapper for densityRealN() using the tol and maxN specification to determine grid size

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• Emin (float) – Lower bound for integration in eV

• mu (float) – Chemical potential in eV

• tol (float, optional) – Convergence tolerance (default: 1e-3)

• T (float, optional) – Temperature in Kelvin (default: 300)

• maxN (int, optional) – Maximum number of integration points (default: 1000)

• debug (bool, optional) – Whether to print per-iteration diagnostics (default: False)

Returns:

Density matrix

Return type:

ndarray

gauNEGF.density.densityRealN(F, S, g, Emin, mu, N=100, T=300, parallel=False, numWorkers=None, showText=True)

Calculate equilibrium density matrix using real-axis integration on a specified grid.

Performs numerical integration along the real energy axis using Gauss-Legendre quadrature. Suitable for equilibrium calculations with energy-dependent self-energies.

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• Emin (float) – Lower bound for integration in eV

• mu (float) – Chemical potential in eV

• N (int, optional) – Number of integration points (default: 100)

• T (float, optional) – Temperature in Kelvin (default: 300)

• parallel (bool, optional) – Whether to use parallel processing (default: False)

• numWorkers (int, optional) – Number of worker processes for parallel mode

• showText (bool, optional) – Whether to print progress messages (default: True)

Returns:

Density matrix

Return type:

ndarray

gauNEGF.density.fermi(E, mu, T)

Calculate Fermi-Dirac distribution.

Parameters:

• E (float) – Energy in eV

• mu (float) – Chemical potential in eV

• T (float) – Temperature in Kelvin

Returns:

Fermi-Dirac occupation number

Return type:

float

gauNEGF.density.getANTPoints(N)

Generate integration points and weights matching ANT.Gaussian implementation.

Follows the IntCompPlane subroutine in device.F90 from ANT.Gaussian package. Uses a modified Gauss-Chebyshev quadrature scheme optimized for transport calculations. _Note: Always generates an even number of points._

Parameters:

N (int) – Number of integration points

Returns:

(points, weights) - Arrays of integration points and weights

Return type:

tuple

gauNEGF.density.getFermi1DContact(gSys, ne, ind=0, tol=0.0001, Eminf=-1000000.0, T=0, maxcycles=1000)

Calculate Fermi energy for a 1D chain contact.

Specialized version of getFermiContact for 1D chain contacts, handling the periodic boundary conditions correctly.

Parameters:

• gSys (surfG object) – Surface Green’s function calculator for the full system

• ne (float) – Target number of electrons per unit cell

• ind (int, optional) – Contact index (0 for left, -1 for right) (default: 0)

• tol (float, optional) – Convergence tolerance (default: 1e-4)

• Eminf (float, optional) – Lower bound for integration (default: -1e6)

• maxcycles (int, optional) – Maximum number of iterations (default: 1000)

Returns:

Fermi energy in eV

Return type:

float

gauNEGF.density.getFermiContact(g, ne, tol=0.0001, Eminf=-1000000.0, maxcycles=1000, T=0, nOrbs=0)

Calculate Fermi energy for a contact.

Determines the Fermi energy for a contact system (Bethe lattice or 1D chain) by matching the electron count.

Parameters:

• g (surfG object) – Surface Green’s function calculator

• ne (float) – Target number of electrons

• tol (float, optional) – Convergence tolerance (default: 1e-4)

• Eminf (float, optional) – Lower bound for integration (default: -1e6)

• maxcycles (int, optional) – Maximum number of iterations (default: 1000)

• nOrbs (int, optional) – Number of orbitals to consider (0 for all) (default: 0)

Returns:

Fermi energy in eV

Return type:

float

gauNEGF.density.integralFit(F, S, g, mu, Eminf=-1000000.0, tol=1e-05, T=0, maxN=1000)

Optimize integration parameters for density calculations.

Determines optimal integration parameters by iteratively testing convergence of the density matrix calculation.

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• mu (float) – Equilibrium contact fermi energy in eV

• Eminf (float) – Lower bound for integration in eV (default: -1e6)

• tol (float, optional) – Convergence tolerance (default: 1e-5)

• T (float) – Temperature in Kelvin for Fermi broadening (default: 0)

• maxN (int, optional) – Max grid points and Emin search iterations(default: 1000)

Returns:

(Emin, N1, N2) - Optimized integration parameters: - Emin: Lower bound for complex contour - N1: Number of complex contour points - N2: Number of real axis points

Return type:

tuple

Notes

The optimization process: 1. Finds Emin by checking DOS convergence 2. Optimizes N1 for complex contour integration 3. Optimizes N2 for real axis integration

gauNEGF.density.integralFitNEGF(F, S, g, fermi, qV, Eminf=-1000000.0, tol=1e-05, T=0, maxGrid=1000)

Determines number of for non-equilibrium density calculations.

Same procedure as integralFit() but applied to densityGrid()

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• mu1 (float) – Left contact Fermi energy in eV

• mu2 (float) – Right contact Fermi energy in eV

• Eminf (float) – Lower bound for integration in eV (default: -1e6)

• tol (float, optional) – Convergence tolerance (default: 1e-5)

• T (float) – Temperature in Kelvin for Fermi broadening (default: 0)

• maxGrid (int, optional) – Maximum number of gridpoints (default: 1000)

Returns:

Number of grid points

Return type:

int

gauNEGF.density.integratePoints(computePointFunc, numPoints, parallel=False, numWorkers=None, chunkSize=None, debug=False)

Perform parallel or serial integration for quantum transport calculations.

This function provides a flexible integration framework that automatically chooses between numpy’s built-in parallelization for matrix operations and process-level parallelization for large workloads.

Parameters:

• computePointFunc (callable) – Function that computes a single integration point. Should take an integer index i and return a matrix/array.

• numPoints (int) – Total number of points to integrate

• parallel (bool, optional) – Whether to force process-level parallel processing (default: False)

• numWorkers (int, optional) – Number of worker processes. If None, automatically determined based on system resources and workload.

• chunkSize (int, optional) – Size of chunks for parallel processing. If None, automatically optimized based on numPoints and numWorkers.

• debug (bool, optional) – Whether to print debug information (default: False)

Returns:

Sum of all computed points

Return type:

ndarray

Notes

• Automatically detects SLURM environment for HPC compatibility

• Falls back to serial processing if parallel execution fails

• Uses numpy’s built-in parallelization for small workloads

• Switches to process-level parallelization for:

  • Large number of points (≥100)

  • Many available cores (≥32)

  • When explicitly requested via parallel=True

gauNEGF.density.integratePointsAdaptiveANT(computePoint, tol=0.001, maxN=1458, debug=False)

Adaptive integration using ANT-modified Gauss-Chebyshev quadrature (IntCompPlane subroutine from ANT.Gaussian package)

Parameters:

• computePoint (callable) – Function that computes a single integration point. Should take a weight and point and return a matrix/array.

• tol (float, optional) – Tolerance for the adaptive integration.

• maxN (int, optional) – Maximum number of points to integrate.

• debug (bool, optional) – Whether to print debug information.

Returns:

Integral of the function.

Return type:

ndarray

Transport Module

Transport calculations for quantum systems using Non-Equilibrium Green’s Functions.

This module provides functions for calculating quantum transport properties: - Coherent transmission through molecular junctions - Spin-dependent transport calculations - Current calculations at finite bias - Density of states calculations

The module supports both energy-independent and energy-dependent self-energies, with implementations for both spin-restricted and spin-unrestricted calculations. Spin-dependent transport follows the formalism described in [1].

References

[1]

Herrmann, C., Solomon, G. C., & Ratner, M. A. J. Chem. Theory Comput. 6, 3078 (2010) DOI: 10.1021/acs.jctc.9b01078

gauNEGF.transport.DOS(Elist, F, S, sig1, sig2)

Calculate density of states with energy-independent self-energies.

Parameters:

• Elist (array_like) – List of energies in eV to calculate DOS at

• F (ndarray) – Fock matrix, NxN

• S (ndarray) – Overlap matrix, NxN

• sig1 (ndarray) – Left contact self-energy, vector (1xN) or matrix (NxN)

• sig2 (ndarray) – Right contact self-energy, vector (1xN) or matrix (NxN)

Returns:

(DOS, DOSList) where: - DOS: Total density of states at each energy - DOSList: Site-resolved DOS at each energy

Return type:

tuple

gauNEGF.transport.DOSE(Elist, F, S, g)

Calculate density of states with energy-dependent self-energies.

Parameters:

• Elist (array_like) – List of energies in eV to calculate DOS at

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

Returns:

(DOS, DOSList) where: - DOS: Total density of states at each energy - DOSList: Site-resolved DOS at each energy

Return type:

tuple

gauNEGF.transport.cohTrans(Elist, F, S, sig1, sig2)

Calculate coherent transmission with energy-independent self-energies.

Parameters:

• Elist (array_like) – List of energies in eV to calculate transmission at

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• sig1 (ndarray) – Left contact self-energy (vector or matrix)

• sig2 (ndarray) – Right contact self-energy (vector or matrix)

Returns:

Transmission values at each energy

Return type:

list

Notes

Supports both vector and matrix self-energies. For vector self-energies, diagonal matrices are constructed internally.

gauNEGF.transport.cohTransE(Elist, F, S, g)

Calculate coherent transmission with energy-dependent self-energies.

Parameters:

• Elist (array_like) – List of energies in eV to calculate transmission at

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

Returns:

Transmission values at each energy

Return type:

list

Notes

Uses the surface Green’s function calculator to compute energy-dependent self-energies at each energy point.

gauNEGF.transport.cohTransSpin(Elist, F, S, sig1, sig2, spin='u')

Calculate spin-dependent coherent transmission with energy-independent self-energies.

Parameters:

• Elist (array_like) – List of energies in eV to calculate transmission at

• F (ndarray) – Fock matrix (2N x 2N for spin-unrestricted)

• S (ndarray) – Overlap matrix (2N x 2N for spin-unrestricted)

• sig1 (ndarray) –

  Left contact self-energy (spin independent vector (1xN) or matrix (NxN),

  spin dependent matrix (2Nx2N))

• sig2 (ndarray) –

  Right contact self-energy (spin independent vector (1xN) or matrix (NxN),

  spin dependent matrix (2Nx2N))

• spin (str, optional) – Spin basis {‘r’, ‘u’, ‘ro’, or ‘g’} (default: ‘u’)

Returns:

(Tr, Tspin) where: - Tr: Total transmission at each energy - Tspin: Array of spin-resolved transmissions [T↑↑, T↑↓, T↓↑, T↓↓]

Return type:

tuple

Notes

For collinear spin calculations (‘u’ or ‘ro’), the matrices are arranged in blocks: [F↑↑ 0 ] [S↑↑ 0 ] [0 F↓↓], [0 S↓↓] For generalized spin basis (‘g’), each orbital contains a 2x2 spinor block: [F↑↑ F↑↓] [S↑↑ S↑↓] [F↓↑ F↓↓], [S↓↑ S↓↓] which are then combined into a 2Nx2N matrix.

gauNEGF.transport.cohTransSpinE(Elist, F, S, g, spin='u')

Calculate spin-dependent coherent transmission with energy-dependent self-energies.

Parameters:

• Elist (array_like) – List of energies in eV to calculate transmission at

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• spin (str, optional) – Spin basis {‘r’, ‘u’, ‘ro’, or ‘g’} (default: ‘u’)

Returns:

(Tr, Tspin) where: - Tr: Total transmission at each energy - Tspin: Array of spin-resolved transmissions [T↑↑, T↑↓, T↓↑, T↓↓]

Return type:

tuple

gauNEGF.transport.current(F, S, sig1, sig2, fermi, qV, T=0, spin='r', dE=0.01)

Calculate coherent current using NEGF with energy-independent self-energies.

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• sig1 (ndarray) – Left contact self-energy (vector or matrix)

• sig2 (ndarray) – Right contact self-energy (vector or matrix)

• fermi (float) – Fermi energy in eV

• qV (float) – Applied bias voltage in eV

• T (float) – Temperature in Kelvin (default: 0)

• spin (str, optional) – Spin configuration (‘r’ for restricted) (default: ‘r’)

• dE (float, optional) – Energy step for integration in eV (default: 0.01)

Returns:

Current in Amperes

Return type:

float

gauNEGF.transport.currentE(F, S, g, fermi, qV, T=0, spin='r', dE=0.01)

Calculate coherent current at T=0K using NEGF with energy-dependent self-energies.

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• g (surfG object) – Surface Green’s function calculator

• fermi (float) – Fermi energy in eV

• qV (float) – Applied bias voltage in eV

• T (float) – Temperature in Kelvin (default: 0)

• spin (str, optional) – Spin configuration (‘r’ for restricted) (default: ‘r’)

• dE (float, optional) – Energy step for integration in eV (default: 0.01)

Returns:

Current in Amperes

Return type:

float

gauNEGF.transport.currentF(fn, dE=0.01, T=0)

Calculate current from saved SCF matrix file.

Parameters:

• fn (str) – Filename of .mat file containing SCF data

• dE (float, optional) – Energy step for integration in eV (default: 0.01)

Returns:

Current in Amperes

Return type:

float

Notes

The .mat file should contain: - F: Fock matrix - S: Overlap matrix - sig1, sig2: Contact self-energies - fermi: Fermi energy - qV: Applied voltage - spin: Spin configuration

gauNEGF.transport.currentSpin(F, S, sig1, sig2, fermi, qV, T=0, spin='r', dE=0.01)

Calculate coherent spin current using NEGF with energy-independent self-energies.

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• sig1 (ndarray) – Left contact self-energy (vector or matrix)

• sig2 (ndarray) – Right contact self-energy (vector or matrix)

• fermi (float) – Fermi energy in eV

• qV (float) – Applied bias voltage in eV

• T (float) – Temperature in Kelvin (default: 0)

• spin (str, optional) – Spin configuration (‘r’ for restricted) (default: ‘r’)

• dE (float, optional) – Energy step for integration in eV (default: 0.01)

Returns:

Spin-currents (in Amperes) [I↑↑, I↑↓, I↓↑, I↓↓]

Return type:

list

Contact Models

Bethe Lattice

Surface Green’s function implementation for Bethe lattice contacts.

This module provides a Bethe lattice implementation for modeling semi-infinite metallic contacts in quantum transport calculations. It supports: - FCC [111] surface geometry with proper orbital symmetries - Slater-Koster parameterization for s, p, and d orbitals - Temperature-dependent calculations - Spin-restricted and unrestricted calculations

The implementation follows the ANT.Gaussian approach [1], using a minimal basis set with s, p, and d orbitals for each contact atom. The Bethe lattice model provides an efficient way to describe bulk metallic electrodes while maintaining proper orbital symmetries and electronic structure. This approach allows for accurate modeling of metal-molecule interfaces without the computational cost of explicit periodic boundary conditions.

References

[1] Jacob, D., & Palacios, J. J. (2011). Critical comparison of electrode models

in density functional theory based quantum transport calculations. The Journal of Chemical Physics, 134(4), 044118. DOI: 10.1063/1.3526044

class gauNEGF.surfGBethe.surfGB(F, S, contacts, bar, latFile='Au', spin='r', eta=1e-09, T=0)

Bases: object

Surface Green’s function calculator for Bethe lattice contacts.

This class implements the Bethe lattice approximation for modeling semi-infinite metallic contacts. It handles: - Contact geometry detection and setup - Slater-Koster parameter management - Energy-dependent self-energy calculations - Temperature effects - Spin configurations

The Bethe lattice model represents the contacts as a semi-infinite tree-like structure with proper coordination number and orbital symmetries matching those of bulk FCC metals. This provides an efficient way to compute surface Green’s functions and self-energies for the electrodes, as demonstrated by Jacob & Palacios in their 2011 paper [1].

Parameters:

• F (ndarray) – Fock matrix

• S (ndarray) – Overlap matrix

• contacts (list of lists) – List of atom indices for each contact

• bar (object) – Gaussian interface object containing basis set information

• latFile (str, optional) – Filename for Bethe lattice parameters (default: ‘Au’)

• spin (str, optional) – Spin configuration (‘r’ for restricted) (default: ‘r’)

• eta (float, optional) – Broadening parameter in eV (default: 1e-9)

• T (float, optional) – Temperature in Kelvin (default: 0)

cVecs

Normal vectors for each contact surface

Type:

list

latVecs

Lattice vectors for each contact

Type:

list

indsLists

Orbital indices for each contact atom

Type:

list

dirLists

Direction vectors for nearest neighbors

Type:

list

nIndLists

Nearest neighbor indices for each atom

Type:

list

gList

surfGBAt objects for each contact

Type:

list

References

[1] Jacob, D., & Palacios, J. J. (2011). Critical comparison of electrode models

in density functional theory based quantum transport calculations. The Journal of Chemical Physics, 134(4), 044118. DOI: 10.1063/1.3526044

constructMat(Mdict, dirCosines)

Construct hopping/overlap matrix using Slater-Koster formalism.

Builds a 9x9 matrix for s, p, and d orbital interactions based on the Slater-Koster two-center approximation. The matrix is first constructed assuming a [0,0,1] bond direction, then rotated to the given direction using direction cosines.

Parameters:

• Mdict (dict) – Dictionary of Slater-Koster parameters (ssσ, spσ, ppσ, etc.)

• dirCosines (ndarray) – Array [l,m,n] of direction cosines for the bond

Returns:

9x9 matrix containing orbital interactions in the rotated frame

Return type:

ndarray

Notes

Matrix blocks: - [0,0]: s-s interaction - [0:4,0:4]: s-p block - [0:4,4:9]: s-d and p-d blocks - [4:9,4:9]: d-d block

genNeighbors(plane_normal, first_neighbor)

Generate 12 nearest neighbor unit vectors for an FCC [111] surface.

Creates a list of unit vectors representing the 12 nearest neighbors in an FCC lattice: - 6 in-plane vectors forming a hexagonal pattern (3 pairs of opposite vectors) - 6 out-of-plane vectors forming triangular patterns (3 pairs of opposite vectors)

Parameters:

• plane_normal (ndarray) – Vector normal to the crystal plane (will be normalized)

• first_neighbor (ndarray) – Vector to one nearest neighbor (will be projected onto plane)

Returns:

12 unit vectors representing nearest neighbor directions

Return type:

list

getSigma(Elist=[None, None], conv=1e-05)

Helper method for getting the left and right contact self-energies

Parameters:

• Elist (tuple, optional) – A list of contact energies for selecting sigma, (default: use contact ermi energy)

• conv (float, optional) – Convergence criterial for the self-energy matrix

Returns:

A tuple of both self-energy matrices (ndarrays)

Return type:

tuple

readBetheParams(filename)

Read Slater-Koster parameters from a .bethe file.

Reads and validates parameters for minimal basis with single s, p, and d orbitals. Parameters are stored in dictionaries for onsite energies, hopping integrals, and overlap matrices.

Parameters:

filename (str) – Name of the .bethe file (without extension)

Raises:

AssertionError – If parameters are missing or invalid

Notes

Parameters are sorted into: - Edict: Onsite energies (converted from Hartrees to eV) - Vdict: Hopping parameters (converted from Hartrees to eV) - Sdict: Overlap parameters

runAllTests()

Run all validation tests for surfGB.

Executes all test methods to validate: - d orbital angular functions - d orbital symmetry - p-d interactions - d-d interactions - General hopping physics

setF(F, muL, muR)

Update Fock matrix and contact chemical potentials.

Sets the Fock matrix and updates the Fermi levels of the left and right contacts if they have changed.

Parameters:

• F (ndarray) – New Fock matrix

• muL (float) – Chemical potential for left contact in eV

• muR (float) – Chemical potential for right contact in eV

sigma(E, i, conv=1e-05)

Calculate self-energy matrix for a specific contact.

Computes the self-energy matrix for contact i by: 1. Calculating surface self-energies for all 9 directions 2. Summing contributions from directions not connected to the device 3. Applying de-orthonormalization if needed 4. Handling spin configurations

Parameters:

• E (float) – Energy point for self-energy calculation (in eV)

• i (int) – Index of the contact to calculate self-energy for

• conv (float, optional) – Convergence criterion for self-energy calculation (default: 1e-5)

Returns:

Self-energy matrix for the specified contact, with dimensions: - (N, N) for restricted calculations - (2N, 2N) for unrestricted or generalized spin calculations

Return type:

ndarray

References

[1] Jacob, D., & Palacios, J. J. (2011). Critical comparison of electrode models

in density functional theory based quantum transport calculations. The Journal of Chemical Physics, 134(4), 044118. DOI: 10.1063/1.3526044

sigmaTot(E, conv=1e-05)

Calculate total self-energy matrix for the extended system.

Computes self-energies for all sites in the extended system (12 neighbors + 1 center). The total self-energy is constructed following the Bethe lattice model described in Jacob & Palacios [1], which provides an efficient representation of bulk metallic electrodes while maintaining proper orbital symmetries.

Parameters:

• E (float) – Energy point for Green’s function calculation (in eV)

• conv (float, optional) – Convergence criterion for self-energy calculation (default: 1e-5)

Returns:

Total self-energy matrix for the extended system

Return type:

ndarray

References

[1] Jacob, D., & Palacios, J. J. (2011). Critical comparison of electrode models

in density functional theory based quantum transport calculations. The Journal of Chemical Physics, 134(4), 044118. DOI: 10.1063/1.3526044

testDDInteraction()

Test d-d orbital interactions.

Validates d-d orbital interactions by checking: - dyz-dyz interaction along x-axis (should be pure delta) - dz2-dz2 interaction along z-axis (should be pure sigma)

testDOrbitalFunctions()

Test d orbital angular functions.

Validates the angular dependence of d orbital interactions by checking: - dxy interaction along x-axis (should be zero) - dx2-y2 interaction along x-axis (should be sqrt(3)/2 * sds) - dz2 interaction along x-axis (should be -1/2 * sds)

testDOrbitalSymmetry()

Test d orbital symmetry properties.

Validates that d orbital interactions respect inversion symmetry by comparing interactions along opposite directions.

testHoppingPhysics()

Test physical properties of hopping matrices.

Validates hopping matrix physics by checking: - s-p hopping antisymmetry - Conservation of total s-p hopping magnitude - Proper angular dependence along principal axes and 45-degree rotations

testPDInteraction()

Test p-d orbital interactions.

Validates p-d orbital interactions by checking: - px-dxy interaction along x-axis (should be zero) - pz-dz2 interaction along z-axis (should be pure sigma)

updateFermi(i, Ef)

Update Fermi energy for a specific contact.

Shifts the Hamiltonian of contact i to align its Fermi level with the specified energy.

Parameters:

• i (int) – Contact index

• Ef (float) – New Fermi energy in eV

class gauNEGF.surfGBethe.surfGBAt(H, Slist, Vlist, eta, T=0)

Bases: object

Atomic-level Bethe lattice Green’s function calculator.

This class implements the surface Green’s function calculation for a single atom in the Bethe lattice, handling: - Onsite and hopping matrix construction - Self-energy calculations for bulk and surface - Temperature effects - Fermi energy optimization

Parameters:

• H (ndarray) – Onsite Hamiltonian matrix (9x9 for minimal basis)

• Slist (list of ndarray) – List of 12 overlap matrices for nearest neighbors

• Vlist (list of ndarray) – List of 12 hopping matrices for nearest neighbors

• eta (float) – Broadening parameter in eV

• T (float, optional) – Temperature in Kelvin (default: 0)

NN

Number of nearest neighbors (fixed to 12 for FCC)

Type:

int

sigmaKprev

Previous bulk self-energy for convergence

Type:

ndarray or None

Eprev

Previous energy point for convergence

Type:

float

fermi

Current Fermi energy

Type:

float

F

Extended Fock matrix including neighbors

Type:

ndarray

S

Extended overlap matrix including neighbors

Type:

ndarray

DOS(E)

Calculate bulk density of states of the Bethe lattice.

Parameters:

E (float) – Energy point for DOS calculation (in eV)

Returns:

Density of states at energy E

Return type:

float

__init__(H, Slist, Vlist, eta, T=0)

Initialize surfGBAt with Hamiltonian and neighbor matrices.

Parameters:

• H (ndarray) – Onsite Hamiltonian matrix (9x9 for minimal basis)

• Slist (list of ndarray) – List of 12 overlap matrices for nearest neighbors

• Vlist (list of ndarray) – List of 12 hopping matrices for nearest neighbors

• eta (float) – Broadening parameter in eV

• T (float, optional) – Temperature in Kelvin (default: 0)

Raises:

AssertionError – If matrix dimensions are incorrect or number of neighbors != 12

calcFermi(ne, fGuess=5, tol=1e-05)

Calculate Fermi energy using bisection method.

Uses getFermiContact from density.py to find the Fermi energy that gives the correct number of electrons.

Parameters:

• ne (float) – Target number of electrons

• fGuess (float, optional) – Initial guess for Fermi energy in eV (default: 5)

• tol (float, optional) – Convergence tolerance (default: 1e-5)

Returns:

Calculated Fermi energy in eV

Return type:

float

Notes

Previous implementation used ANT.Gaussian approach with complex contour integration. Current version uses simpler bisection method from density.py.

setF(F, mu1, mu2)

Empty function for compatibility with density.py methods.

Bethe lattice bulk properties are intrinsic (dependent on TB parameters).

Parameters:

• F (ndarray) – Fock matrix (unused)

• mu1 (float) – First chemical potential (unused)

• mu2 (float) – Second chemical potential (unused)

sigma(E, inds=None, conv=1e-05, mix=0.5)

Calculate surface self-energies for an FCC lattice.

Computes self-energies for atoms at the surface with the geometry: [3x plane dir] - o - [3x plane dir]

/|

[3x out of plane dir]

Uses a self-consistent iteration scheme with mixing to solve the Dyson equation. The implementation follows the Bethe lattice approach described in Jacob & Palacios (2011), where the self-energy is computed recursively for a semi-infinite tree-like structure that preserves the proper coordination number and orbital symmetries of bulk FCC metals.

Parameters:

• E (float) – Energy point for Green’s function calculation (in eV)

• inds (list or int, optional) – Indices of the sigma matrix to return. If None, returns full list (default: None)

• conv (float, optional) – Convergence criterion for Dyson equation (default: 1e-5)

• mix (float, optional) – Mixing factor for Dyson equation (default: 0.5)

Returns:

List of self-energy matrices for the surface atom. If inds is specified, returns only the requested matrices.

Return type:

list

Notes

First calculates bulk self-energies using sigmaK, then iterates to find surface self-energies for the 9 surface directions. The recursive method ensures proper treatment of the metal-molecule interface while maintaining computational efficiency.

References

[1] Jacob, D., & Palacios, J. J. (2011). Critical comparison of electrode models

in density functional theory based quantum transport calculations. The Journal of Chemical Physics, 134(4), 044118. DOI: 10.1063/1.3526044

sigmaK(E, conv=1e-05, mix=0.5)

Calculate bulk self-energies for all 12 lattice directions.

Computes self-energies for an FCC lattice with the following geometry:

[3x out of plane dir]

|/

[3x plane dir] - o - [3x plane dir]

/|

[3x out of plane dir]

Uses a self-consistent iteration scheme with mixing to solve the Dyson equation.

Parameters:

• E (float) – Energy point for Green’s function calculation (in eV)

• conv (float, optional) – Convergence criterion for Dyson equation (default: 1e-5)

• mix (float, optional) – Mixing factor for Dyson equation (default: 0.5)

Returns:

Array of 12 self-energy matrices (9x9 each) in order by lattice direction

Return type:

ndarray

Notes

Uses previous solution as initial guess when energy point is close to previous calculation to improve convergence.

sigmaTot(E, conv=1e-05)

updateH(fermi=None)

Update Hamiltonian and extended matrices.

Updates onsite and hopping matrices, as well as extended lattice matrices. The extended matrices H0x and S0x include 13 sites total (12 neighbor sites followed by 1 onsite term). These are stored as F and S for compatibility with density.py functions.

Parameters:

fermi (float, optional) – New Fermi energy setpoint in eV (default: None)

Notes

When fermi is provided and different from current value: - Shifts onsite energies by the Fermi level difference - Updates hopping matrices with overlap contributions - Rebuilds extended matrices for the full system

1D Chain

Surface Green’s function implementation for 1D chain contacts.

This module provides a 1D chain implementation for modeling semi-infinite contacts in quantum transport calculations. It supports three usage patterns:

1. Fully automatic extraction from Fock matrix: surfG1D(F, S, [[contact1], [contact2]], [[contact1connection], [contact2connection]]) - All parameters extracted from F/S using contact and connection indices

2. Fock matrix with custom coupling: surfG1D(F, S, [[contact1], [contact2]], [tau1, tau2], [stau1, stau2]) - Contact parameters from F/S, but with custom coupling matrices

3. Fully specified contacts: surfG1D(F, S, [[contact1], [contact2]], [tau1, tau2], [stau1, stau2],

   [alpha1, alpha2], [salpha1, salpha2], [beta1, beta2], [sbeta1, sbeta2])

   • All contact parameters specified manually

The implementation uses an iterative scheme to calculate surface Green’s functions for 1D chain contacts, with support for both manual parameter specification and automatic extraction from DFT calculations.

class gauNEGF.surfG1D.surfG(Fock, Overlap, indsList, taus=None, staus=None, alphas=None, aOverlaps=None, betas=None, bOverlaps=None, eta=1e-09)

Bases: object

Surface Green’s function calculator for 1D chain contacts.

This class implements the surface Green’s function calculation for 1D chain contacts. It supports three usage patterns:

1. Fully automatic extraction from Fock matrix: - Provide contact indices and connection indices - All parameters extracted from F/S matrices Example: surfG1D(F, S, [[c1], [c2]], [[c1conn], [c2conn]])

2. Fock matrix with custom coupling: - Provide contact indices and coupling matrices - Onsite contact parameters from F/S, coupling specified manually Example: surfG1D(F, S, [[c1], [c2]], [tau1, tau2], [stau1, stau2])

3. Fully specified contacts: - All contact parameters provided manually Example: surfG1D(F, S, [[c1], [c2]], [tau1, tau2], [stau1, stau2],

   [alpha1, alpha2], [salpha1, salpha2], [beta1, beta2], [sbeta1, sbeta2])

Parameters:

• Fock (ndarray) – Fock matrix for the extended system

• Overlap (ndarray) – Overlap matrix for the extended system

• indsList (list of lists) – Lists of orbital indices for each contact region

• taus (list or None, optional) – Either coupling matrices or connection indices (default: None) - If indices: [[contact1connection], [contact2connection]] - If matrices: [tau1, tau2]

• staus (list or None, optional) – Overlap matrices for coupling, required if taus are matrices (default: None)

• alphas (list of ndarray or None, optional) – On-site energies for contacts, required for pattern (c) (default: None)

• aOverlaps (list of ndarray or None, optional) – On-site overlap matrices for contacts, required for pattern (c) (default: None)

• betas (list of ndarray or None, optional) – Hopping matrices between contact unit cells, required for pattern (c) (default: None)

• bOverlaps (list of ndarray or None, optional) – Overlap matrices between contact unit cells, required for pattern (c) (default: None)

• eta (float, optional) – Broadening parameter in eV (default: 1e-9)

F

Fock matrix

Type:

ndarray

S

Overlap matrix

Type:

ndarray

X

Inverse square root of overlap matrix for orthogonalization

Type:

ndarray

tauList

Contact coupling matrices

Type:

list

stauList

Contact coupling overlap matrices

Type:

list

aList

On-site energy matrices for contacts

Type:

list

aSList

On-site overlap matrices for contacts

Type:

list

bList

Hopping matrices between contact unit cells

Type:

list

bSList

Overlap matrices between contact unit cells

Type:

list

gPrev

Previous surface Green’s functions for convergence

Type:

list

__init__(Fock, Overlap, indsList, taus=None, staus=None, alphas=None, aOverlaps=None, betas=None, bOverlaps=None, eta=1e-09)

Initialize the surface Green’s function calculator.

The initialization follows one of three patterns: a) Fully automatic: Only provide Fock, Overlap, indsList, and connection indices in taus b) Custom coupling: Provide Fock, Overlap, indsList, coupling matrices in taus, and staus c) Fully specified: Provide all parameters including alphas, aOverlaps, betas, bOverlaps

Parameters:

• Fock (ndarray) – Fock matrix for the extended system

• Overlap (ndarray) – Overlap matrix for the extended system

• indsList (list of lists) – Lists of orbital indices for each contact region

• taus (list or None, optional) – Either coupling matrices or connection indices (default: None) - If indices: [[contact1connection], [contact2connection]] - If matrices: [tau1, tau2]

• staus (list or None, optional) – Overlap matrices for coupling, required if taus are matrices (default: None)

• alphas (list of ndarray or None, optional) – On-site energies for contacts, required for pattern (c) (default: None)

• aOverlaps (list of ndarray or None, optional) – On-site overlap matrices for contacts, required for pattern (c) (default: None)

• betas (list of ndarray or None, optional) – Hopping matrices between contact unit cells, required for pattern (c) (default: None)

• bOverlaps (list of ndarray or None, optional) – Overlap matrices between contact unit cells, required for pattern (c) (default: None)

• eta (float, optional) – Broadening parameter in eV (default: 1e-9)

Notes

The initialization will raise an error if: - The parameters don’t match one of the three usage patterns - taus contains matrices but staus is None - alphas is provided but aOverlaps is None - betas is provided but bOverlaps is None

denFunc(E, ind=None, mu=None, T=300)

Calculate density matrix contribution at a single energy point.

Computes the contribution to the density matrix at energy E for specified orbitals. This is used by integration routines to calculate the total density matrix.

Parameters:

• E (float) – Energy point in eV

• ind (int or None, optional) – Contact index for partial density calculation (default: None)

• mu (float or None, optional) – Chemical potential in eV (default: None)

• T (float, optional) – Temperature in Kelvin (default: 300)

Returns:

Density matrix contribution at energy E

Return type:

ndarray

densityGrid(Emin, Emax, ind=None, dE=0.001, mu=None, T=300)

Calculate density matrix contribution on a grid of energy points.

Computes the contribution to the density matrix for each energy point in a grid from Emin to Emax. This is used for real-axis integration of the density matrix.

Parameters:

• Emin (float) – Minimum energy for integration in eV

• Emax (float) – Maximum energy for integration in eV

• ind (int or None, optional) – Contact index for partial density calculation (default: None)

• dE (float, optional) – Energy step size in eV (default: 0.001)

• mu (float or None, optional) – Chemical potential in eV (default: None)

• T (float, optional) – Temperature in Kelvin (default: 300)

Returns:

Array of density matrix contributions at each energy point

Return type:

ndarray

g(E, i, conv=1e-05, relFactor=0.1)

Calculate surface Green’s function for a contact.

Uses an iterative scheme to calculate the surface Green’s function for contact i at energy E. The iteration continues until the change in the Green’s function is below the convergence criterion.

Parameters:

• E (float) – Energy point in eV

• i (int) – Contact index

• conv (float, optional) – Convergence criterion for iteration (default: 1e-5)

• relFactor (float, optional) – Relaxation factor for iteration mixing (default: 0.1)

Returns:

Surface Green’s function matrix for contact i

Return type:

ndarray

Notes

The method uses the previous solution as an initial guess to improve convergence. For the first calculation at a given energy, it uses zeros as the initial guess. The relaxation factor controls mixing between iterations to help convergence.

setContacts(alphas=None, aOverlaps=None, betas=None, bOverlaps=None)

Update contact parameters for the 1D chain.

This method is used internally during initialization and can be called later to update contact parameters. It follows the same patterns as initialization:

1. If self.contactFromFock is True (patterns a and b): - Parameters are extracted from F/S matrices - Any provided parameters are ignored

2. If self.contactFromFock is False (pattern c): - All parameters must be provided together - Partial updates are not supported

Parameters:

• alphas (list of ndarray or None, optional) – On-site energies for contacts (default: None)

• aOverlaps (list of ndarray or None, optional) – On-site overlap matrices for contacts (default: None)

• betas (list of ndarray or None, optional) – Hopping matrices between contact unit cells (default: None)

• bOverlaps (list of ndarray or None, optional) – Overlap matrices between contact unit cells (default: None)

Notes

When using pattern (c), all parameters must be provided together. Partial updates (providing some parameters but not others) are not supported and will raise an error.

setF(F, mu1=None, mu2=None)

Update the Fock matrix and contact chemical potentials.

This method updates the system’s Fock matrix and optionally shifts the contact chemical potentials. If the contacts are extracted from the Fock matrix, their parameters are automatically updated.

Parameters:

• F (ndarray) – New Fock matrix for the system

• mu1 (float or None, optional) – Chemical potential for first contact in eV (default: None)

• mu2 (float or None, optional) – Chemical potential for second contact in eV (default: None)

Notes

If chemical potentials are provided, the corresponding contact parameters are shifted to align with the new potentials.

sigma(E, i, conv=1e-05)

Calculate self-energy matrix for a contact.

Computes the self-energy matrix for contact i at energy E using the surface Green’s function. The self-energy represents the effect of the semi-infinite contact on the device region.

Parameters:

• E (float) – Energy point in eV

• i (int) – Contact index

• conv (float, optional) – Convergence criterion for surface Green’s function (default: 1e-5)

Returns:

Self-energy matrix for contact i

Return type:

ndarray

sigmaTot(E, conv=1e-05)

Calculate total self-energy matrix from all contacts.

Computes the total self-energy matrix at energy E by summing contributions from all contacts. This represents the combined effect of all semi-infinite contacts on the device region.

Parameters:

• E (float) – Energy point in eV

• conv (float, optional) – Convergence criterion for surface Green’s functions (default: 1e-5)

Returns:

Total self-energy matrix from all contacts

Return type:

ndarray

Constant Self Energy

class gauNEGF.surfGTester.surfGTest(Fock, Overlap, indsList, sig1=None, sig2=None)

Bases: object

setF(F, mu1, mu2)

sigma(E, i, conv=0.001)

sigmaTot(E, conv=0.001)

Utilities

Matrix Tools

Matrix manipulation utilities for quantum transport calculations.

This module provides helper functions for handling matrices in quantum transport calculations, with a focus on: - Self-energy matrix construction - Density and Fock matrix manipulation - Spin treatment (restricted, unrestricted, and generalized) - Integration with Gaussian’s matrix formats

The functions handle three types of spin treatments: - ‘r’: Restricted (same orbitals for alpha and beta electrons) - ‘ro’/’u’: Unrestricted (separate alpha and beta orbitals) - ‘g’: Generalized (non-collinear spin treatment)

gauNEGF.matTools.formSigma(inds, V, nsto, S=0)

Form a self-energy matrix for specified orbitals.

Creates a self-energy matrix of size nsto x nsto with values V at the specified orbital indices. Can handle both scalar and matrix-valued self-energies.

Parameters:

• inds (list of int) – Orbital indices where self-energy should be applied

• V (complex or ndarray) – Self-energy value(s) to insert - If scalar: Same value used for all specified orbitals - If matrix: Must match size of indices

• nsto (int) – Total number of orbitals (size of resulting matrix)

• S (ndarray or int, optional) – Overlap matrix for broadening term. If 0, identity used (default: 0)

Returns:

Complex self-energy matrix of size nsto x nsto

Return type:

ndarray

gauNEGF.matTools.getDen(bar, spin)

Build density matrix from Gaussian checkpoint file.

Extracts the density matrix from a Gaussian checkpoint file based on the specified spin treatment. Handles restricted, unrestricted, and generalized cases.

Parameters:

• bar (QCBinAr) – Gaussian interface object

• spin (str) – Spin treatment to use: - ‘r’: Restricted (same orbitals for alpha/beta) - ‘ro’/’u’: Unrestricted (separate alpha/beta) - ‘g’: Generalized (non-collinear)

Returns:

Density matrix in appropriate format for spin treatment: - Restricted: Single block - Unrestricted: Two blocks (alpha/beta) - Generalized: Single block with complex entries

Return type:

ndarray

Raises:

ValueError – If spin treatment is not recognized

gauNEGF.matTools.getEnergies(bar, spin)

Get orbital energies from Gaussian checkpoint file.

Extracts orbital energies from a Gaussian checkpoint file based on the specified spin treatment. Converts energies from Hartrees to eV.

Parameters:

• bar (QCBinAr) – Gaussian interface object

• spin (str) – Spin treatment to use: - ‘r’: Restricted (same orbitals for alpha/beta) - ‘ro’/’u’: Unrestricted (separate alpha/beta) - ‘g’: Generalized (non-collinear)

Returns:

Array of orbital energies in eV, sorted in ascending order. Format depends on spin treatment: - Restricted: Alternating alpha/beta pairs - Unrestricted: Alternating alpha/beta pairs - Generalized: Single set of energies

Return type:

ndarray

Raises:

ValueError – If spin treatment is not recognized

gauNEGF.matTools.getFock(bar, spin)

Build Fock matrix from Gaussian checkpoint file.

Extracts the Fock matrix and orbital indices from a Gaussian checkpoint file based on the specified spin treatment. Handles restricted, unrestricted, and generalized cases.

Parameters:

• bar (QCBinAr) – Gaussian interface object

• spin (str) – Spin treatment to use: - ‘r’: Restricted (same orbitals for alpha/beta) - ‘ro’/’u’: Unrestricted (separate alpha/beta) - ‘g’: Generalized (non-collinear)

Returns:

(Fock, locs) where: - Fock: ndarray, Fock matrix in appropriate format for spin treatment - locs: ndarray, Orbital indices with positive for alpha/paired and

negative for beta orbitals

Return type:

tuple

Raises:

ValueError – If spin treatment is not recognized

gauNEGF.matTools.storeDen(bar, P, spin)

Store density matrix in Gaussian checkpoint format.

Converts the density matrix to the appropriate format based on spin treatment and stores it in the Gaussian checkpoint file. Handles compression and proper typing of matrices.

Parameters:

• bar (QCBinAr) – Gaussian interface object

• P (ndarray) – Density matrix to store

• spin (str) – Spin treatment to use: - ‘r’: Restricted (same orbitals for alpha/beta) - ‘ro’/’u’: Unrestricted (separate alpha/beta) - ‘g’: Generalized (non-collinear)

Notes

For restricted calculations, the density is divided by 2 to account for double occupation. For generalized calculations, complex matrices are used.

Raises:

ValueError – If spin treatment is not recognized