The Jordan-Wigner and Bravyi-Kitaev Transforms

Ladder operators and the canonical anticommutation relations

A system of $N$ fermionic modes is described by a set of fermionic annihilation operators $\{a_p\}_{p=0}^{N-1}$ satisfying the canonical anticommutation relations $$\begin{aligned} \{a_p, a_q\} &= 0, \label{eq:car1} \\ \{a_p, a^\dagger_q\} &= \delta_{pq}, \label{eq:car2} \end{aligned}$$ where $\{A, B\} := AB + BA$. The adjoint $a^\dagger_p$ of an annihilation operator $a_p$ is called a creation operator, and we refer to creation and annihilation operators as fermionic ladder operators. In a finite-dimensional vector space the anticommutation relations have the following consequences:

• The operators $\{a^\dagger_p a_p\}_{p=0}^{N-1}$ commute with each other and have eigenvalues 0 and 1. These are called the occupation number operators.

• There is a normalized vector $\lvert{\text{vac}}\rangle$, called the vacuum state, which is a mutual 0-eigenvector of all the $a^\dagger_p a_p$.

• If $\lvert{\psi}\rangle$ is a 0-eigenvector of $a_p^\dagger a_p$, then $a_p^\dagger\lvert{\psi}\rangle$ is a 1-eigenvector of $a_p^\dagger a_p$. This explains why we say that $a_p^\dagger$ creates a fermion in mode $p$.

• If $\lvert{\psi}\rangle$ is a 1-eigenvector of $a_p^\dagger a_p$, then $a_p\lvert{\psi}\rangle$ is a 0-eigenvector of $a_p^\dagger a_p$. This explains why we say that $a_p$ annihilates a fermion in mode $p$.

• $a_p^2 = 0$ for all $p$. One cannot create or annihilate a fermion in the same mode twice.

• The set of $2^N$ vectors $$\lvert n_0, \ldots, n_{N-1} \rangle := (a^\dagger_0)^{n_0} \cdots (a^\dagger_{N-1})^{n_{N-1}} \lvert{\text{vac}}\rangle, \qquad n_0, \ldots, n_{N-1} \in \{0, 1\}$$ are orthonormal. We can assume they form a basis for the entire vector space.

• The annihilation operators $a_p$ act on this basis as follows: $$\begin{aligned} a_p \lvert n_0, \ldots, n_{p-1}, 1, n_{p+1}, \ldots, n_{N-1} \rangle &= (-1)^{\sum_{q=0}^{p-1} n_q} \lvert n_0, \ldots, n_{p-1}, 0, n_{p+1}, \ldots, n_{N-1} \rangle \,, \\ a_p \lvert n_0, \ldots, n_{p-1}, 0, n_{p+1}, \ldots, n_{N-1} \rangle &= 0 \,.\end{aligned}$$

See here for a derivation and discussion of these consequences.

Mapping fermions to qubits with transforms

To simulate a system of fermions on a quantum computer, we must choose a representation of the ladder operators on the Hilbert space of the qubits. In other words, we must designate a set of qubit operators (matrices) which satisfy the canonical anticommutation relations. Qubit operators are written in terms of the Pauli matrices $X$, $Y$, and $Z$. In OpenFermion a representation is specified by a transform function which maps fermionic operators (typically instances of FermionOperator) to qubit operators (instances of QubitOperator). In this demo we will discuss the Jordan-Wigner and Bravyi-Kitaev transforms, which are implemented by the functions jordan_wigner and bravyi_kitaev.

The Jordan-Wigner Transform

Under the Jordan-Wigner Transform (JWT), the annihilation operators are mapped to qubit operators as follows: $$\begin{aligned} a_p &\mapsto \frac{1}{2} (X_p + \mathrm{i}Y_p) Z_1 \cdots Z_{p - 1} \\ &= (\lvert{0}\rangle\langle{1}\rvert)_p Z_1 \cdots Z_{p - 1} \\ &=: \tilde{a}_p. \end{aligned}$$

This operator has the following action on a computational basis vector $\lvert z_0, \ldots, z_{N-1} \rangle$: $$\begin{aligned} \tilde{a}_p \lvert z_0 \ldots, z_{p-1}, 1, z_{p+1}, \ldots, z_{N-1} \rangle &= (-1)^{\sum_{q=0}^{p-1} z_q} \lvert z_0 \ldots, z_{p-1}, 0, z_{p+1}, \ldots, z_{N-1} \rangle \\ \tilde{a}_p \lvert z_0 \ldots, z_{p-1}, 0, z_{p+1}, \ldots, z_{N-1} \rangle &= 0. \end{aligned}$$

Note that $\lvert n_0, \ldots, n_{N-1} \rangle$ is a basis vector in the Hilbert space of fermions, while $\lvert z_0, \ldots, z_{N-1} \rangle$ is a basis vector in the Hilbert space of qubits. Similarly, in OpenFermion $a_p$ is a FermionOperator while $\tilde{a}_p$ is a QubitOperator.

Let's instantiate some FermionOperators, map them to QubitOperators using the JWT, and check that the resulting operators satisfy the expected relations.

The parity transform

By comparing the action of $\tilde{a}_p$ on $\lvert z_0, \ldots, z_{N-1} \rangle$ in the JWT with the action of $a_p$ on $\lvert n_0, \ldots, n_{N-1} \rangle$ (described in the first section of this demo), we can see that the JWT is associated with a particular mapping of bitstrings $e: \{0, 1\}^N \to \{0, 1\}^N$, namely, the identity map $e(x) = x$. In other words, under the JWT, the fermionic basis vector $\lvert n_0, \ldots, n_{N-1} \rangle$ is represented by the computational basis vector $\lvert z_0, \ldots, z_{N-1} \rangle$, where $z_p = n_p$ for all $p$. We can write this as $$\lvert x \rangle \mapsto \lvert e(x) \rangle,$$ where the vector on the left is fermionic and the vector on the right is qubit. We call the mapping $e$ an encoder.

There are other transforms which are associated with different encoders. To see why we might be interested in these other transforms, observe that under the JWT, $\tilde{a}_p$ acts not only on qubit $p$ but also on qubits $0, \ldots, p-1$. This means that fermionic operators with low weight can get mapped to qubit operators with high weight, where by weight we mean the number of modes or qubits an operators acts on. There are some disadvantages to having high-weight operators; for instance, they may require more gates to simulate and are more expensive to measure on some near-term hardware platforms. In the worst case, the annihilation operator on the last mode will map to an operator which acts on all the qubits. To emphasize this point let's apply the JWT to the annihilation operator on mode 99:

The purpose of the string of Pauli $Z$'s is to introduce the phase factor $(-1)^{\sum_{q=0}^{p-1} n_q}$ when acting on a computational basis state; when $e$ is the identity encoder, the modulo-2 sum $\sum_{q=0}^{p-1} n_q$ is computed as $\sum_{q=0}^{p-1} z_q$, which requires reading $p$ bits and leads to a Pauli $Z$ string with weight $p$. A simple solution to this problem is to consider instead the encoder defined by $$e(x)_p = \sum_{q=0}^p x_q \quad (\text{mod 2}),$$ which is associated with the mapping of basis vectors $\lvert n_0, \ldots, n_{N-1} \rangle \mapsto \lvert z_0, \ldots, z_{N-1} \rangle,$ where $z_p = \sum_{q=0}^p n_q$ (again addition is modulo 2). With this encoding, we can compute the sum $\sum_{q=0}^{p-1} n_q$ by reading just one bit because this is the value stored by $z_{p-1}$. The associated transform is called the parity transform because the $p$-th qubit is storing the parity (modulo-2 sum) of modes $0, \ldots, p$. Under the parity transform, annihilation operators are mapped as follows: $$\begin{aligned} a_p &\mapsto \frac{1}{2} (X_p Z_{p - 1} + \mathrm{i}Y_p) X_{p + 1} \cdots X_{N} \\ &= \frac{1}{4} [(X_p + \mathrm{i} Y_p) (I + Z_{p - 1}) - (X_p - \mathrm{i} Y_p) (I - Z_{p - 1})] X_{p + 1} \cdots X_{N} \\ &= [(\lvert{0}\rangle\langle{1}\rvert)_p (\lvert{0}\rangle\langle{0}\rvert)_{p - 1} - (\lvert{0}\rangle\langle{1}\rvert)_p (\lvert{1}\rangle\langle{1}\rvert)_{p - 1}] X_{p + 1} \cdots X_{N} \\ \end{aligned}$$

The term in brackets in the last line means "if $z_p = n_p$ then annihilate in mode $p$; otherwise, create in mode $p$ and attach a minus sign". The value stored by $z_{p-1}$ contains the information needed to determine whether a minus sign should be attached or not. However, now there is a string of Pauli $X$'s acting on modes $p+1, \ldots, N-1$ and hence using the parity transform also yields operators with high weight. These Pauli $X$'s perform the necessary update to $z_{p+1}, \ldots, z_{N-1}$ which is needed if the value of $n_{p}$ changes. In the worst case, the annihilation operator on the first mode will map to an operator which acts on all the qubits.

Since the parity transform does not offer any advantages over the JWT, OpenFermion does not include a standalone function to perform it. However, there is functionality for defining new transforms by specifying an encoder and decoder pair, also known as a binary code (in our examples the decoder is simply the inverse mapping), and the binary code which defines the parity transform is included in the library as an example. See examples/binary_code_transforms_demo.ipynb for a demonstration of this functionality and how it can be used to reduce the qubit resources required for certain applications.

Let's use this functionality to map our previously instantiated FermionOperators to QubitOperators using the parity transform with 10 total modes and check that the resulting operators satisfy the expected relations.

Now let's map one of the FermionOperators again but with the total number of modes set to 100.

Note that with the JWT, it is not necessary to specify the total number of modes in the system because $\tilde{a}_p$ only acts on qubits $0, \ldots, p$ and not any higher ones.

The Bravyi-Kitaev transform

The discussion above suggests that we can think of the action of a transformed annihilation operator $\tilde{a}_p$ on a computational basis vector $\lvert z \rangle$ as a 4-step classical algorithm:

1. Check if $n_p = 0$. If so, then output the zero vector. Otherwise,
2. Update the bit stored by $z_p$.
3. Update the rest of the bits $z_q$, $q \neq p$.
4. Multiply by the parity $\sum_{q=0}^{p-1} n_p$.

Under the JWT, Steps 1, 2, and 3 are represented by the operator $(\lvert{0}\rangle\langle{1}\rvert)_p$ and Step 4 is accomplished by the operator $Z_{0} \cdots Z_{p-1}$ (Step 3 actually requires no action). Under the parity transform, Steps 1, 2, and 4 are represented by the operator $(\lvert{0}\rangle\langle{1}\rvert)_p (\lvert{0}\rangle\langle{0}\rvert)_{p - 1} - (\lvert{0}\rangle\langle{1}\rvert)_p (\lvert{1}\rangle\langle{1}\rvert)_{p - 1}$ and Step 3 is accomplished by the operator $X_{p+1} \cdots X_{N-1}$.

To obtain a simpler description of these and other transforms (with an aim at generalizing), it is better to put aside the ladder operators and work with an alternative set of $2N$ operators defined by $$c_p = a_p + a_p^\dagger\,, \qquad d_p = -\mathrm{i} (a_p - a_p^\dagger)\,.$$ These operators are known as Majorana operators. Note that if we describe how Majorana operators should be transformed, then we also know how the annihilation operators should be transformed, since $$a_p = \frac{1}{2} (c_p + \mathrm{i} d_p).$$

For simplicity, let's consider just the $c_p$; the $d_p$ are treated similarly. The action of $c_p$ on a fermionic basis vector is given by $$c_p \lvert n_0, \ldots, n_{p-1}, n_p, n_{p+1}, \ldots, n_{N-1} \rangle = (-1)^{\sum_{q=0}^{p-1} n_q} \lvert n_0, \ldots, n_{p-1}, 1 - n_p, n_{p+1}, \ldots, n_{N-1} \rangle$$

In words, $c_p$ flips the occupation of mode $p$ and multiplies by the ever-present parity factor. If we transform $c_p$ to a qubit operator $\tilde{c}_p$, we should be able to describe the action of $\tilde{c}_p$ on a computational basis vector $\lvert z \rangle$ with a 2-step classical algorithm:

1. Update the string $z$ to a new string $z'$.
2. Multiply by the parity $\sum_{q=0}^{p-1} n_q$.

Step 1 amounts to flipping some bits, so it will be performed by some Pauli $X$'s, and Step 2 will be performed by some Pauli $Z$'s. So $\tilde{c}_p$ should take the form $$\tilde{c}_p = X_{U(p)} Z_{P(p - 1)},$$ where $U(j)$ is the set of bits that need to be updated upon flipping $n_j$, and $P(j)$ is a set of bits that stores the sum $\sum_{q=0}^{j} n_q$ (let's define $P(-1)$ to be the empty set). Let's see how this looks under the JWT and parity transforms.

For the JWT, $U(j) = \{j\}$ and $P(j) = \{0, \ldots, j\}$, whereas for the parity transform, $U(j) = \{j, \ldots, N-1\}$ and $P(j) = \{j\}$. The size of these sets can be as large as $N$, the total number of modes. These sets are determined by the encoding function $e$.

It is possible to pick a clever encoder with the property that these sets have size $O(\log N)$. The corresponding transform will map annihilation operators to qubit operators with weight $O(\log N)$, which is much smaller than the $\Omega(N)$ weight associated with the JWT and parity transforms. This fact was noticed by Bravyi and Kitaev, and later Havlíček and others pointed out that the encoder which achieves this is implemented by a classical data structure called a Fenwick tree. The transforms described in these two papers actually correspond to different variants of the Fenwick tree data structure and give different results when the total number of modes is not a power of 2. OpenFermion implements the one from the first paper as bravyi_kitaev and the one from the second paper as bravyi_kitaev_tree. Generally, the first one (bravyi_kitaev) is preferred because it results in operators with lower weight and is faster to compute.

Let's transform our previously instantiated Majorana operator using the Bravyi-Kitaev transform.

The advantage of the Bravyi-Kitaev transform is not apparent in a system with so few modes. Let's look at a system with 100 modes.

Now let's go back to a system with 10 modes and check that the Bravyi-Kitaev transformed operators satisfy the expected relations.