1. Hamiltonian operators

In this tutorial, we symbolically represent Hamiltonian operators and evaluate commutators between them.

# import the fh_comm package; see the README for installation instructions
import fh_comm as fhc

1.1. Elementary operators

As first step, we define a hopping operator \(h_{ij\sigma}^{} = a_{i\sigma}^\dagger a_{j\sigma}^{} + a_{j\sigma}^\dagger a_{i\sigma}^{}\) between two lattice sites.

# Syntax: HoppingOp(i, j, s, coeff) for lattice coordinates `i` and `j`,
# spin `s = 0` or `s = 1` and a real coefficient `coeff`
h1 = fhc.HoppingOp((2,), (3,), 0, 0.25)
print(h1)

(0.25) h_{(2,), (3,), up}

Let’s define another hopping operator:

h2 = fhc.HoppingOp((3,), (7,), 0, 1.0)
print(h2)

h_{(3,), (7,), up}

We can now symbolically evaluate the commutator of these operators, resulting in an “antisymmetric” hopping operator of the form \(\tilde{h}_{ij\sigma}^{} = a_{i\sigma}^\dagger a_{j\sigma}^{} - a_{j\sigma}^\dagger a_{i\sigma}^{}\):

hc = fhc.commutator(h1, h2)
print(hc)

(0.25) g_{(2,), (7,), up}

Note that hopping operators with different spin or disjoint support (lattice sites they act on) commute. In this case, the returned commutator is the ZeroOp:

# third argument specifies spin-down
h3 = fhc.HoppingOp((3,), (7,), 1, 1.0)
print(fhc.commutator(h1, h3))

<zero op>

Technically, the “support” also includes the spin as last coordinate entry:

h1.support()

[(2, 0), (3, 0)]

The “fermi weight” is the maximum number of fermionic creation and annihilation operators multiplied together in an operator:

h1.fermi_weight

2

Besides the (anti-)symmetric hopping operators, fh_comm also implements number operators \(n_{i\sigma}^{}\):

# Syntax: NumberOp(i, s, coeff) for lattice coordinate `i`,
# spin `s = 0` or `s = 1` and a real coefficient `coeff`
n = fhc.NumberOp((3,), 0, 3.0)
print(n)

(3.0) n_{(3,), up}

print(fhc.commutator(h1, n))

(0.75) g_{(2,), (3,), up}

We can obtain (an upper bound on) the spectral norm of an operator as follows:

h1.norm_bound()

0.25

For elementary and quadratic operators, fh_comm computes the exact spectral norm.

1.2. Products and sums

fh_comm facilitates the representation of products and sums of the elementary operators. A ProductOp collects all the coefficients into a single overal prefactor:

prod = fhc.ProductOp((h1, h3, n), 1.0)
print(prod)

(0.75) (h_{(2,), (3,), up}) @ (h_{(3,), (7,), dn}) @ (n_{(3,), up})

This product operator contains 6 creation and annihilation operators:

prod.fermi_weight

6

Finally, we can define sums of operators via SumOp. Note that the ZeroOp is automatically ignored here.

sumop = fhc.SumOp([prod, fhc.ZeroOp(), h2])
print(sumop)

h_{(3,), (7,), up} + (0.75) (h_{(2,), (3,), up}) @ (h_{(3,), (7,), dn}) @ (n_{(3,), up})

Commutators can be evaluated for products and sums, too:

print(fhc.commutator(prod, sumop))

(0.75) (h_{(2,), (3,), up}) @ (h_{(3,), (7,), dn}) @ (g_{(3,), (7,), up} + (-0.75) (g_{(2,), (3,), up}) @ (h_{(3,), (7,), dn}) @ (n_{(3,), up})) + (0.75) (g_{(2,), (7,), up} + (0.75) (h_{(2,), (3,), up}) @ (h_{(3,), (7,), dn}) @ (g_{(2,), (3,), up})) @ (h_{(3,), (7,), dn}) @ (n_{(3,), up})

1.3. Field operators and matrix representations

fh_comm also implements a lower-level representation (see the class FieldOp) directly based on fermionic creation, annihilation and number operators, mainly used internally for constructing sparse matrix representations.

In the following example, ad_{i, s} stands for \(a_{i\sigma}^\dagger\) and a_{i, s} for \(a_{i\sigma}^{}\).

h1fop = h1.as_field_operator()
print(h1fop)
print(type(h1fop))

(0.25) ad_{(2,), up} a_{(3,), up} + (0.25) ad_{(3,), up} a_{(2,), up}
<class 'fh_comm.field_ops.FieldOp'>

as_compact_matrix() constructs the sparse matrix representation of a fermionic operator based on Jordan-Wigner transformation. The corresponding fermionic modes correspond to the lattice sites which the operator acts on non-trivially.

h1fop.as_compact_matrix()

<4x4 sparse matrix of type '<class 'numpy.float64'>'
        with 2 stored elements in Compressed Sparse Row format>

h1fop.as_compact_matrix().todense()

matrix([[0.  , 0.  , 0.  , 0.  ],
        [0.  , 0.  , 0.25, 0.  ],
        [0.  , 0.25, 0.  , 0.  ],
        [0.  , 0.  , 0.  , 0.  ]])