Bose-Hubbard Hamiltonians

Section author: Josh Izaac <josh@xanadu.ai>

The Bose-Hubbard model describes the dynamics of multiple bosons on a lattice composed of orthogonal positions states — referred to as vertices or nodes — with the overall energy of the system dependent on the location of the bosons, surrounding potentials, and possible interactions.

Bose-Hubbard Hamiltonians have the following form:

\[\hat{H} = -J\sum_{i=1}^N\sum_{j=1}^N A_{ij} \ad_i\a_j + \frac{1}{2}U\sum_{i=1}^N \ad_i \a_i (\ad_i \a_i - 1) - \mu \sum_{i=1}^N \ad_i \a_i + V \sum_{i=1}^N\sum_{j=1}^N A_{ij} \ad_i \a_i \ad_j \a_j.\]

where:

  • \(A\) is a real symmetric matrix of ones and zeros defining the adjacency of each pairwise combination of nodes \((i,j)\) in the \(N\)-node system.
  • \(J\) represents the transfer integral or hopping term of the boson between nodes.
  • \(U\) is the on-site interaction potential, the interaction strength between bosons on the same node. This can be attractive (\(U<0\)) or repulsive (\(U>0\)).
  • \(\mu\) is the chemical potential, dependent on the number of bosons per node.
  • \(V\) is the dipole-dipole or nearest-neighbour interaction term, the interaction strength between bosons on adjacent nodes. Like the on-site interaction, this can be attractive or repulsive.

CV decomposition

As the Bose-Hubbard Hamiltonian is time-independent, it suffices to find a continuous-variable (CV) gate decomposition for the unitary operator \(\hat{U}=e^{-i\hat{H}t}\) [1]. We can rewrite the Hamiltonian in the following form:

\[\hat{H} = -J\sum_{i=1}^N\sum_{j=1}^N A_{ij} \ad_i\a_j + \sum_{\ell=1}^N \left(\frac{1}{2}U \hat{n}_\ell^2 - \left(\frac{1}{2}U+\mu\right) \hat{n}_\ell\right) + V \sum_{i=1}^N\sum_{j=1}^N \hat{n}_i\hat{n}_j,\]

where \(\hat{n}_i=\ad_i\a_i\) is the bosonic number operator. Taking the matrix exponential, and applying the Lie product formula,

\[\begin{split}e^{-iHt} = \lim_{k\rightarrow\infty}\left[\prod_{\substack{i,j\\i\sim j}}\exp\left({i\frac{ J t}{k}(\ad_i\a_j + \ad_j\a_i)}\right)\exp\left({-i\frac{V t}{k}\hat{n}_i\hat{n}_j}\right)\prod_{\ell}\exp\left(-i\frac{Ut}{2k}\hat{n}_\ell^2\right)\exp\left(i\frac{(U+2\mu)t}{2k}\hat{n}_\ell\right)\right]^k,\end{split}\]

where \(i\sim j\) indicates we are only summing over adjacent nodes (i.e., those where \(A_{ij}=1\)). Truncating \(k\) to a reasonable value, we end up with the approximation

\[\begin{split}e^{-iHt} = \left[\prod_{\substack{i,j\\i\sim j}}\exp\left({i\frac{ J t}{k}(\ad_i\a_j + \ad_j\a_i)}\right)\exp\left({-i\frac{V t}{k}\hat{n}_i\hat{n}_j}\right)\prod_{\ell}\exp\left(-i\frac{Ut}{2k}\hat{n}_\ell^2\right)\exp\left(i\frac{(U+2\mu)t}{2k}\hat{n}_\ell\right)\right]^k + \mathcal{O}(t^2/k).\end{split}\]

Comparing these individual bracketed operators with the known CV gate set (see here for more details) we see that they correspond to a beamsplitter, cross-Kerr gate, Kerr gate, and phase-space rotation respectively:

  • \(\exp\left({i\frac{ J t}{k}(\ad_i\a_j + \ad_j\a_i)}\right)\equiv BS(\theta, \phi)\) where \(\theta=Jt/k\), \(\phi=\pi/2\),
  • \(\exp\left({-i\frac{V t}{k}\hat{n}_i\hat{n}_j}\right)\equiv CK(\kappa)\) where \(\kappa=-Vt/k\),
  • \(\exp\left(-i\frac{Ut}{2k}\hat{n}_\ell^2\right)\equiv K(\kappa)\) where \(\kappa=-Ut/2k\),
  • \(\exp\left(i\frac{(U+2\mu)t}{2k}\hat{n}_\ell\right)\equiv R(r)\) where \(r=(U+2\mu)t/2k\).

Decomposition

A Bose-Hubbard Hamiltonian can be implemented to arbitrary error via a decomposition of beamsplitters, cross-Kerr interactions, Kerr gates, and phase-space rotations.

Tip

Implemented in SFOpenBoson as a quantum operation by sfopenboson.ops.BoseHubbardPropagation