Designing Custom Hoppings

In the previous section we defined custom layer geometries. The layer tells MoirePy where the atoms sit; the hopping parameters tell it how strongly they couple. By default you pass a constant float and every bond of that type gets the same amplitude. In this section we replace those constants with callables so the coupling varies bond-by-bond, driven by the actual distance between atoms.

1. The Six Hopping Slots

generate_hamiltonian accepts six parameters:

Parameter	Description
`tll`	Intra-layer hopping inside the lower layer
`tuu`	Intra-layer hopping inside the upper layer
`tlu`	Inter-layer hopping from lower to upper
`tul`	Inter-layer hopping from upper to lower
`tlself`	On-site energy for atoms in the lower layer
`tuself`	On-site energy for atoms in the upper layer

Each slot accepts either:

A float or int: every bond of that type gets the same constant amplitude.
Or a callable (a function for example): MoirePy calls it once per hopping channel and expects an array back.

You can freely mix: pass a float for tll and a function for tul. All six slots follow the same rule and default to zero (0) if not specified.

2. Function Signatures

The signature differs slightly depending on whether the parameter is a pair-hopping term (tll, tuu, tlu, tul) or an on-site term (tlself, tuself).

Pair-hopping (`tll`, `tuu`, `tlu`, `tul`)

def my_hopping(pos_i, pos_j, R, type_i, type_j, lattice, **kwargs):
    ...
    return values

Argument	Type	Description
`pos_i`	`ndarray (N, 2)`	Coordinates of the source atoms
`pos_j`	`ndarray (N, 2)`	Coordinates of the destination atoms. In case of PBC, this is the mapped coordinate of the atom inside the supercell.
`R`	`ndarray (N, 2)`	Supercell shift vectors for the point `j`. Add this to `pos_j` to get the actual coordinate of the atom.
`type_i`	`list[str]` of length N	Atom-type labels of source atoms (e.g. `"A"`, `"B"`)
`type_j`	`list[str]` of length N	Atom-type labels of destination atoms
`lattice`	`BilayerMoireLattice`	The full lattice object; has `theta`, `mlv1`, `mlv2`, etc.
`**kwargs`	any	Anything you pass via `extra_inputs`

N is the total number of bonds of that hopping type in the supercell.

Note

The shape of the returned value should be (N, k, k) where k is the number of orbitals per atom. For single-orbital (k=1) systems (N,) shape doesn't raise error.

On-site (`tlself`, `tuself`)

def my_onsite(pos_i, type_i, lattice, **kwargs):
    ...
    return values

R is absent here because on-site terms do not hop anywhere.

Note

This example only uses pos_i and pos_j to compute bond length. The remaining arguments are there when you need them: type_i/type_j lets you distinguish A->B from A->A bonds; R gives the periodic image offset of each hop; lattice exposes lattice.theta, lattice.mlv1, and the rest of the geometry.

3. The Hopping Model

We use a Gaussian decay for the inter-layer hopping:

\[t(r) = t_0 \, \exp\!\left(-\left(\frac{r}{r_c}\right)^2\right)\]

\(t_0\) is the maximum hopping amplitude at zero separation.
\(r_c\) is the decay length. Bonds much longer than \(r_c\) contribute almost nothing.

Orbital overlap decays rapidly with distance, and the Gaussian captures that behavior with just two parameters. We keep tll and tuu as plain floats (constant nearest-neighbour hopping within each layer) and put the Gaussian only on tul and tlu, which control the inter-layer coupling.

4. Setup

import numpy as np
import matplotlib.pyplot as plt
from moirepy import BilayerMoireLattice, SquareLayer

ll1, ll2, ul1, ul2 = 2, 3, 3, 2

lattice = BilayerMoireLattice(
    SquareLayer, ll1, ll2, ul1, ul2,
    n1=1, n2=1, pbc=True,
)
lattice.generate_connections(inter_layer_radius=3.0)

# Output:
# twist angle = 0.3948 rad (22.6199 deg)
# 13 cells in upper lattice
# 13 cells in lower lattice

Warning

Set inter_layer_radius comfortably above your largest r_c. Here r_c stays at or below 1.0 while the radius is 3.0, so every bond that matters falls within the search radius, and the Gaussian handles the suppression of distant ones.

5. Defining the Custom Function

def gaussian_inter(pos_i, pos_j, R, type_i, type_j, lattice, t0, rc):
    """Gaussian-decay inter-layer hopping."""
    # R maps pos_j to its true physical coordinate inside the supercell
    r = np.linalg.norm(pos_i - (pos_j + R), axis=1)  
    return t0 * np.exp(-(r / rc) ** 2)

pos_i - (pos_j + R) has shape (N, 2), so np.linalg.norm(..., axis=1) gives a length-N array of distances, one per bond. The function returns a length-N array, exactly what MoirePy expects.

6. Building the Hamiltonian

Pass t0 and rc through extra_inputs. Every key there becomes a keyword argument to your callable.

INTRA_T = -1.0  # constant nearest-neighbour hopping within each layer
T0      =  0.5  # inter-layer amplitude at r=0
RC      =  1.0  # Gaussian decay length

ham = lattice.generate_hamiltonian(
    tll=INTRA_T,           # float: same hop for every intra-lower bond
    tuu=INTRA_T,           # float: same hop for every intra-upper bond
    tul=gaussian_inter,    # callable
    tlu=gaussian_inter,    # callable (tlu = tul dagger for a Hermitian system)
    tlself=0.0,            # float: zero on-site energy
    tuself=0.0,            # float: zero on-site energy
    extra_inputs={"t0": T0, "rc": RC},
).toarray()

print(f"Hamiltonian shape: {ham.shape}")

# Output:
# Hamiltonian shape: (26, 26)

Hermiticity

MoirePy does not enforce tlu = tul* automatically. For a physical Hermitian system, make sure your tlu callable returns the complex conjugate of what tul returns, or pass the same real-valued function for both as we do here.

This setup can be naturally generalized to multiorbital and non-Hermitian systems through a careful choice of hopping amplitudes.

Summary

All six hopping slots accept a float, int, or any callable; you can mix and match freely.
Pair-hopping callables receive (pos_i, pos_j, R, type_i, type_j, lattice, **kwargs); on-site callables receive (pos_i, type_i, lattice, **kwargs).
Pass extra parameters to your callable via the extra_inputs dict.
Set inter_layer_radius larger than your largest r_c; let the hopping function handle suppression.
For a Hermitian system, ensure tlu is the complex conjugate of tul.

Next Steps

K-Space and Band Structures: Fold this Hamiltonian into momentum space and compute a band structure.
OBC vs PBC: See how boundary conditions change the spectrum.
Tutorials and Replicated Papers: See full physical results with realistic hop models.