Gaussian Hamiltonians¶
Section author: Josh Izaac <josh@xanadu.ai>
Note
This section assumes a familiarity with the operators and conventions used in quantum optics and continuous-variable (CV) quantum computation. For an introduction to these concepts, see the Strawberry Fields documentation.
Quadrature operators¶
A Gaussian Hamiltonian is a Hamiltonian that can be written as a quadratic polynomial of the quadrature operators \(\q_i\) (canonical position, also written as \(\x_i\)) and \(\p_j\) (canonical momentum) [2]:
where:
- \(\r=(\q_1,\dots,\q_{N},\p_1,\dots,\p_N)\) is a vector of the canonical position and momentum operators acting on qumode \(i\), in \(\q\p\)-ordering,
- \(A\in\mathbb{R}^{2N\times 2N}\) is a symmetric matrix containing the quadratic coefficients of terms of the form \(\hat{x}_i\hat{y}_j\),
- \(\mathbf{d}\in\mathbb{R}^{2N}\) is a vector containing the linear coefficients corresponding to terms of the form \(\hat{x}_i\).
For example, the following are all quadratic Hamiltonians:
- \(\hat{H} = \q_0^2 -q_0\)
- \(\hat{H} = \q_0 \q_1 - \p_0\p_1\)
while \(\hat{H}=\q_0^2\p_1\), with a third-order term in the quadrature operators, is not.
Ladder operators¶
In terms of the bosonic annihilation and creation operators
this corresponds to a Hamiltonian of the form [2]
where \(\av = (\a_1, \a_2,\dots,\a_N)\), \(\mathbf{w}\in\mathbb{C}^N\) and \(F,G\in\mathbb{C}^{N\times N}\). This corresponds to a linear unitary Bogoliubov transform of the annihilation and creation operators:
where \(AB^T=BA^T\) and \(AA^\dagger = BB^\dagger+\I\).
Time propagation¶
Once the matrix of quadratic coefficients \(A\) and vector of linear coefficients \(\mathbf{d}\) have been found for the Gaussian Hamiltonian, it can be shown that the Hamiltonian can always be written in the following form, up to a local phase factor [3]:
where \(\mathbf{d}'=-A^{-1}\mathbf{d}\). Let’s consider the case of zero linear coefficients and non-zero linear coefficients separately.
Zero linear coefficients¶
In the Heisenberg picture, with \(\hat{H}=\frac{1}{2}\r^T A\r\), the time-evolution of the quadrature operators must satisfy the Heisenberg equations of motion:
where
is the symplectic matrix, coming from the definition of the canonical commutation relation \([\r_i,\r_j]=i\hbar \Omega_{ij}\).
Solving this differential equation, we find that the symplectic Gaussian transformation describing the time-evolution of the Hamiltonian \(\hat{H}\) acting on the quadrature operators, for time \(t\), is given by:
Non-zero linear coefficients¶
If we have non-zero linear coefficients, we need to take this into account by performing the required displacement operation. Consider the following Gaussian Hamiltonian:
where \(\mathbf{d}'=-A^{-1}\mathbf{d}\), as before. This corresponds to the action of a Weyl operator or displacement of the form \(\hat{D}(\mathbf{s})\) with \(\mathbf{s}=-\mathbf{d}'/\sqrt{2\hbar}\):
Calculating the time-evolution operator,
In order to write this as a symplectic matrix transformation, we need to move all displacement operators to the left. To do this, we can post-multiply by \(\I=e^{i\hat{H}(0)t}e^{-i\hat{H}(0)t}\):
Finally, we can rewrite this as a symplectic transformation, by making the substitution \(e^{-i\hat{H}(0)t}\rightarrow e^{\Omega A t}\) and by noting that the bracketed term is simply a displacement by \(-\mathbf{s}\), evolved under \(\hat{H}(0)\) for time \(t\):
Definition
For a quadratic Hamiltonian of the form \(\hat{H} = \frac{1}{2}\r^T A\r + \r^T \mathbf{d}\), the symplectic transformation \(S\in\mathbb{R}^{2N\times 2N}\) characterizing the time-evolution unitary operator \(\hat{U}(t) = e^{-i\hat{H}t/\hbar}\) is given by
where \(\Omega\) is the symplectic matrix, \(\hat{D}\) the displacement operation, and \(\mathbf{s} = -A^{-1}\mathbf{d}/\sqrt{2\hbar}\).
Tip
Implemented in SFOpenBoson as a quantum operation by sfopenboson.ops.GaussianPropagation
Warning
In the case where the quadratic coefficient matrix \(A\) is singular, for example the Hamltonian \(\hat{H}=\frac{1}{2}p_0^2+q_0\), in order to calculate \(A^{-1}\) to determine the resulting displacement, the matrix \(A\) is perturbed by \(\epsilon\ll 1\):