Protocols¶

This page documents the Python protocols defined and used in ffsim.

In Python, a protocol defines an interface that a class can satisfy implicitly, without explicit inheritance. You may already be familiar with some standard protocols. For example, any class that defines a __len__ method implements the Sized protocol, and the built-in function len dispatches to it.

ffsim defines its own set of protocols that follow the same pattern. Each protocol specifies a single “dunder”-style method (a method whose name begins and ends with underscores, like _apply_unitary_). A class implements the protocol by defining this method. For each protocol, ffsim also provides a corresponding module-level function that dispatches to the protocol method on the given object.

The protocol-based design means that if you define your own class and implement the appropriate methods, it will work with ffsim’s functions without needing to inherit from any ffsim base class.

The table below summarizes all of ffsim’s protocols. The rest of this page describes each one in detail.

Protocol	Method	Function	Description
SupportsApplyUnitary	`_apply_unitary_`	ffsim.apply_unitary	Apply a unitary transformation to a state vector.
SupportsLinearOperator	`_linear_operator_`	ffsim.linear_operator	Convert to a SciPy `LinearOperator`.
SupportsFermionOperator	`_fermion_operator_`	ffsim.fermion_operator	Convert to a FermionOperator.
SupportsTrace	`_trace_`	ffsim.trace	Compute the trace of a linear operator.
SupportsDiagonal	`_diag_`	ffsim.diag	Return the diagonal entries of a linear operator.
SupportsApproximateEquality	`_approx_eq_`	ffsim.approx_eq	Compare objects with numerical tolerance.

[1]:

import numpy as np
import scipy.sparse.linalg

import ffsim

# Set up a small system for demonstrations
norb = 4
nelec = (2, 2)
rng = np.random.default_rng(12345)

SupportsApplyUnitary¶

The SupportsApplyUnitary protocol is for objects that represent a unitary transformation and can apply it to a state vector. A class implements this protocol by defining the _apply_unitary_ method, which has the signature

def _apply_unitary_(
    self, vec: np.ndarray, norb: int, nelec: int | tuple[int, int], copy: bool
) -> np.ndarray:
    ...

The function ffsim.apply_unitary(vec, obj, norb, nelec, copy=True) dispatches to this method.

Classes that implement this protocol include the variational ansatz operators such as UCJOpSpinBalanced and the UCCSD operator classes.

[2]:

# Create a UCJOpSpinBalanced, which implements SupportsApplyUnitary
operator = ffsim.random.random_ucj_op_spin_balanced(norb, n_reps=2, seed=rng)

# Apply the unitary to the Hartree-Fock state using the protocol function
vec = ffsim.hartree_fock_state(norb, nelec)
result = ffsim.apply_unitary(vec, operator, norb=norb, nelec=nelec)
result

[2]:

array([ 0.04007362-0.05264454j, -0.02721016-0.12343819j,
        0.04930984-0.07652777j, -0.04194643+0.07713257j,
       -0.03589734+0.10488243j,  0.10943586-0.006303j  ,
       -0.02721016-0.12343819j, -0.06439712-0.2474382j ,
        0.02990384-0.03644825j, -0.22082698-0.08395184j,
        0.06341983-0.04759223j,  0.17155003+0.05574341j,
        0.04930984-0.07652777j,  0.02990384-0.03644825j,
       -0.13034281+0.08251763j, -0.08186874+0.20012938j,
        0.07125724+0.26673363j, -0.14560141-0.10587217j,
       -0.04194643+0.07713257j, -0.22082698-0.08395184j,
       -0.08186874+0.20012938j, -0.19607996+0.21057831j,
        0.21819935-0.01279336j, -0.01216177+0.04060845j,
       -0.03589734+0.10488243j,  0.06341983-0.04759223j,
        0.07125724+0.26673363j,  0.21819935-0.01279336j,
       -0.16261186+0.07414021j, -0.04053588+0.1794134j ,
        0.10943586-0.006303j  ,  0.17155003+0.05574341j,
       -0.14560141-0.10587217j, -0.01216177+0.04060845j,
       -0.04053588+0.1794134j , -0.08791435-0.01156778j])

SupportsLinearOperator¶

The SupportsLinearOperator protocol is for objects that can be converted to a SciPy LinearOperator. A class implements this protocol by defining the _linear_operator_ method, which has the signature

def _linear_operator_(
    self, norb: int, nelec: int | tuple[int, int]
) -> scipy.sparse.linalg.LinearOperator:
    ...

The function ffsim.linear_operator(obj, norb, nelec) dispatches to this method. It also has special handling for FermionOperator instances, which are converted to a LinearOperator directly.

The LinearOperator representation enables matrix-vector multiplication, eigenvalue computation, and time evolution. The Hamiltonian classes (MolecularHamiltonian, DiagonalCoulombHamiltonian, DoubleFactorizedHamiltonian, SingleFactorizedHamiltonian) all implement this protocol.

[3]:

# Create a MolecularHamiltonian, which implements SupportsLinearOperator
mol_hamiltonian = ffsim.random.random_molecular_hamiltonian(norb, seed=rng, dtype=float)

# Convert to a LinearOperator using the protocol function
linop = ffsim.linear_operator(mol_hamiltonian, norb=norb, nelec=nelec)

# Use the LinearOperator to compute the ground state energy
eigs, _ = scipy.sparse.linalg.eigsh(linop, k=1, which="SA")
eigs[0]

[3]:

np.float64(-36.53619903781154)

SupportsFermionOperator¶

The SupportsFermionOperator protocol is for objects that can be converted to a FermionOperator, which is ffsim’s symbolic representation of a fermionic operator as a linear combination of products of creation and annihilation operators. A class implements this protocol by defining the _fermion_operator_ method, which has the signature

def _fermion_operator_(self) -> FermionOperator:
    ...

The function ffsim.fermion_operator(obj) dispatches to this method.

This protocol is implemented by the Hamiltonian classes.

[4]:

# Convert the molecular Hamiltonian to a FermionOperator
fermion_op = ffsim.fermion_operator(mol_hamiltonian)

# The result is a FermionOperator, a dict-like mapping from terms to coefficients
print(f"Number of terms: {len(fermion_op)}")

Number of terms: 1057

SupportsDiagonal¶

The SupportsDiagonal protocol is for linear operators whose diagonal entries can be returned. A class implements this protocol by defining the _diag_ method, which has the signature

def _diag_(self, norb: int, nelec: int | tuple[int, int]) -> np.ndarray:
    ...

The function ffsim.diag(obj, norb, nelec) dispatches to this method.

The diagonal is useful for constructing preconditioners in iterative eigensolvers like Davidson’s method. The Hamiltonian classes (MolecularHamiltonian, DiagonalCoulombHamiltonian, DoubleFactorizedHamiltonian) implement this protocol.

[5]:

# Get the diagonal of the molecular Hamiltonian
diagonal = ffsim.diag(mol_hamiltonian, norb=norb, nelec=nelec)
print(f"Diagonal shape: {diagonal.shape}")
print(f"Diagonal entries: {diagonal}")

Diagonal shape: (36,)
Diagonal entries: [148.55257853 136.67339982  92.05299049  84.69780382  38.91095487
  47.10339618 136.67339982 154.38312715 105.29562943  80.4305581
  30.17662077  67.95796813  92.05299049 105.29562943 111.99593791
  24.56673251  30.10060138  63.41486035  84.69780382  80.4305581
  24.56673251  98.90176349  41.87149829  57.6758726   38.91095487
  30.17662077  30.10060138  41.87149829  40.62903928  51.9663252
  47.10339618  67.95796813  63.41486035  57.6758726   51.9663252
  92.89251718]

SupportsTrace¶

The SupportsTrace protocol is for linear operators whose trace can be computed. A class implements this protocol by defining the _trace_ method, which has the signature:

def _trace_(self, norb: int, nelec: int | tuple[int, int]) -> float:
    ...

The function ffsim.trace(obj, norb, nelec) dispatches to this method. This function first checks for the _trace_ method and uses that if it exists. If the _trace_ method does not exist, it then checks for a _diag_ method, since the trace can be computed by summing the diagonal. Thus, classes that implement the SupportsDiagonal protocol automatically implement SupportsTrace as well.

One practical use of the trace is to pass it to scipy.sparse.linalg.expm_multiply for Hamiltonian time evolution, which uses the trace to improve its performance. Without it, SciPy issues a warning.

[6]:

# Compute the trace of the molecular Hamiltonian
trace = ffsim.trace(mol_hamiltonian, norb=norb, nelec=nelec)
print(f"Trace: {trace}")

# The trace is just the sum of the diagonal
print(f"Sum of diagonal: {np.sum(diagonal)}")

# Use the trace for efficient time evolution
time = 1.0
linop = ffsim.linear_operator(mol_hamiltonian, norb=norb, nelec=nelec)
evolved_vec = scipy.sparse.linalg.expm_multiply(
    -1j * time * linop,
    ffsim.hartree_fock_state(norb, nelec),
    traceA=-1j * time * trace,
)
evolved_vec

Trace: 2553.1453874075833
Sum of diagonal: 2553.1453874075833

[6]:

array([ 0.19181745-0.13514659j,  0.04791743+0.13143759j,
        0.07135074+0.2036616j , -0.11691946-0.10604848j,
        0.14771837+0.06653658j, -0.10777344+0.03926832j,
        0.04791743+0.13143759j,  0.3815037 -0.16765241j,
        0.12202809+0.02339165j,  0.16799309-0.17841154j,
        0.04098119+0.00198756j,  0.20070653+0.01560637j,
        0.07135074+0.2036616j ,  0.12202809+0.02339165j,
       -0.13877315-0.08873788j,  0.14875463-0.03728076j,
        0.03076148-0.10516136j, -0.00493455+0.03374819j,
       -0.11691946-0.10604848j,  0.16799309-0.17841154j,
        0.14875463-0.03728076j,  0.09707893-0.0119299j ,
        0.00801057-0.04815717j,  0.10232012-0.21822743j,
        0.14771837+0.06653658j,  0.04098119+0.00198756j,
        0.03076148-0.10516136j,  0.00801057-0.04815717j,
       -0.00050482+0.0428504j ,  0.02724199+0.1352515j ,
       -0.10777344+0.03926832j,  0.20070653+0.01560637j,
       -0.00493455+0.03374819j,  0.10232012-0.21822743j,
        0.02724199+0.1352515j ,  0.06175497+0.01483006j])

SupportsApproximateEquality¶

The SupportsApproximateEquality protocol is for objects that can be compared approximately, accounting for floating-point imprecision. This is analogous to using numpy.allclose for arrays. A class implements this protocol by defining the _approx_eq_ method, which has the signature

def _approx_eq_(self, other: Any, rtol: float, atol: float) -> bool:
    ...

The function ffsim.approx_eq(obj, other, rtol=1e-5, atol=1e-8) dispatches to this method. It tries obj._approx_eq_(other) first, then other._approx_eq_(obj), and falls back to obj == other. The tolerance parameters rtol and atol follow the same convention as numpy.allclose.

This protocol is widely implemented across ffsim: Hamiltonians, ansatz operators, and FermionOperator all support it.

[7]:

# Create a second Hamiltonian that differs by a small amount
mol_hamiltonian_2 = ffsim.MolecularHamiltonian(
    one_body_tensor=mol_hamiltonian.one_body_tensor + 1e-9,
    two_body_tensor=mol_hamiltonian.two_body_tensor,
    constant=mol_hamiltonian.constant,
)

# They are approximately equal with the default tolerance
print(ffsim.approx_eq(mol_hamiltonian, mol_hamiltonian_2))

# But not with a very tight tolerance
print(ffsim.approx_eq(mol_hamiltonian, mol_hamiltonian_2, rtol=1e-12, atol=1e-12))

True
False