Skip to content

Harrison Tight-Binding (universal parameters)

Dependency-free (numpy only) sp3 tight-binding Hamiltonian builder for real atoms and crystals -- no PySCF/OpenFermion, no SCF/DFT. See vhd_tb for the material-specific alternative when this module's accuracy isn't enough.

Source

Walter A. Harrison, Electronic Structure and the Properties of Solids: The Physics of the Chemical Bond. Originally published by W. H. Freeman, 1980; reprinted by Dover Publications (Dover Books on Physics), 1989, ISBN 0-486-66021-4. Atomic term values (ELEMENTS) and the universal eta coefficients (ETA) are transcribed from that book's Solid State Table, cross-checked against jarvist/HarrisonSolidStateTable.jl, an independent Julia implementation of the same table.

The method

Harrison's tight-binding model builds a solid's electronic Hamiltonian from two ingredients only, both universal (materials-independent functional form):

  1. Atomic term values -- the free-atom s and p orbital energies (on-site Hamiltonian diagonal), tabulated per element.
  2. A universal bond-scaling law for the off-diagonal (hopping) matrix elements between neighboring atoms' orbitals:

$\(V_{ll'm} = \eta_{ll'm} \cdot \frac{\hbar^2}{m_e d^2}\)$

where \(d\) is the bond length and the four dimensionless \(\eta\) coefficients are the same for every element pair -- only \(d\) and the atomic term values change between materials. This is what makes the method "universal": no fitting per material.

\(\hbar^2/m_e = 7.62\ \text{eV·Å}^2\).

coefficient value
\(\eta_{ss\sigma}\) -1.40
\(\eta_{sp\sigma}\) +1.84
\(\eta_{pp\sigma}\) +3.24
\(\eta_{pp\pi}\) -0.81

Off-diagonal sp3 matrix elements follow the standard Slater-Koster (1954) table for an (s, px, py, pz) basis and a bond of direction cosines \((l, m, n)\):

\[E(s,s) = V_{ss\sigma}, \quad E(s,x) = l\,V_{sp\sigma}, \quad E(x,s) = -l\,V_{sp\sigma}$$ $$E(x,x) = l^2 V_{pp\sigma} + (1-l^2)V_{pp\pi}, \quad E(x,y) = lm\,(V_{pp\sigma}-V_{pp\pi})\]

(and cyclic permutations for y, z). sp3_bond_block implements this; sp3_dimer_hamiltonian builds a 2-atom cluster from it, and zincblende_hamiltonian sums it with Bloch phases over the 4 nearest-neighbor bonds to build the full periodic crystal Hamiltonian.

Accuracy

This is a universal model -- one parameter table for every material, no per-material fitting, no d-orbitals. That buys zero setup cost per new material at the price of accuracy: gaps typically come out ~2-3x off from experiment, and for indirect-gap materials it can misplace the conduction-band minimum entirely (see vhd_tb for why and the fix). Validation numbers against real experimental gaps (GaAs, Si, Ge) are tracked in Dense-Evolution-Discovery.

harrison_tb

Harrison empirical tight-binding parameters -- builds an sp3 tight-binding Hamiltonian for a cluster of real atoms directly from published atomic term values and the universal bond-scaling law, with no SCF/DFT and no external quantum-chemistry dependency (PySCF, OpenFermion).

Source: Walter A. Harrison, "Electronic Structure and the Properties of Solids" (Dover reprint). Atomic term values (ELEMENTS) and the universal eta coefficients (ETA) below are transcribed from that book's Solid State Table, cross-checked against the numbers in jarvist/HarrisonSolidStateTable.jl (github.com/jarvist/HarrisonSolidStateTable.jl), a Julia implementation of the same table. Only elements with both an s and a p term value in that table are included here (the "simple atom" sp3 entries Harrison uses for tetrahedral semiconductors); d-block entries are left out since Harrison's own table flags them "not well checked".

Sign convention: ELEMENTS stores true orbital energies (negative, eV -- bound states below vacuum), i.e. the negative of the magnitudes printed in Harrison's table.

Two-center bond integrals follow Harrison's universal scaling law: V_ll'm = eta_ll'm * hbar^2 / (m_e * d^2) with d the bond length in Angstrom and hbar^2/m_e = 7.62 eV*Angstrom^2 (the standard constant quoted alongside this law). eta_ssσ/spσ/ppσ/ppπ are dimensionless and materials-independent -- the same four numbers apply to every element pair.

Slater-Koster sp3 matrix elements (sp3_bond_block) follow the standard 1954 Slater-Koster table for an (s, px, py, pz) basis.

zincblende_hamiltonian builds the periodic Bloch Hamiltonian for a two-atom zinc-blende basis (nearest-neighbor sp3, 4 bonds per atom) -- validated against real GaAs (a=5.6533 Angstrom): computed direct gap at Gamma is 2.91 eV vs. the experimental 1.42 eV, roughly 2x too large. This is a known, documented limitation of Harrison's universal (materials-independent) parameter set on polar/ionic compound semiconductors -- not a bug here -- since it uses no per-material fitting and omits d-orbitals. Useful as a fast, dependency-free qualitative estimate; not a substitute for this project's DFT-derived GaAs parameters where quantitative accuracy matters.

hopping_integral

hopping_integral(eta, d_angstrom)

V_ll'm = eta * hbar^2/(m_e d^2) [eV], for bond length d in Angstrom.

Source code in dense_evolution/harrison_tb.py
def hopping_integral(eta, d_angstrom):
    """V_ll'm = eta * hbar^2/(m_e d^2) [eV], for bond length d in Angstrom."""
    if d_angstrom <= 0:
        raise ValueError(f"bond length must be positive, got {d_angstrom}")
    return eta * HBAR2_OVER_M_EV_ANG2 / d_angstrom ** 2

sp3_bond_block

sp3_bond_block(l, m, n, d_angstrom, eta=ETA)

4x4 Slater-Koster hopping block for a bond from atom A to atom B along direction cosines (l, m, n) (unit vector, l^2+m^2+n^2 = 1) and bond length d_angstrom.

Basis order: s, px, py, pz. Standard Slater-Koster (1954) table.

Source code in dense_evolution/harrison_tb.py
def sp3_bond_block(l, m, n, d_angstrom, eta=ETA):
    """
    4x4 Slater-Koster hopping block <A, {s,px,py,pz}| H |B, {s,px,py,pz}>
    for a bond from atom A to atom B along direction cosines (l, m, n)
    (unit vector, l^2+m^2+n^2 = 1) and bond length d_angstrom.

    Basis order: s, px, py, pz. Standard Slater-Koster (1954) table.
    """
    norm = l * l + m * m + n * n
    if not np.isclose(norm, 1.0, atol=1e-6):
        raise ValueError(f"(l, m, n) must be a unit vector, got norm={norm}")

    Vssσ = hopping_integral(eta['ss_sigma'], d_angstrom)
    Vspσ = hopping_integral(eta['sp_sigma'], d_angstrom)
    Vppσ = hopping_integral(eta['pp_sigma'], d_angstrom)
    Vppπ = hopping_integral(eta['pp_pi'], d_angstrom)

    block = np.zeros((4, 4), dtype=np.complex128)
    block[0, 0] = Vssσ
    block[0, 1], block[0, 2], block[0, 3] = l * Vspσ, m * Vspσ, n * Vspσ
    block[1, 0], block[2, 0], block[3, 0] = -l * Vspσ, -m * Vspσ, -n * Vspσ

    block[1, 1] = l * l * Vppσ + (1 - l * l) * Vppπ
    block[2, 2] = m * m * Vppσ + (1 - m * m) * Vppπ
    block[3, 3] = n * n * Vppσ + (1 - n * n) * Vppπ

    block[1, 2] = block[2, 1] = l * m * (Vppσ - Vppπ)
    block[2, 3] = block[3, 2] = m * n * (Vppσ - Vppπ)
    block[1, 3] = block[3, 1] = l * n * (Vppσ - Vppπ)
    return block

sp3_dimer_hamiltonian

sp3_dimer_hamiltonian(
    element_a,
    element_b,
    bond_length_angstrom,
    direction=(0.0, 0.0, 1.0),
    eta=ETA,
)

8x8 sp3 tight-binding Hamiltonian for a 2-atom A-B cluster (one bond), basis order [A:s,px,py,pz, B:s,px,py,pz]. On-site blocks are each atom's diagonal (eps_s, eps_p, eps_p, eps_p); the A-B off-diagonal block is sp3_bond_block along direction (unit vector, default: bond along z), Hermitian-conjugated into the B-A block.

This is a minimal, directly checkable sanity case (bonding/antibonding sp3 splitting), not a periodic-solid band structure.

Source code in dense_evolution/harrison_tb.py
def sp3_dimer_hamiltonian(element_a, element_b, bond_length_angstrom,
                           direction=(0.0, 0.0, 1.0), eta=ETA):
    """
    8x8 sp3 tight-binding Hamiltonian for a 2-atom A-B cluster (one bond),
    basis order [A:s,px,py,pz, B:s,px,py,pz]. On-site blocks are each
    atom's diagonal (eps_s, eps_p, eps_p, eps_p); the A-B off-diagonal
    block is sp3_bond_block along `direction` (unit vector, default: bond
    along z), Hermitian-conjugated into the B-A block.

    This is a minimal, directly checkable sanity case (bonding/antibonding
    sp3 splitting), not a periodic-solid band structure.
    """
    for name in (element_a, element_b):
        if name not in ELEMENTS:
            raise ValueError(f"no Harrison sp term values for element {name!r}; "
                              f"available: {sorted(ELEMENTS)}")
    l, m, n = direction
    a, b = ELEMENTS[element_a], ELEMENTS[element_b]

    H = np.zeros((8, 8), dtype=np.complex128)
    H[0, 0] = a['eps_s']
    H[1, 1] = H[2, 2] = H[3, 3] = a['eps_p']
    H[4, 4] = b['eps_s']
    H[5, 5] = H[6, 6] = H[7, 7] = b['eps_p']

    hop = sp3_bond_block(l, m, n, bond_length_angstrom, eta=eta)
    H[0:4, 4:8] = hop
    H[4:8, 0:4] = hop.conj().T
    return H

zincblende_hamiltonian

zincblende_hamiltonian(
    k, cation, anion, lattice_constant_angstrom, eta=ETA
)

8x8 Bloch Hamiltonian for a zinc-blende crystal's two-atom basis (cation at (0,0,0), anion at (1/4,1/4,1/4) of the conventional cubic cell), sp3 nearest-neighbor tight-binding, at crystal momentum k (Cartesian, 1/Angstrom -- e.g. Gamma=(0,0,0)).

Basis order [cation:s,px,py,pz, anion:s,px,py,pz]. The four cation->anion nearest-neighbor bonds are the standard zinc-blende tetrahedral set (lattice_constant/4)*(1,1,1), (1,-1,-1), (-1,1,-1), (-1,-1,1); each contributes sp3_bond_block(...) weighted by its Bloch phase exp(i k . d), summed into the off-diagonal block.

Source code in dense_evolution/harrison_tb.py
def zincblende_hamiltonian(k, cation, anion, lattice_constant_angstrom, eta=ETA):
    """
    8x8 Bloch Hamiltonian for a zinc-blende crystal's two-atom basis
    (cation at (0,0,0), anion at (1/4,1/4,1/4) of the conventional cubic
    cell), sp3 nearest-neighbor tight-binding, at crystal momentum k
    (Cartesian, 1/Angstrom -- e.g. Gamma=(0,0,0)).

    Basis order [cation:s,px,py,pz, anion:s,px,py,pz]. The four
    cation->anion nearest-neighbor bonds are the standard zinc-blende
    tetrahedral set (lattice_constant/4)*(1,1,1), (1,-1,-1), (-1,1,-1),
    (-1,-1,1); each contributes sp3_bond_block(...) weighted by its
    Bloch phase exp(i k . d), summed into the off-diagonal block.
    """
    for name in (cation, anion):
        if name not in ELEMENTS:
            raise ValueError(f"no Harrison sp term values for element {name!r}; "
                              f"available: {sorted(ELEMENTS)}")
    k = np.asarray(k, dtype=float)
    d_vectors = (lattice_constant_angstrom / 4) * np.array([
        [1, 1, 1], [1, -1, -1], [-1, 1, -1], [-1, -1, 1],
    ], dtype=float)
    bond_length = np.linalg.norm(d_vectors[0])

    T = np.zeros((4, 4), dtype=np.complex128)
    for d in d_vectors:
        l, m, n = d / bond_length
        phase = np.exp(1j * np.dot(k, d))
        T += phase * sp3_bond_block(l, m, n, bond_length, eta=eta)

    c, a = ELEMENTS[cation], ELEMENTS[anion]
    H = np.zeros((8, 8), dtype=np.complex128)
    H[0, 0] = c['eps_s']
    H[1, 1] = H[2, 2] = H[3, 3] = c['eps_p']
    H[4, 4] = a['eps_s']
    H[5, 5] = H[6, 6] = H[7, 7] = a['eps_p']
    H[0:4, 4:8] = T
    H[4:8, 0:4] = T.conj().T
    return H

See also: vhd_tb for material-specific fitted parameters when the universal table's ~2-3x gap error isn't good enough.