A Method to Compute the Schrieffer–Wolff Generator for Analysis of Quantum Memory

Abstract

Quantum illumination uses entangled light that consists of signal and idler modes to achieve higher detection rate of a low-reflective object in noisy environments. The best performance of quantum illumination can be achieved by measuring the returned signal mode together with the idler mode. Thus, it is necessary to prepare a quantum memory that can keep the idler mode ideal. To send a signal towards a long-distance target, entangled light in the microwave regime is used. There was a recent demonstration of a microwave quantum memory using microwave cavities coupled with a transmon qubit. We propose an ordering of bosonic operators to efficiently compute the Schrieffer–Wolff transformation generator to analyze the quantum memory. Our proposed method is applicable to a wide class of systems described by bosonic operators whose interaction part represents a definite number of transfer in quanta.

1. Introduction

Quantum memories are required to store and retrieve quantum states with high fidelity. To synchronize various events, quantum memories are essential for quantum information networks, including quantum computation [1], quantum communication [2], and quantum illumination [3]. Quantum illumination (QI), a target detection scheme using quantum entangled light with signal and idler modes, has as its objective enhancing the detection rate of a target with low-reflectivity in a highly noisy environment [3]. In QI, the signal mode is sent to the target while the idler mode is retained. Although the noisy environment destroys the entanglement between the signal and idler modes, we can take quantum advantage over the classical limit by jointly measuring the returned signal mode and the idler mode when the signal arrives [4,5,6,7]. During this process, it is highly appreciable to keep the idler mode in an ideal quantum memory. This was investigated using various systems, such as a microwave cavity [8], mechanical oscillators [9], or spin ensembles [10,11].

Here, we focus on quantum memories using microwave cavities that can have high-quality factors and allow continuous-variable quantum information processes. By coupling a microwave cavity to a transmon qubit, it is able to write arbitrary states on the cavity and infer information about the cavity [8,12]. It is based on the cross Kerr effect, where the energy gap of neighboring levels of the cavity (transmon qubit) depends on the excitations of the transmon qubit (cavity). The anharmonicity of the transmon qubit gives rise to the cross Kerr effect through coupling of the qubit and cavity [13,14].

In dealing with such systems, it is crucial to understand how the coupling affects the energy structure. The Schrieffer–Wolff transformation computes this shift in energy structure by using a basis change unitary to remove the coupling [15,16]. For multiple bosonic modes containing nonlinear terms, it is complicated to find the exact generator of the unitary for the Schrieffer–Wolff transformation. Here, we propose a systematic approach to find the generator and compute the energy corrections induced by this transformation. An ordering of operators, which we call computational ordering, greatly simplifies the commutation structure of operators making it suitable in calculating the Schrieffer–Wolff transformation generator.

2. Equivalent Circuit Analysis of a Quantum Memory

The quantum memory demonstrated in Ref. [8] couples two microwave cavities through a transmon qubit. One cavity is used as a memory (storage) to store quantum states and the other cavity is used as a readout port whose response depends on the memory-cavity state through the transmon qubit. Such a system can be described by an equivalent circuit depicted in Figure 1b. Two LC oscillators represent the microwave cavities, while the middle oscillator represents the transmon qubit. The oscillators are labeled as s, t, and r for storage, transmon, and readout, respectively, as in Ref. [8]. The transmon qubit is coupled to both cavities by capacitors. The Hamiltonian describing this system is

$$\frac{\widehat{\mathcal{H}}}{\hslash }=\sum _{i=s,t,r}{\omega }_i{\widehat{a}}_i^{†}{\widehat{a}}_i-\frac{E_C}{2\hslash }{\widehat{a}}_t^{†}{\widehat{a}}_t^{†}{\widehat{a}}_t{\widehat{a}}_t+g_1({\widehat{a}}_s^{†}{\widehat{a}}_t+{\widehat{a}}_s{\widehat{a}}_t^{†})+g_2({\widehat{a}}_t^{†}{\widehat{a}}_r+{\widehat{a}}_t{\widehat{a}}_r^{†}),$$

(1)

where ${\widehat{a}}_s$, ${\widehat{a}}_t$, ${\widehat{a}}_r$ are bosonic annihilation operators corresponding to each oscillator mode. A detailed derivation of this Hamiltonian and expressions of ${\omega }_s,{\omega }_t,\dots$ in terms of $L_i$, $C_i$, $C_{c1}$, $C_{c2}$$(i=s,t,r)$ are given in Appendix A. We assume that the system is in the dispersive regime, where the couplings $g_1$ and $g_2$ are much smaller than the detunings $|{\omega }_s-{\omega }_t|$ and $|{\omega }_t-{\omega }_r|$.

The elimination of capacitive couplings in Equation (1) gives rise to cross Kerr effects among each cavity and the transmon. This is done by diagonalizing the Hamiltonian. There are two ways in achieving this; namely, second-order perturbation and the Schrieffer–Wolff transformation [16]. These methods were previously applied to systems of transmon qubits coupled with LC oscillators in evaluating the energy structure [13,14]. For the Schrieffer–Wolff transformation, one must find an operator $\widehat{S}$, the Schrieffer–Wolff generator, which is an off-diagonal operator satisfying a given commutator equation. It is complicated to determine the operator $\widehat{S}$ and compute various commutators to obtain the second-order energy corrections. Thus, we introduce a method that simplifies the computation and apply it to Equation (1). After the computation, we truncate the transmon qubit to the lowest two levels to obtain a Jaynes–Cummings-like Hamiltonian [13]. Truncation of the transmon qubit should be done after diagonalization since virtual excitations of the transmon need to be considered.

2.1. Computational Ordering for Schrieffer–Wolff Transformation

We propose an ordering of bosonic operators, which gives direct computation of the Schrieffer–Wolff generator and second-order energy corrections. For a short recall of the Schrieffer–Wolff transformation, let $\widehat{\mathcal{H}}$ be the Hamiltonian in interest. We separate the Hamiltonian into diagonal and off-diagonal parts, $\widehat{H_0}$, $\widehat{V}$ respectively, so $\widehat{\mathcal{H}}={\widehat{H}}_0+\widehat{V}$. The Schrieffer–Wolff generator $\widehat{S}$ is defined as the off-diagonal operator which satisfies $[\widehat{S},\widehat{H_0}]=-\widehat{V}$. Then,

$$e^{\widehat{S}}\widehat{\mathcal{H}}{e}^{-\widehat{S}}={\widehat{H}}_0+\frac{1}{2}[\widehat{S},\widehat{V}]+\dots .$$

(2)

The energy separations of ${\widehat{H}}_0$ must be larger than $\widehat{V}$, so that $\widehat{S}$ becomes small and a perturbative approach is applicable [17]. This condition becomes evident when we write down the explicit form of $\widehat{S}$ in Equation (5). The second-order energy corrections to ${\widehat{H}}_0$ are given as the diagonal part of $\frac{1}{2}[\widehat{S},\widehat{V}]$.

We consider a system of N bosonic modes, where ${\widehat{a}}_i$ is the annihilation operator of the i-th mode and satisfies $[{\widehat{a}}_i,{\widehat{a}}_j^{†}]={\delta }_{ij}$. We propose an ordering of operators

$${\widehat{a}}^{†\mathbf{n}}f({\widehat{a}}^{†}\widehat{a}){\widehat{a}}^{\mathbf{m}}$$

(3)

to efficiently compute the Schrieffer–Wolff transformation generator and second-order energy corrections to ${\widehat{H}}_0$. In Equation (3), ${\widehat{a}}^{†}\widehat{a}=({\widehat{a}}_1^{†}{\widehat{a}}_1,\dots ,{\widehat{a}}_N^{†}{\widehat{a}}_N)$, $\mathbf{n}=(n_1,\dots ,n_N)$, $\mathbf{m}=(m_1,\dots ,m_N)$ are N tuples and ${\widehat{a}}^{†\mathbf{n}}={\widehat{a}}_1^{†n_1}\dots {\widehat{a}}_N^{†n_N}$, ${\widehat{a}}^{\mathbf{m}}={\widehat{a}}_1^{m_1}\dots {\widehat{a}}_N^{m_N}$. To ensure that f is unique, we require n, m to have disjoint support, i.e. $\mathbf{n}·\mathbf{m}:=(n_1{m}_1,\dots ,n_N{m}_N)=(0,\dots ,0)$. For example, the operator ${\widehat{a}}_1^{†}{\widehat{a}}_1^{†}{\widehat{a}}_1^{†}{\widehat{a}}_1{\widehat{a}}_1$ will be written as

$$\begin{array}{cc}\hfill {\widehat{a}}_1^{†}{\widehat{a}}_1^{†}{\widehat{a}}_1^{†}{\widehat{a}}_1{\widehat{a}}_1& ={\widehat{a}}_1^{†}({\widehat{a}}_1^{†}{\widehat{a}}_1{\widehat{a}}_1^{†}{\widehat{a}}_1-{\widehat{a}}_1^{†}{\widehat{a}}_1)={\widehat{a}}_1^{†}f({\widehat{a}}^{†}\widehat{a}),\hfill \\ \hfill f(x_1,\dots ,x_N)& =x_1^2-x_1.\hfill \end{array}$$

The main motivation of this ordering is that diagonal operators in the Fock basis correspond to functions defined on ${\mathbb{N}}_0^N$, with ${\mathbb{N}}_0=\{0,1,2,3,\dots \}$, and functions are in general easier to manipulate than operators. The computational ordering is then equivalent to writing a given operator in terms of number operators as much as possible. Explicit expressions for writing normal-ordered or antinormal-ordered operators in this ordering are given in Appendix B. Note that operators that have $\mathbf{n}=\mathbf{m}=0$ are exactly the diagonal operators in the Fock basis.

We write the Hamiltonian $\widehat{\mathcal{H}}$ in this ordering:

$$\widehat{\mathcal{H}}=f({\widehat{a}}^{†}\widehat{a})+{\sum }_{\mathbf{n},\mathbf{m}}{\widehat{a}}^{†\mathrm{n}}{g}_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}\widehat{a}){\widehat{a}}^{\mathrm{m}}.$$

(4)

The sum, here and henceforth, is over all n, m satisfying $\mathbf{n}·\mathbf{m}=(0,\dots ,0)$ and n, m are not both 0. This automatically splits the Hamiltonian into diagonal and off-diagonal parts. The hermitian condition on $\widehat{\mathcal{H}}$ forces f to be real valued and $g_{\mathbf{n}\mathbf{m}}^{*}=g_{\mathbf{m}\mathbf{n}}$, $z^{*}$ being the complex conjugate of z. The main results are

$$\begin{array}{cc}\hfill \widehat{S}& =\mathbf{\sum }_{\mathbf{n},\mathbf{m}}{\widehat{a}}^{†\mathbf{n}}\frac{g_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}\widehat{a})}{f({\widehat{a}}^{†}\widehat{a}+\mathbf{n})-f({\widehat{a}}^{†}\widehat{a}+\mathbf{m})}{\widehat{a}}^{\mathbf{m}},\hfill \end{array}$$

(5)

$$\begin{array}{cc}\hfill \frac{1}{2}{[\widehat{S},\widehat{V}]}_{(d)}& =\mathbf{\sum }_{\mathbf{n},\mathbf{m}}\frac{{({\widehat{a}}^{†}\widehat{a})}^{\underline{\mathbf{n}}}{({\widehat{a}}^{†}\widehat{a}-\mathbf{n}+\mathbf{m})}^{\underline{\mathbf{m}}}{|g_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}a-\mathbf{n})|}^2}{f({\widehat{a}}^{†}\widehat{a})-f({\widehat{a}}^{†}\widehat{a}-\mathbf{n}+\mathbf{m})}.\hfill \end{array}$$

(6)

where $x^{\underline{n}}=x(x-1)\dots (x-n+1)$ is the falling factorial and the falling factorial of tuples is defined element-wise. The subscript ${}_{(d)}$ means to take the diagonal part, so $\frac{1}{2}{[\widehat{S},\widehat{V}]}_{(d)}$ is the second-order correction to energy. Since $f,g_{\mathbf{n}\mathbf{m}}$ are essentially functions defined on ${\mathbb{N}}_0^N$ as noted before, the computation of Equations (5) and (6) is straightforward. In the end, the original Hamiltonian is transformed via the Schrieffer–Wolff transformation as

$$e^{\widehat{S}}\widehat{\mathcal{H}}{e}^{-\widehat{S}}={\widehat{H}}_0+\frac{1}{2}{[\widehat{S},\widehat{V}]}_{(d)}+\cdots ={\widehat{H}}_0+{\widehat{H}}^{(2)}+\dots$$

(7)

where the omitted terms are off-diagonal terms of second-order in $\widehat{V}$ or diagonal terms of third-order in $\widehat{V}$. The superscript ${}^{(2)}$ indicates that the term is second-order in $\widehat{V}$.

To obtain the main results Equations (5) and (6), we need a computational lemma.

Lemma 1.

The commutators of ${\widehat{a}}^{†\mathbf{n}}$, ${\widehat{a}}^{\mathbf{m}}$ with $f({\widehat{a}}^{†}\widehat{a})$ are as follows.

$$\begin{array}{cc}\hfill [{\widehat{a}}^{†\mathbf{n}},f({\widehat{a}}^{†}\widehat{a})]& ={\widehat{a}}^{†\mathbf{n}}(f({\widehat{a}}^{†}\widehat{a})-f({\widehat{a}}^{†}\widehat{a}+\mathbf{n})),\hfill \end{array}$$

(8)

$$\begin{array}{cc}\hfill {\widehat{a}}^{\mathbf{m}},f({\widehat{a}}^{†}\widehat{a})]& =(f({\widehat{a}}^{†}\widehat{a}+\mathbf{m})-f({\widehat{a}}^{†}\widehat{a})){\widehat{a}}^{\mathbf{m}}.\hfill \end{array}$$

(9)

Proof. 

It suffices to check on number states $|\mathbf{k}\rangle =|k_1,\dots ,k_N\rangle$. One can verify

$$\begin{array}{cc}\hfill [{\widehat{a}}^{†\mathbf{n}},f({\widehat{a}}^{†}\widehat{a})]|\mathbf{k}\rangle & ={\widehat{a}}^{†\mathbf{n}}f({\widehat{a}}^{†}\widehat{a})|\mathbf{k}\rangle -f({\widehat{a}}^{†}\widehat{a}){\widehat{a}}^{†\mathbf{n}}|\mathbf{k}\rangle \hfill \\ \hfill & =f(\mathbf{k}){\widehat{a}}^{†\mathbf{n}}|\mathbf{k}\rangle -{\left(\prod _{i=1}^N{(k_i+1)}^{\overline{{\mathbf{n}}_i}}\right)}^{1/2}f({\widehat{a}}^{†}\widehat{a})|\mathbf{k}+\mathbf{n}\rangle \hfill \\ \hfill & ={\left(\prod _{i=1}^N{(k_i+1)}^{\overline{n_i}}\right)}^{1/2}(f(\mathbf{k})-f(\mathbf{k}+\mathbf{n}))|\mathbf{k}+\mathbf{n}\rangle \hfill \\ \hfill & ={\widehat{a}}^{†\mathbf{n}}(f(\mathbf{k})-f(\mathbf{k}+\mathbf{n}))|\mathbf{k}\rangle \hfill \\ \hfill & ={\widehat{a}}^{†\mathbf{n}}(f({\widehat{a}}^{†}\widehat{a})-f({\widehat{a}}^{†}\widehat{a}+\mathbf{n}))|\mathbf{k}\rangle .\hfill \end{array}$$

Here, $x^{\overline{n}}=x(x+1)\dots (x+n-1)$ is the rising factorial. The commutator with ${\widehat{a}}^{\mathbf{m}}$ follows from taking the adjoint. □

Now one can compute the commutator of $\widehat{S}$ with ${\widehat{H}}_0$. Write $\widehat{S}$ as

$$\widehat{S}=\sum _{\mathbf{n},\mathbf{m}}{\widehat{a}}^{†\mathbf{n}}{h}_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}\widehat{a}){\widehat{a}}^{\mathbf{m}},$$

(10)

with $h_{\mathbf{n}\mathbf{m}}^{*}=-h_{\mathbf{m}\mathbf{n}}$ so that $\widehat{S}$ is antihermitian. Then, one has

$$\begin{array}{cc}\hfill [\widehat{S},{\widehat{H}}_0]& =\sum _{\mathbf{n},\mathbf{m}}[{\widehat{a}}^{†\mathbf{n}}{h}_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}\widehat{a}){\widehat{a}}^{\mathbf{m}},{\widehat{H}}_0]\hfill \\ \hfill & =\sum _{\mathbf{n},\mathbf{m}}{\widehat{a}}^{†\mathbf{n}}{h}_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}\widehat{a})(f({\widehat{a}}^{†}\widehat{a}+\mathbf{m})-f({\widehat{a}}^{†}\widehat{a}+\mathbf{n})){\widehat{a}}^{\mathbf{m}}.\hfill \end{array}$$

The choice of $\widehat{S}$ as in Equation (5) yields $[\widehat{S},\widehat{H_0}]=-\widehat{V}$, i.e., we take $h_{\mathbf{n}\mathbf{m}}$ as

$$h_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}\widehat{a})=\frac{g_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}\widehat{a})}{f({\widehat{a}}^{†}\widehat{a}+\mathbf{n})-f({\widehat{a}}^{†}\widehat{a}+\mathbf{m})}.$$

(11)

The conditions on f, $g_{\mathbf{n}\mathbf{m}}$ ensure that $h_{\mathbf{n}\mathbf{m}}^{*}=-h_{\mathbf{m}\mathbf{n}}$ holds. This is well-defined as long as the diagonal part is nondegenerate, which is true when considering low excitations of transmons.

Using the generator $\widehat{S}$ defined as Equation (5), we can compute the correction to energies as the diagonal part of $\frac{1}{2}[\widehat{S},\widehat{V}]$. The only term in $\widehat{V}$ that gives a diagonal contribution with the ${\widehat{a}}^{†\mathbf{n}}{h}_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}\widehat{a}){\widehat{a}}^{\mathbf{m}}$ term in $\widehat{S}$ is ${\widehat{a}}^{†\mathbf{m}}{g}_{\mathbf{m}\mathbf{n}}({\widehat{a}}^{†}\widehat{a}){\widehat{a}}^{\mathbf{n}}$. A pictorial representation of this statement is given in Figure 2. Hence,

$$\begin{array}{cc}\hfill \frac{1}{2}{[\widehat{S},\widehat{V}]}_{(d)}& =\frac{1}{2}\sum _{\mathbf{n},\mathbf{m}}[{\widehat{a}}^{†\mathbf{n}}{h}_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}\widehat{a}){\widehat{a}}^{\mathbf{m}},{\widehat{a}}^{†\mathbf{m}}{g}_{\mathbf{m}\mathbf{n}}({\widehat{a}}^{†}\widehat{a}){\widehat{a}}^{\mathbf{n}}]\hfill \end{array}$$

(12)

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill & =\frac{1}{2}\sum _{\mathbf{n},\mathbf{m}}\left\{\frac{{({\widehat{a}}^{†}\widehat{a})}^{\underline{\mathbf{n}}}{({\widehat{a}}^{†}\widehat{a}-\mathbf{n}+\mathbf{m})}^{\underline{\mathbf{m}}}{|g_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}a-\mathbf{n})|}^2}{f({\widehat{a}}^{†}\widehat{a})-f({\widehat{a}}^{†}\widehat{a}-\mathbf{n}+\mathbf{m})}\right.\hfill \\ \hfill & +\left.\frac{{({\widehat{a}}^{†}\widehat{a})}^{\underline{\mathbf{m}}}{({\widehat{a}}^{†}\widehat{a}+\mathbf{n}-\mathbf{m})}^{\underline{\mathbf{n}}}{|g_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}a-\mathbf{m})|}^2}{f({\widehat{a}}^{†}\widehat{a})-f({\widehat{a}}^{†}\widehat{a}+\mathbf{n}-\mathbf{m})}\right\}\hfill \end{array}\end{array}$$

(13)

$$\begin{array}{cc}\hfill & =\sum _{\mathbf{n},\mathbf{m}}\frac{{({\widehat{a}}^{†}\widehat{a})}^{\underline{\mathbf{n}}}{({\widehat{a}}^{†}\widehat{a}-\mathbf{n}+\mathbf{m})}^{\underline{\mathbf{m}}}{|g_{\mathbf{n}\mathbf{m}}({\widehat{a}}^{†}a-\mathbf{n})|}^2}{f({\widehat{a}}^{†}\widehat{a})-f({\widehat{a}}^{†}\widehat{a}-\mathbf{n}+\mathbf{m})}.\hfill \end{array}$$

(14)

The calculation of the commutator can be done by using Lemma 1 and results in Appendix B. Note that the summand in Equation (13) is symmetric under change of n, m, which leads to Equation (14). This result is equivalent to nondegenerate second-order perturbation energy correction,

$$E_{\mathbf{k}}^{(2)}=\sum _{{\mathbf{k}}^{\prime }\ne \mathbf{k}}\frac{{|\langle {\mathbf{k}}^{\prime }|\widehat{V}|\mathbf{k}\rangle |}^2}{E_{\mathbf{k}}-E_{{\mathbf{k}}^{\prime }}},$$

(15)

computed in our proposed ordering.

In most cases, interaction terms are of form ${\widehat{a}}_i^{†}{g}_{ij}({\widehat{a}}^{†}\widehat{a}){\widehat{a}}_j$, which represent a single transfer of quantum excitations. If we restrict the interaction to only these terms, our main results Equations (5) and (6) are simplified to

$$\begin{array}{cc}\hfill \widehat{S}& =\sum _{i\ne j}{\widehat{a}}_i^{†}\frac{g_{ij}({\widehat{a}}^{†}\widehat{a})}{f({\widehat{a}}^{†}\widehat{a}+e_i)-f({\widehat{a}}^{†}\widehat{a}+e_j)}{\widehat{a}}_j,\hfill \end{array}$$

(16)

$$\begin{array}{cc}\hfill \frac{1}{2}{[\widehat{S},\widehat{V}]}_{(d)}& =\sum _{i\ne j}\frac{{\widehat{a}}_i^{†}{\widehat{a}}_i({\widehat{a}}_j^{†}{\widehat{a}}_j+1){|g_{ij}({\widehat{a}}^{†}\widehat{a}-e_i)|}^2}{f({\widehat{a}}^{†}\widehat{a})-f({\widehat{a}}^{†}\widehat{a}-e_i+e_j)},\hfill \end{array}$$

(17)

where $e_i$ is the N tuple, which has 1 as its i-th component, and all other elements are 0 and $i,j\in \{1,2,\dots ,N\}$.

2.2. Application to Analyzing the Circuit Hamiltonian

We return to the diagonalization of the circuit Hamiltonian Equation (1). To apply the previous formalism, define functions $f,g_{12},g_{23}$ as

$$\begin{array}{cc}\hfill f(n,m,\ell )& :={\omega }_{s}n+{\omega }_{t}m+{\omega }_r\ell -\frac{E_C}{2\hslash }m(m-1),\hfill \end{array}$$

(18)

$$\begin{array}{cc}\hfill g_{12}(n,m,\ell )& =g_{21}(n,m,\ell ):=g_1,\hfill \end{array}$$

(19)

$$\begin{array}{cc}\hfill g_{23}(n,m,\ell )& =g_{32}(n,m,\ell ):=g_2.\hfill \end{array}$$

(20)

These functions give a full description of the Hamiltonian Equation (1). $n,m,\ell$ correspond to ${\widehat{a}}_s^{†}{\widehat{a}}_s$, ${\widehat{a}}_t^{†}{\widehat{a}}_t$, ${\widehat{a}}_r^{†}{\widehat{a}}_r$, respectively. The second-order energy corrections can be directly computed by our main result Equation (17).

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill \frac{{\widehat{H}}^{(2)}}{\hslash }& =|g_1{|}^2\left(\frac{n(m+1)}{f(n,m,\ell )-f(n-1,m+1,\ell )}+\frac{m(n+1)}{f(n,m,\ell )-f(n+1,m-1,\ell )}\right)\hfill \\ \hfill & +|g_2{|}^2\left(\frac{m(\ell +1)}{f(n,m,\ell )-f(n,m-1,\ell +1)}+\frac{\ell (m+1)}{f(n,m,\ell )-f(n,m+1,\ell -1)}\right)\hfill \end{array}\end{array}$$

(21)

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill & =|g_1{|}^2\left(\frac{{\widehat{a}}_s^{†}{\widehat{a}}_s({\widehat{a}}_t^{†}{\widehat{a}}_t+1)}{{\Delta }_{st}+{\widehat{a}}_t^{†}{\widehat{a}}_t{E}_C/\hslash }-\frac{{\widehat{a}}_t^{†}{\widehat{a}}_t({\widehat{a}}_s^{†}{\widehat{a}}_s+1)}{{\Delta }_{st}+({\widehat{a}}_t^{†}{\widehat{a}}_t-1)E_C/\hslash }\right)\hfill \\ \hfill & +|g_2{|}^2\left(\frac{{\widehat{a}}_r^{†}{\widehat{a}}_r({\widehat{a}}_t^{†}{\widehat{a}}_t+1)}{{\Delta }_{rt}+{\widehat{a}}_t^{†}{\widehat{a}}_t{E}_C/\hslash }-\frac{{\widehat{a}}_t^{†}{\widehat{a}}_t({\widehat{a}}_r^{†}{\widehat{a}}_r+1)}{{\Delta }_{rt}+({\widehat{a}}_t^{†}{\widehat{a}}_t-1)E_C/\hslash }\right),\hfill \end{array}\end{array}$$

(22)

with ${\Delta }_{st}:={\omega }_s-{\omega }_t$, ${\Delta }_{rt}:={\omega }_r-{\omega }_t$. To read off shifts in frequency, cross Kerr coefficients, and anharmonicities, we must put Equation (22) in normal order:

$$\frac{{\widehat{H}}^{(2)}}{\hslash }={\delta }_s{\widehat{a}}_s^{†}{\widehat{a}}_s+{\delta }_t{\widehat{a}}_t^{†}{\widehat{a}}_t+{\delta }_r{\widehat{a}}_r^{†}{\widehat{a}}_r+\frac{{\delta }_K}{2\hslash }{\widehat{a}}_t^{†}{\widehat{a}}_t^{†}{\widehat{a}}_t{\widehat{a}}_t+{\chi }_{st}{\widehat{a}}_s^{†}{\widehat{a}}_s{\widehat{a}}_t^{†}{\widehat{a}}_t+{\chi }_{rt}{\widehat{a}}_r^{†}{\widehat{a}}_r{\widehat{a}}_t^{†}{\widehat{a}}_t+\dots .$$

(23)

Using the result from Equation (A16), the shifts are given as

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill \frac{{\widehat{H}}^{(2)}}{\hslash }& =\frac{|g_1{|}^2}{{\Delta }_{st}}{\widehat{a}}_s^{†}{\widehat{a}}_s+\frac{|g_2{|}^2}{{\Delta }_{rt}}{\widehat{a}}_r^{†}{\widehat{a}}_r+\sum _{k=1}^{\infty }\frac{{(-1)}^k}{E_C}\left(\frac{|g_1{|}^2}{{({\Delta }_{st}/E_C)}^{\overline{k}}}+\frac{|g_2{|}^2}{{({\Delta }_{rt}/E_C)}^{\overline{k}}}\right){\widehat{a}}_t^{†k}{\widehat{a}}_t^k\hfill \\ \hfill & +\sum _{k=1}^{\infty }\frac{{(-1)}^k{|g_1|}^2(k+1)}{E_C{({\Delta }_{st}/E_C)}^{\overline{k+1}}}{\widehat{a}}_s^{†}{\widehat{a}}_t^{†k}{\widehat{a}}_t^k{\widehat{a}}_s+\sum _{k=1}^{\infty }\frac{{(-1)}^k{|g_2|}^2(k+1)}{E_C{({\Delta }_{rt}/E_C)}^{\overline{k+1}}}{\widehat{a}}_r^{†}{\widehat{a}}_t^{†k}{\widehat{a}}_t^k{\widehat{a}}_r,\hfill \end{array}\end{array}$$

(24)

where $x^{\overline{n}}=x(x+1)\dots (x+n-1)$ is the rising factorial and factors of ℏ were omitted in the right-hand side for simplicity. Restoring these factors are done by replacing $E_C$ with $E_C/\hslash$. The shifts in physical quantities are found by simply reading off the coefficients of Equation (24):

$$\begin{array}{cc}\hfill {\delta }_s& =\frac{|g_1{|}^2}{{\Delta }_{st}},{\delta }_r=\frac{|g_2{|}^2}{{\Delta }_{rt}},{\delta }_t=-\frac{|g_1{|}^2}{{\Delta }_{st}}-\frac{|g_2{|}^2}{{\Delta }_{rt}},\hfill \end{array}$$

(25)

$$\begin{array}{cc}\hfill {\delta }_K& =\frac{2|g_1{|}^2{E}_C}{{\Delta }_{st}({\Delta }_{st}+E_C/\hslash )}+\frac{2|g_2{|}^2{E}_C}{{\Delta }_{rt}({\Delta }_{rt}+E_C/\hslash )},\hfill \end{array}$$

(26)

$$\begin{array}{cc}\hfill {\chi }_{st}& =-\frac{2|g_1{|}^2{E}_C/\hslash }{{\Delta }_{st}({\Delta }_{st}+E_C/\hslash )},{\chi }_{rt}=-\frac{2|g_2{|}^2{E}_C/\hslash }{{\Delta }_{rt}({\Delta }_{rt}+E_C/\hslash )}.\hfill \end{array}$$

(27)

Hence, the total transformed Hamiltonian can be written as

$$\begin{array}{cc}\hfill e^{\widehat{S}}\frac{\widehat{\mathcal{H}}}{\hslash }& e^{-\widehat{S}}=\frac{1}{\hslash }\left({\widehat{H}}_0+{\widehat{H}}^{(2)}+\dots \right)\hfill \\ \hfill & =\sum _{i=s,t,r}({\omega }_i+{\delta }_i){\widehat{a}}_i^{†}{\widehat{a}}_i-\frac{E_C-{\delta }_K}{2\hslash }{\widehat{a}}_t^{†}{\widehat{a}}_t^{†}{\widehat{a}}_t{\widehat{a}}_t+{\chi }_{st}{\widehat{a}}_s^{†}{\widehat{a}}_s{\widehat{a}}_t^{†}{\widehat{a}}_t+{\chi }_{rt}{\widehat{a}}_r^{†}{\widehat{a}}_r{\widehat{a}}_t^{†}{\widehat{a}}_t+\dots .\hfill \end{array}$$

(28)

This extends the results using Bogoliubov approach to diagonalize the Hamiltonian of a coupled single LC oscillator and transmon [18] in the sense that the frequency shift of the transmon qubit is the sum of contributions from coupling to each LC oscillator. Such a system is described by a Hamiltonian

$$\widehat{\mathcal{H}}=\hslash {\omega }_1{\widehat{a}}^{†}\widehat{a}+\hslash {\omega }_2{\widehat{b}}^{†}\widehat{b}-\frac{E_C}{2}{\widehat{b}}^{†}{\widehat{b}}^{†}\widehat{b}\widehat{b}+\hslash g({\widehat{a}}^{†}\widehat{b}+\widehat{a}{\widehat{b}}^{†}).$$

(29)

Elimination of the $\hslash g({\widehat{a}}^{†}\widehat{b}+\widehat{a}{\widehat{b}}^{†})$ term gives rise to cross Kerr coefficient between $\widehat{a},\widehat{b}$,

$$\chi =-\frac{{2|g|}^2{E}_C/\hslash }{\Delta (\Delta +E_C/\hslash )},$$

(30)

where $\Delta :={\omega }_1-{\omega }_2$, which highly resembles the results in Equation (27). Note that there is a sign difference in the definition of $\Delta$ compared with that of Ref. [18].

To obtain a form similar to that given in Refs. [13,14], we truncate the transmon Hilbert space to the first two levels. Such truncation is done by replacing ${\widehat{a}}_t^{†}{\widehat{a}}_t$ with $({\sigma }_z+1)/2$ in Equation (22). The result is

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill \frac{{\widehat{H}}^{(2)}}{\hslash }& =\frac{|g_1{|}^2}{{\Delta }_{st}+E_C/\hslash }{\widehat{a}}_s^{†}{\widehat{a}}_s-\left(\frac{|g_1{|}^2}{{\Delta }_{st}}+\frac{|g_2{|}^2}{{\Delta }_{rt}}\right)\frac{{\sigma }_z}{2}+\frac{|g_2{|}^2}{{\Delta }_{rt}+E_C/\hslash }{\widehat{a}}_r^{†}{\widehat{a}}_r\hfill \\ \hfill & -\frac{|g_1{|}^2{E}_C/\hslash }{{\Delta }_{st}({\Delta }_{st}+E_C/\hslash )}{\widehat{a}}_s^{†}{\widehat{a}}_s{\sigma }_z-\frac{|g_2{|}^2{E}_C/\hslash }{{\Delta }_{rt}({\Delta }_{rt}+E_C/\hslash )}{\widehat{a}}_r^{†}{\widehat{a}}_r{\sigma }_z\hfill \end{array}\end{array}$$

(31)

up to an overall constant. Again, the contributions from each oscillator-transmon coupling stated in Ref. [13] are added independently. The frequency shifts of LC oscillators seem to be different compared with Equation (25), but this is due to a subtle difference of physical interpretation. The coefficients ${\delta }_s$, ${\delta }_r$ in Equation (25) are the frequency shifts when the transmon is in the ground state, while the coefficients of ${\widehat{a}}_s^{†}{\widehat{a}}_s$, ${\widehat{a}}_r^{†}{\widehat{a}}_r$ in Equation (31) are the average of the frequency shifts when the transmon is in the ground state and excited state. With the ${\widehat{a}}_s^{†}{\widehat{a}}_s{\sigma }_z$, ${\widehat{a}}_r^{†}{\widehat{a}}_r{\sigma }_z$ terms in consideration, both Equations (22) and (31) give the same energy spectrum when considering up to the first excitation of the transmon.

3. Discussion

We proposed an ordering of bosonic operators to efficiently compute the Schrieffer–Wolff transformation generator and energy corrections. This formalism was applied to a system with a transmon coupled to two different LC oscillators to model a quantum memory and readout device demonstrated in Ref. [8]. We solved the normal ordering problem for an operator that appears in the second-order energy correction, so that shifts in physical parameters such as frequency, anharmonicity, and cross Kerr coefficients can be directly read off from the normal-ordered form.

Our proposed method can be directly applied to systems consisting of LC circuits coupled with multiple transmons, and even to systems that have nonlinear couplings provided that the couplings represent a definite number of excitations or de-excitations in the Fock basis. It is possible to generalize this method to incorporate fermionic operators in this formalism, which can be used to reproduce the results of the original application of Schrieffer–Wolff transformation to the Anderson impurity model as in Appendix C. With such a general method, one can analyze a wide range of time-independent systems.

Quantum illumination, an example of quantum information technology, uses entangled light to achieve higher detection rate of a target with low-reflectivity. The idler mode, a part of the entangled light, should be stored in a quantum memory for ideal operation. Our method was used to analyze a demonstrated quantum memory and can be used to analyze other systems operating in various quantum technologies. For further research, it is required to find methods to store and release the propagating idler mode efficiently [19], leading to applications of quantum memories to quantum illumination.