Many-body quantum vacuum fluctuation engines

The minimal quantum system featuring an entangled ground state consists of two coupled qubits. As such, we consider two qubits, denoted by $\mathrm{A}$ and $\mathrm{B}$, coupled to one another via the horizontal components of their spins. Denoting by $\omega_{j}$ qubit $j$’s transition frequency, and by $g$ the interqubit coupling strength, the local and interaction Hamiltonians read as (we set $\hbar=1$)

where we have introduced the raising and lowering operators for qubit $j$, $\sigma_{j}^{+}=\ket{1_{j}}\!\bra{0_{j}}$ and $\sigma_{j}^{-}=\ket{0_{j}}\!\bra{1_{j}}$, as well as its spin’s $x$ component, $\sigma_{j}^{x}=\sigma_{j}^{+}+\sigma_{j}^{-}$. We note that $\braket{0_{j}|\sigma_{j}^{x}|0_{j}}=\braket{1_{j}|\sigma_{j}^{x}|1_{j}}=0$, so the interaction terms can be turned off with no expected energy cost when the system is in a separable state.

As a result of the interaction between qubits $\mathrm{A}$ and $\mathrm{B}$, the eigenstates of the total Hamiltonian $H=H_{\mathrm{loc}}+H_{\mathrm{int}}$ are entangled; we write them as

where the angles $\varphi$ and $\psi$ satisfy

with $\gamma=g/(\omega_{\mathrm{A}}+\omega_{\mathrm{B}})$ and $\delta=(\omega_{\mathrm{A}}-\omega_{\mathrm{B}})/(\omega_{\mathrm{A}}+\omega_{\mathrm{B}})$, $0<\delta<1$. The corresponding eigenenergies are given by

with $E_{\varphi}^{-}<E_{\psi}^{-}<E_{\psi}^{+}<E_{\varphi}^{+}$. Using the notations in Table 1, we write $\ket{\tilde{0}}=\ket{\varphi^{-}}$ and $\ket{0_{\mathrm{loc}}}=\ket{00}$. The definition of the local Hamiltonian in Eq. (9) sets the local ground-state energy at zero: $E_{0_{\mathrm{loc}}}=0$; the local entanglement gap is then given by $\Delta=-E_{\tilde{0}}=-E_{\varphi}^{-}$.

The engine cycle described in Sec. II applied to the two-qubit system considered here is depicted in Fig. 2. One should note that a projective measurement on one the qubits, say $\mathrm{A}$, in the local eigenbasis $\{\ket{0},\ket{1}\}$ suffices to fully characterize the state of the two-qubit system: “0” corresponds to the two-qubit state $\ket{00}$, and similarly “1” corresponds to $\ket{11}$. In the latter case, we turn off interactions with no energetic cost since $\braket{11|H_{\mathrm{int}}|11}=0$, and we apply local pulses to each qubit so as to bring the two-qubit system to the local ground state $\ket{00}$. In doing so, work is extracted as both qubits emit a stimulated photon at their resonant frequency. Using Eqs. (14) and (15), we find the probability for each outcome of the local measurement:

Once the two-qubit system is in its local ground state, we turn interactions back on with no energetic cost given that $\braket{00|H_{\mathrm{int}}|00}=0$, and we couple the system to a cold bath so that it relaxes to its entangled ground state. The cycle can then resume. We note that to carry out the projective measurements in the local basis, an external system must couple strongly and quickly. This process is modeled in Appendix A and we show meter coupling around a factor of $10$ larger than the system energies is sufficient to carry out an approximate projective measurement.

Using the general results in Table 1, one can readily obtain the two-qubit engine work and efficiency. The average work output is given by

while the quantum heat necessary for the engine to operate reads as

This yields the efficiency:

The engine’s work output and efficiency are plotted as functions of $\gamma$ in Fig. 3. In the weak coupling limit where $\gamma\ll 1$, that is, $g\ll\omega_{\mathrm{A}}+\omega_{\mathrm{B}}$, the interaction Hamiltonian is negligible as compared with the local Hamiltonian. This means that the many-body ground state is almost equal to the local ground state: $\ket{\tilde{0}}\simeq\ket{0_{\mathrm{loc}}}=\ket{00}$. As a consequence, the probability that the local measurement yields the excited state $\ket{11}$ is vanishingly small, and almost every realization of the cycle results in no work extracted. The fact that the many-body and local ground states almost coincide also implies that the local entanglement gap closes: $\Delta\simeq 0$. Coincidentally, the work output and local entanglement gap have the same asymptotic behavior as $\gamma$ vanishes. The efficiency then tends to $1/2$—its maximum value—in the weak coupling limit.

As $\gamma$ increases, the work output increases and the efficiency decreases. The trade-off between these quantities is described through the relation

which is plotted in the inset of Fig. 3.

In the deep strong coupling limit where $\gamma\gg 1$, or equivalently $g\gg\omega_{\mathrm{A}}+\omega_{\mathrm{B}}$, the interacting ground state becomes a Bell state: $\ket{\tilde{0}}\simeq(\ket{00}-\ket{11})/\sqrt{2}$. As such, the local ground state $\ket{00}$ and the excited state $\ket{11}$ are equiprobable outcomes of the local measurement. Work extraction then occurs for half of the cycle’s realizations, so the average work output is $(\omega_{\mathrm{A}}+\omega_{\mathrm{B}})/2$. This is the upper bound for the engine’s average work output. Conversely, since the local entanglement gap grows linearly with $\gamma$ in the deep strong coupling limit, the efficiency vanishes.

Finally, fluctuations of the work output between different realizations of the cycle are described through the standard deviation $\sigma$. We obtain (see Table 1)

In the weak coupling regime ($\gamma\ll 1$), the standard deviation scales linearly with $\gamma$ while the average extracted work is quadratic, so fluctuations swamp the average in this limit. In contrast, both the work output’s average and standard deviation saturate to the value $(\omega_{\mathrm{A}}+\omega_{\mathrm{B}})/2$ in the deep strong coupling regime ($\gamma\gg 1$).

To compute the power generated through a realization of the cycle, it is necessary to model its dynamics, which is discussed in Appendix A.

We extend the working principle of the two-qubit engine to $N$ coupled qubits. We will focus here on the case of a closed chain with periodic boundary conditions (Appendix B provides some details about the case of an open chain). We then consider a closed chain of $N$ qubits with transition frequencies $\omega$ coupled to their nearest neighbors with the coupling amplitude $g$. The local and interaction Hamiltonians are then written as

with $\sigma_{N+1}^{x}=\sigma_{1}^{x}$. The procedure considered previously to extract work using local operations can be applied here as well: The chain, initially in its entangled ground state, is projected onto the local eigenbasis. To do so, local measurements must be performed on all but one of the qubits within the chain as can be seen in Eq. (35) below. The interaction Hamiltonian is then turned off at no energetic cost. Each qubit found in an excited state is flipped using a local pulse, which transfers energy to the field driving the pulse. With the chain now in its local ground state, interactions are turned back on at no energetic cost, and the system is coupled to a cold environment which relaxes it to its many-body ground state. The cycle can then restart.

The total Hamiltonian $H=H_{\mathrm{loc}}+H_{\mathrm{int}}$ is analogous to the transverse-field Ising model, where $\omega$ and $g$ would respectively represent the intensity of the transverse magnetic field and the exchange parameter characterizing the interaction between neighboring spins. This connection to quantum magnetism provides helpful tools to assess the qubit chain engine’s performance. The work output is given by the expectation value of the local Hamiltonian in the ground state (see Table 1). In the Ising model language, this corresponds to the energy from the transverse magnetic field which is proportional to the transverse magnetization. To further obtain the engine’s efficiency, it is necessary to calculate the local entanglement gap which is readily deduced from the ground-state energy (see Table 1). One can then fully characterize the performance of a many-qubit engine from the transverse magnetization and ground-state energy in the corresponding Ising model.

Furthermore, the one-dimensional transverse-field Ising model considered here undergoes a quantum phase transition [29, 30, 31, 32]. With the conventions used to define the Hamiltonian in Eqs. (26) and (27), this is a transition from a diamagnetic phase ($g<\omega$) to an antiferromagnetic phase ($g>\omega$) as depicted in Fig. 4. The magnetic order in both of these phases correspond to different behaviors of the longitudinal magnetization, but other quantities also bear the signature of the phase transition. This includes the transverse magnetization and ground-state energy from which the engine’s work output and efficiency are inferred.

An analytical solution for the one-dimensional transverse-field Ising model can be obtained using a Jordan–Wigner transformation which maps spins to fermions [29, 30, 31, 32]. We then introduce the fermionic operators

where the first term to the right-hand side of the equation above ensures that fermionic anticommutation relations are satisfied. We next move to momentum space [30, 31, 32]:

and obtain

The Hamiltonian is then diagonalized by performing the adequate Bogoliubov transformation: $\zeta_{p}=u_{p}c_{p}+\mathrm{i}v_{p}c_{-p}^{\dagger}$, with¹¹ 1 If $p$ is a multiple of $\pi$, the Bogoliubov transformation is unnecessary; we then have $u_{p}=1$, $v_{p}=0$, and $\Omega_{p}=\omega+g\cos p$.

The fact that $u_{p}^{2}+v_{p}^{2}=1$ ensures that fermionic commutation relations are satisfied for the quasiparticle operators: $\{\zeta_{p},\zeta_{q}^{\dagger}\}=\delta_{pq}$, $\{\zeta_{p},\zeta_{q}\}=0$. The transverse-field Ising Hamiltonian can then be rewritten as a free fermion Hamiltonian:

It must be noted that the set of possible momentum values in the equation above depend on the number of excitations [31], or, equivalently, the number of Jordan–Wigner fermions, $n=\sum_{j}\sigma_{j}^{+}\sigma_{j}^{-}=\sum_{j}c_{j}^{\dagger}c_{j}$. Namely, $p=(2m-1)\pi/N$ for an even number of excitations, while $p=2m\pi/N$ for an odd number of excitations, where $m$ is an integer ranging from $-\lfloor(N-1)/2\rfloor$ to $\lfloor N/2\rfloor$. More details about the derivation of Eq. (34) are given in Appendix C.

The ground state of $H$ is the quasiparticle vacuum, that is, the state $\ket{\tilde{0}}$ such that $\zeta_{p}\ket{\tilde{0}}=0$ for all $p$. We find that it can be written as

where the local ground state $\ket{0_{\mathrm{loc}}}=\otimes_{j}\ket{0}$ corresponds to the vacuum state for Jordan–Wigner fermions, that is, the state where all qubits are in their individual ground state. The entangled ground state $\ket{\tilde{0}}$ in Eq. (35) above is a linear combination of states with an even number of excitations. As a consequence, the state of the $N$-qubit chain in the local basis can be deduced from local measurements on $N-1$ of the qubits within the chain. This is similar to the two-qubit case where it was sufficient to measure one of the qubits as depicted in Fig. 2.

One can then read the ground-state energy $E_{\tilde{0}}$ from Eq. (34). Since the local Hamiltonian in Eq. (26) has been defined such that its ground-state energy $E_{0_{\mathrm{loc}}}$ is set at zero, we find that the local entanglement gap is given by (see Table 1)

As the entangled ground state $\ket{\tilde{0}}$ in Eq. (35) is a linear combination of states with an even number of excitations, the values of the momentum $p$ in the summation of Eq. (36) above are such that: $p=(2m-1)\pi/N$, with $m$ an integer between $-\lfloor(N-1)/2\rfloor$ and $\lfloor N/2\rfloor$.

We now calculate the engine’s work output which, as shown in Table 1, is equal to the expectation value for the local Hamiltonian in the interacting ground state. We then rewrite the local Hamiltonian in terms of quasiparticle operators:

where we have expressed the number operator in momentum space: $n=\sum_{j}\sigma_{j}^{+}\sigma_{j}^{-}=\sum_{p}c_{p}^{\dagger}c_{p}$, and inverted the Bogoliubov transformation. We then obtain

The engine’s efficiency can now be deduced from the work output and the local entanglement gap using the relation $\eta=W/(W+\Delta)$ from Table 1.

Analytical results can be obtained in the thermodynamic limit where $N\to\infty$: The summations in Eqs. (36) and (38) become elliptic integrals, and we find

where $K(x)$ and $E(x)$ are the complete elliptic integrals first and second kinds:

Both $\Delta$ and $W$ above are nonanalytic at $g=\omega$: the second derivative of $\Delta$ and the first derivative of $W$ diverge at this point. As a consequence, the engine’s work output and efficiency both exhibit a vertical tangent at the critical point $g=\omega$. Nevertheless, these quantities have well-defined critical values:

One can gain insights about the scaling of work and efficiency with the number of qubits by analyzing the limiting cases $g\ll\omega$ (weak coupling limit, diamagnetic phase) and $g\gg\omega$ (deep strong coupling limit, antiferromagnetic phase). In the weak coupling limit $g\ll\omega$, we find that the local entanglement gap and the work output behave in a similar way, scaling linearly with the number of qubits. Indeed, for $N>2$, we have²² 2 We refer the reader to Sec. III.1 for the two-qubit case.

As a consequence, the efficiency is independent of the number of qubits and is simply given by: $\eta\simeq 1/2$.

In the deep strong coupling limit $g\gg\omega$, we find that the local entanglement gap and work output again scale linearly with the number of qubits:

Again, the efficiency is independent on the number of qubits here: $\eta\simeq\omega/g$. The discrepancy observed between even and odd values of $N$ above is a consequence of geometrical frustration in the antiferromagnetic phase: When $g\gg\omega$, the interaction Hamiltonian dominates which favors the antialignment of neighboring spins. For an even number of spins, it is possible to find a configuration where each spin is antialigned with its neighbors, while such a configuration does not exist for an odd number of spins. This results in different structures for the antiferromagnetic ground state in each case [31]. This difference manisfests in the quasiparticle picture of Eq. (34) through the fact that one of the quasiparticle energies, namely, $\Omega_{\pi}=\omega-g$, is negative for $g>\omega$ when $N$ is odd. In contrast, all quasiparticle energies remain positive when $N$ is even; this is because $p=\pi$ is not a relevant value for the momentum in this case.

As shown in Fig. 5, the observations made above concerning the engine’s work output and efficiency still hold in the general case: The work output scales linearly with the number of qubits, while the efficiency is essentially constant. One should also note that both density plots in Fig. 5 exhibit a horizontal stripe pattern which corresponds to the dependence of work and efficiency on the parity of the number of qubits within the chain. More precisely, the different results obtained for an even or odd number of qubits are a consequence of geometric frustration in the antiferromagnetic phase as explained previously. According to Eqs. (46) and (47), such a difference is of order $1/N$.

Fig. 5 also highlights the impact of the quantum phase transition occurring at $g=\omega$ on the engine’s performance. One can clearly distinguish two regimes for the work output (upper panel): The work output vanishes for $g=0$ and remains small in the diamagnetic phase ($g<\omega$). It then abruptly increases around the critical point, with a vertical tangent at $g=\omega$. In the antiferromagnetic phase ($g>\omega$), the work output plateaus as it approaches its maximum value given in Eq. (47). As for the efficiency, it is essentially constant in the diamagnetic phase ($g<\omega$) with $\eta\approx 1/2$. Similar to the work output, the quantum phase transition is marked by an abrupt increase of the efficiency, with a vertical tangent at $g=\omega$. The efficiency peaks for $g\gtrsim\omega$, and then decreases in the antiferromagnetic phase ($g>\omega$) to eventually vanish in the deep strong coupling limit. We note that the observation of a clear trade-off between work and efficiency made in the case of two coupled qubits (see Fig. 3) cannot be generalized to a longer chain: For a larger number of coupled qubits, maximum efficiency is reached at nonzero coupling. Interestingly, we find that the work output at maximum efficiency is not zero in this case, and it is in fact relatively close to its maximum value. This is a consequence of the sharp increase of the efficiency around $g=\omega$ caused by the divergence of its first derivative at the critical point. Our results corroborate recent analyses which showed that the thermodynamic performance of various quantum machines operating close to a quantum critical point bears the signature of the corresponding phase transition [33, 34, 35, 36, 37]. In particular, an enhancement of efficiency around the critical point has been noted in Refs. [33, 34, 35].

We conclude this section by addressing the question of fluctuations. The work output standard deviation $\sigma$ is calculated using the relation in Table 1 and the expression for the local Hamiltonian in Eq (37):

Once again, we scrutinize the weak and deep strong coupling limits to obtain insights about the scaling of the standard deviation with the number of qubits within the chain. We find

which suggests that $\sigma$ scales like $\sqrt{N}$. In fact, the results in Eqs. (49) and (50) above accurately approximate the standard deviation within the whole diamagnetic ($g<\omega$) and antiferromagnetic ($g>\omega$) phase respectively. In the thermodynamic limit $N\to\infty$, the approximation matches the exact result which is obtained by replacing the summation in Eq. (48) by an integral:

We note that the only effect of the phase transition on the standard deviation is a slope discontinuity at the critical point. Further, we observe that $\sigma$ scales as $\sqrt{N}$ which shows that the fluctuations of the work output become negligible in comparison to its average value for $N\to\infty$.

In the previous section, we showed the physics of the qubit chain could be mapped to free fermions and solved. Here, we consider the bosonic case, and we begin with the simplest version of two coupled harmonic oscillators:

where $x_{j}$ and $p_{j}$ respectively denote the position and momentum for oscillator $j$, $k_{0}$ is the force constant for both oscillators, and $g$ is the coupling constant. The Hamiltonian in Eq. (52) above can be written as $H=H_{\mathrm{loc}}+H_{\mathrm{int}}$, where the local and interaction Hamiltonians are given by

We denote by $\ket{n_{1},n_{2}}$ the local eigenstates, where $n_{j}$ is the energy quantum number for local oscillator $j$. Obviously, the local ground state is $\ket{0_{\mathrm{loc}}}=\ket{0,0}$, and its energy is given by the local oscillator frequency: $E_{0_{\mathrm{loc}}}=\omega=\sqrt{k_{0}+g}$. Furthermore, we note that the local eigenstates satisfy the condition $\braket{n_{1},n_{2}|H_{\mathrm{int}}|n_{1},n_{2}}=0$ because of parity, indicating that there is no expected energy in the coupling term.

Introducing the sum and difference variables, $x_{\pm}=(x_{1}\pm x_{2})/\sqrt{2}$, and quantizing the total Hamiltonian $H=H_{\mathrm{loc}}+H_{\mathrm{int}}$ in the standard fashion, it takes the form of two independent oscillators:

with the natural frequencies $\omega_{+}=\sqrt{k_{0}}$ and $\omega_{-}=\sqrt{k_{0}+2g}$. In this new basis, it is clear that the ground-state energy is $E_{\tilde{0}}=(\omega_{+}+\omega_{-})/2$. The associated ground-state wave function is given by

Following the engine cycle considered in previous sections, we now make local energy measurements of each oscillator, projecting each of them onto local energy eigenstates with energies $E_{n_{1},n_{2}}=(n_{1}+n_{2}+1)\omega$. If oscillator $j$ is found in an excited state $\ket{n_{j}}$ with $n_{j}\neq 0$, the excess energy can be transferred to a battery with local operations which lower the oscillator to its ground state. Once in the local ground state $\ket{0_{\mathrm{loc}}}=\ket{0,0}$, the interaction with the cold environment relaxes the coupled oscillators to the entangled ground state $\ket{\tilde{0}}$, and the cycle resets.

The engine’s work and efficiency can be obtained using the relations in Table 1. The average work output of the engine corresponds to the energy transferred from the local oscillators to the idealized energy source during local operations. It is given by

The Gaussian nature of the integrals permits an exact calculation of the expectation value:

Thus, the total work is

The amount of quantum heat given to the system by the measurement is the difference of the energy before and after the measurement, given, on average, by

The efficiency of the engine, defined as $\eta=W/Q$, is controlled by the local entanglement gap:

As for fluctuations, we calculate the work output standard deviation $\sigma$, and we find (see Table 1)

The work output and efficiency are plotted in contour plots for varying $k_{0}$ and $g$ in Fig. 6. Interestingly, both quantities increase as $k_{0}\to 0$ and $g\to\infty$. In this limit, $W\simeq g/4\sqrt{k_{0}}$ while $\eta\simeq 1-(4-2\sqrt{2})\sqrt{k_{0}/g}$. We see that a perfect efficiency engine is possible in this situation with the caveat that the standard deviation diverges since $\sigma\simeq g/2\sqrt{2k_{0}}$.

It is also of interest to find the explicit probabilities of finding the each excitation $\ket{n_{1},n_{2}}$ as the result of a local measurement on the system in its entangled ground state:

We consider the state overlap

where $H_{n}(z)$ are the Hermite polynomials. These state overlaps can be found with a joint generating function,

Using the identity $\sum_{n=0}^{\infty}H_{n}(x)t^{n}/n!=\mathrm{e}^{2xt-t^{2}}$, we find the result

where

We now generalize the previous treatment to an array of $N$ oscillators, linearly coupled to any other oscillator in the array. We follow the analysis of a related problem in Ref. [38]. The system is described by the Hamiltonian,

where the sums all range from $1$ to $N$, and $K$ is a real symmetric matrix with positive eigenvalues. The Hamiltonian in Eq. (69) above is broken into a local and interacting parts as follows

The matrix $K$ can be diagonalized by an orthogonal matrix $O$ as $K=O^{\mathrm{T}}K_{\mathrm{D}}O$, and the diagonal matrix $K_{\mathrm{D}}$ has positive eigenvalues $\{k_{m}\}$, so $\Omega_{m}=\sqrt{k_{m}}$ is the natural frequency of normal mode $m$. The normalized many-body ground-state wave function is then given by

where $\Omega=O^{\mathrm{T}}K_{\mathrm{D}}^{1/2}O$ is the square root of $K$, and the vector $x$ has components $x_{j}$. This is generally an entangled state.

The engine protocol described in the previous sections can be generalized here: Local energy measurements are performed on each local oscillator, with resulting energy

where $\sqrt{K_{jj}}$ is the natural frequency of local oscillator $j$. The globally entangled nature of the ground state permits the possibility of finding the local oscillators in locally excited states [25]. The interaction Hamiltonian in Eq. (71) is then turned off with no energetic cost. Any excess energy for a local oscillator found in an excited state is transferred to a battery with local operations, extracting its energy and lowering the oscillator to its ground state. Once the oscillator array is in its local ground state $\ket{0_{\mathrm{loc}}}=\ket{n_{1}=0,\dots,n_{N}=0}$, interactions are turned back on with no energetic cost. The system is finally let to relax to its many-body ground state, closing the cycle.

As shown in Table 1, the extracted work is given by

where the local ground-state energy is

The expectation value in Eq. (74) may be calculated using multidimensional Gaussian integration to find

which is our main result. In the special case where $K$ is diagonal, corresponding to a situation without entanglement in the ground state, $\Omega_{jj}=\sqrt{K_{jj}}$, and we recover the vacuum energy of the local oscillators, resulting in impossible work extraction.

The quantum heat is given by the mean energy difference before and after the measurement (see Table 1),

which is shown to be greater than or equal to 0 in Appendix D. The efficiency is given by

where $\Delta=E_{0_{\mathrm{loc}}}-E_{\tilde{0}}$ is again the local entanglement gap. A generalized expression for the generating function of all local energies can be found of Gaussian form in the generating variables. The probabilities are given explicitly in terms of a matrix Hafnian in Ref. [39].

The work output standard deviation $\sigma$ can also be calculated exactly using similar techniques. Using the relation in Table 1, we eventually find

In the case of two oscillators discussed in the previous section, we have

Inserting these formulas into Eqs. (76) and (79) respectively recovers the expression for $\braket{H_{\mathrm{loc}}}_{\tilde{0}}$ from Eq. (58) and the expression for $\sigma$ from Eq. (62).

We consider a simple model of a linear chain of masses on springs with nearest-neighbor coupling to see how the work and efficiency scale with $N$. The symmetric, tridiagonal coupling matrix is given by $K=k_{0}(2I-T)$, where $T$ has $1$s on the first off-diagonals and $0$s on the main diagonal. In this case, the local expected energy takes the form

The tridiagonal coupling matrix $K$ can be diagonalized exactly [40], with eigenvalues

where $m=1,\dots,N$ delineates the mode number. The corresponding natural frequencies are

Thus, we can express the local energy as

The second sum in Eq. (84) above can be calculated exactly:

but the first one does not have a closed form solution. Therefore, we examine the limit $N\to\infty$ and approximate the sum as an integral. We find the asymptotic behavior: