Differential Parametric Formalism for the Evolution of Gaussian States: Nonunitary Evolution and Invariant States

The study of Gaussian states has been of an essential interest in the last decades. This type of states, associated to classical random fields, were considered as a possibility to connect covariance matrices of the states as quantum density matrices and, with this definition, to study the quantum–classical relation of randomness with the quantization procedure [1, 2]. The problems of new developments of foundations of quantum mechanics and applications of new results in quantum information and quantum probabilities, as well as in areas like mathematical finance and economics, attract the attention of the researchers; they are intensely discussed in the literature [3, 4, 5]. The important role in this development is played by discussing the problems which appeared from the very beginning of quantum mechanics, like the notion of quantum system states and interpretation of the states associated in conventional formulation of quantum mechanics with Hilbert space vectors and density operators, using the quasiprobability distributions and the probability distributions containing the complete information on quantum states. There exists increasing interest to quantum foundations since a deeper understanding the essence and formalism of quantum theory is needed for the development of quantum technologies and possibilities to extend the applications of quantum formalism in physics to all other areas of science like economy, finance, and social disciplines.

Some examples of Gaussian states of quantum fields as the coherent, squeezed, and thermal light states are regularly used in the theoretical and experimental framework of quantum mechanics, optics, information, and computing. The use of these states in quantum information has been of particular importance [6, 7, 8]. One can list some of the most recent applications of the use of Gaussian systems – it has been demonstrated [9] that it is not possible to distill more entanglement from a bipartite Gaussian state, using local Gaussian transformations. In [10], several properties of the purity of Gaussian states were found. The connection between the symplectic invariants of bipartite Gaussian states, the von Neumann entropy, and the mutual information was established in [11]. The extremality of entanglement measures and secret key rates for Gaussian states was observed in [12]. It was shown [13] that Gaussian attacks are characterized by an optimum efficiency against eavesdrop protocols. Quantum illumination of a target using Gaussian light states was explored by Tan et al [14]. A quantum discord for systems of continuous variables, such as Gaussian states, was implemented in [15]. In [16], an invariant describing the nonclassicality in a two-mode Gaussian state was reported. The entanglement of $m$ modes with other $n$ modes of a Gaussian multipartite system was treated in [17]. The linear response for systems close to steady states under Gaussian processes was obtained in [18]. The optimal measurement of the fidelity of multimode Gaussian states has been studied in [19]. On the other hand, the study of Gaussian wave packets by nonlinear differential equations, as the Riccati equation, has been studied in [20, 21, 22]. Several coherent states have been defined by the use of quadratic operators [23]. The behavior of different quantities as covariances in thermal relaxation phenomena has been also studied in [24].

The aim of this work is to present a new way to characterize the dynamics of Gaussian states using the differential equations for the parameters, which determine their continuous variable density matrices. The proposed method makes use of the integrals of motion of such systems, and it can be used to clarify new aspects of multimode Gaussian quantum states, such as an explicit form of nonunitary evolution of the states of subsystems and the existence of invariant states with constant covariance matrices and mean values.

The time evolution of a quantum system was first established by Schrödinger [25]. The dynamics of the system given by a Hamiltonian operator $\hat{H}$ for a pure state $|\psi(t)\rangle$ must follow the Schrödinger equation

$$\hat{H}|\psi(t)\rangle=i\hbar\frac{\partial}{\partial t}|\psi(t)\rangle\,;$$

this expression corresponds to a second-order differential equation in the position representation. In the case of an arbitrary state represented by the density matrix $\hat{\rho}(t)$, which can be not pure [26, 27], the evolution is determined by the von Neumann equation

$$i\hbar\frac{\partial}{\partial t}\hat{\rho}(t)=[\hat{H},\hat{\rho}(t)]\,;$$

the general solution of this equation is given by the unitary transform $\hat{U}(t)$, i.e., $|\psi(t)\rangle=\hat{U}(t)|\psi(0)\rangle$ or $\hat{\rho}(t)=\hat{U}(t)\hat{\rho}(0)\hat{U}^{\dagger}(t)$, where $|\psi(0)\rangle$ and $\hat{\rho}(0)$ describe the system at time $t=0$. Also it is a common knowledge that, when the system interacts with an environment, its dynamics is described by the master equation [28, 30, 29].

The Gaussian states can be determined by their covariance matrix $\boldsymbol{\sigma}$ and mean values $\langle q_{j}\rangle$ and $\langle p_{j}\rangle$. This property also implies that the evolution of a Gaussian state can be obtained, if the time dependence of these parameters is known. In this work, we review the differential equations to which the covariance matrix and the mean values satisfy [31, 32]; employing these results we can define differential equations for the density matrix parameters of a general multimode state satisfying these equations, and then use the equations to discuss some physical characteristics of the unitary and nonunitary evolutions of Gaussian states.

This paper is organized as follows.

In section 2, the evolution of non-pure Gaussian states for a one-dimensional quadratic Hamiltonian is presented. To obtain this evolution, we make use of the derivatives of the covariance matrix, the mean values, and the parameters of the density operator; also we define and obtain invariant states for this system. The generalization of these results to the case of a multidimensional quadratic system is explored in section 3. As examples of the application of the general results to the nonunitary evolution of the subsystems of a two-mode state, as well as the definition of invariant and quasi-invariant states is done in section 4. Also in section 5, we obtain new invariant states for the frequency converter and quasi-invariant states for the parametric amplifier. The detection of this invariant states using the quantum tomographic representation of the states is discussed for single-mode Gaussian states in section 5 and for the bipartite system in section 6; in these sections, the correspondence between the time-independent states and thermal density matrices is mentioned. Finally, we give our conclusions.

In this section, we analyze some properties of the one-dimensional quadratic Hamiltonian. Particularly, we are interested in the invariant operators, which in the quadratic case happen to be linear in the quadrature operators $\hat{p}$ and $\hat{q}$.

The most general (in a unit system where $\hbar=m=1$), one-dimensional quantum quadratic Hamiltonian can be obtained in terms of the quadrature operators $\hat{p}$ and $\hat{q}$ as follows:

$$\hat{H}=(\hat{p},\hat{q})\left(\begin{array}[]{cc}\omega_{1}(t)&\omega_{2}(t)\\ \omega_{2}(t)&\omega_{3}(t)\end{array}\right)\left(\begin{array}[]{cc}\hat{p}\\ \hat{q}\end{array}\right)+(\hat{p},\hat{q})\left(\begin{array}[]{cc}\delta_{1}\\ \delta_{2}\end{array}\right)\,,$$

(1)

where the parameters $\omega_{1}(t)$, $\omega_{2}(t)$, $\omega_{3}(t)$, and $\delta_{1,2}$ are real functions of time. The dynamics associated to this Hamiltonian can be solved by different methods. One of them is the method of time-dependent invariants (integrals of motion) [33, 34]. These invariants are quantum operators $\hat{R}(t)$, whose total time derivative is equal to zero $\dfrac{d\hat{R}(t)}{dt}=0$. In the quadratic case, it is known that there exist invariants linearly depending on the quadrature operators $\hat{p}$ and $\hat{q}$, i.e., $\hat{R}(t)=\lambda_{1}(t)\hat{p}+\lambda_{2}(t)\hat{q}+\lambda_{3}(t)$.

By substituting this expression into the von Neumann equation, which determines the dynamics of $\hat{R}$, i.e., $\dfrac{d\hat{R}(t)}{dt}=\dfrac{i}{\hbar}[\hat{H}(t),\hat{R}(t)]+\dfrac{\partial\hat{R}(t)}{\partial t}=0$, one can show that $\hat{R}(t)$ is an invariant operator, if the following differential equations are satisfied:

$\displaystyle\dot{\lambda}_{1}=2(\omega_{2}\lambda_{1}-\omega_{1}\lambda_{2})\,,\qquad\dot{\lambda}_{2}=2(-\omega_{2}\lambda_{2}+\omega_{3}\lambda_{1})\,,\qquad\dot{\lambda}_{3}=\delta_{2}\lambda_{1}-\delta_{1}\lambda_{2}\,.$

(2)

We point out that parameters $\lambda_{1,2}$ satisfy the classical Hamilton equations

$\displaystyle\dot{p}=-2(\omega_{2}p+\omega_{3}q)-\delta_{2}\,,\qquad\dot{q}=2(\omega_{2}q+\omega_{1}p)+\delta_{1}\,.$

(3)

with $\delta_{1}=\delta_{2}=0$. To show this, one can see that the differential equations for $\lambda_{1,2}$ correspond to the classical equations with the substitution $\lambda_{1}\rightarrow q$ and $\lambda_{2}\rightarrow-p$; in other words, they correspond to the time inversion of the classical equations. In the case of the differential equation for $\lambda_{3}$, one can show, in view of the Hamilton equations, that it corresponds to the classical Lagrangian (with $\delta_{1,2}\neq 0$) plus the time variation of the function $pq$ of the system, that is,

$$\dot{\lambda}_{3}=-2\mathcal{L}+\dot{L}\,,$$

(4)

where $\mathcal{L}=p\dot{q}-H$ and $L=pq$. From these identifications, one can conclude that the classical dynamics given by the Hamilton equation or the equation of motion can lead to the solution of quantum dynamics given by the Hamiltonian operator of equation (2). For example, one can derive the propagator of the system $G(x,x^{\prime},t)=\langle x|\hat{U}(t)|x^{\prime}\rangle$ using each one of the solutions of the classical problem [35].

In this section, we review the equations determining the evolution of the covariance matrix and mean values $\langle p_{j}\rangle$ and $\langle q_{j}\rangle$ for an arbitrary system under the evolution of a quadratic Hamiltonian; also we mention the connection and dynamics of the continuous density matrix parameters. To obtain these properties, we use, as in the one-dimensional case, the invariant operators defined in [33, 34].

In the case of an $N$-dimensional quadratic system, the time evolution is characterized by the Hamiltonian

$$\hat{H}=\tilde{\mathbf{r}}\mathbf{B}(t)\mathbf{r}+\tilde{\boldsymbol{\Delta}}(t)\mathbf{r}\,,$$

(26)

where the tilde corresponds to the transposition operation, the vector $\tilde{\mathbf{r}}=\left(\hat{p}_{1},\hat{q}_{1},\hat{p}_{2},\hat{q}_{2},\ldots,\hat{p}_{N},\hat{q}_{N}\right)$ corresponds to the vector of quadrature operators. The time dependence of this Hamiltonian is contained in the matrices

$$\mathbf{B}(t)=\left(\begin{array}[]{cccc}\omega_{1,1}(t)&\omega_{1,2}(t)&\cdots&\omega_{1,2N}(t)\\ \omega_{1,2}(t)&\omega_{2,2}(t)&\cdots&\omega_{2,2N}(t)\\ \vdots&\vdots&\ddots&\vdots\\ \omega_{1,2N}(t)&\omega_{2,2N}(t)&\cdots&\omega_{2N,2N}(t)\end{array}\right),\qquad\boldsymbol{\Delta}(t)=\left(\begin{array}[]{cccc}\delta_{1}(t)\\ \delta_{2}(t)\\ \vdots\\ \delta_{2N}(t)\end{array}\right)\,,$$

(27)

where $\mathbf{B}(t)$ is a real and symmetric matrix and $\boldsymbol{\Delta}(t)$ is a real vector. As in the one-dimensional case, there exist $2N$ linear time-dependent operators $\hat{R}_{j}$ ($j=1,\ldots,2N$), whose time derivatives are equal to zero $\left(\dfrac{d\hat{R}_{j}}{dt}=0.\right)$. These operators can be arranged on a vector as follows:

$$\boldsymbol{\mathbf{R}}=\boldsymbol{\Lambda}(t)\mathbf{r}+\boldsymbol{\Gamma}(t)\,,$$

(28)

with the matrix $\boldsymbol{\Lambda}(t)$ and the vector $\tilde{\boldsymbol{\Gamma}}(t)=\left(\gamma_{1}(t),\gamma_{2}(t),\ldots,\gamma_{2N}(t)\right)$. By taking the time-derivative of the operator $\boldsymbol{\mathbf{R}}_{j}=\hat{R}_{j}$ and equating it to zero $\left(\dfrac{d\hat{R}_{j}}{dt}=0.\right)$, one can demonstrate that $\boldsymbol{\Lambda}(t)$ and $\boldsymbol{\Gamma}(t)$ must satisfy the following differential equations:

$$\dot{\boldsymbol{\Lambda}}(t)=2i\boldsymbol{\Lambda}(t)\mathbf{D}\mathbf{B}(t),\quad\dot{\boldsymbol{\Gamma}}(t)=i\boldsymbol{\Lambda}(t)\mathbf{D}\boldsymbol{\Delta}(t),\quad{\rm with}\quad\mathbf{D}_{j,k}=[\mathbf{r}_{j},\mathbf{r}_{k}]\,.$$

(29)

The solution to these differential equations with the initial conditions $\hat{R}_{j}(0)=\mathbf{r}_{j}$ provides, as a result, the invariant operators $\tilde{\boldsymbol{R}}=\left(\hat{P}_{1},\hat{Q}_{1},\hat{P}_{2},\hat{Q}_{2},\ldots,\hat{P}_{N},\hat{Q}_{N}\right)$, satisfying the standard commutation rules for the operators at zero time: $\mathbf{r}$, i.e., $[\mathbf{R}_{j},\mathbf{R}_{k}]=[\mathbf{r}_{j},\mathbf{r}_{k}]=\mathbf{D}_{j,k}$. This property leads us to the conclusion that the matrix $\boldsymbol{\Lambda}(t)$ must be symplectic and satisfies the equation

$$\boldsymbol{\Lambda}(t)\mathbf{D}\tilde{\boldsymbol{\Lambda}}(t)=\mathbf{D},$$

(30)

with $\mathbf{D}_{j,k}=[\mathbf{r}_{j},\mathbf{r}_{k}]$; this relation can be then used to obtain the inverse of $\boldsymbol{\Lambda}(t)$, which results in the expression

$$\boldsymbol{\Lambda}^{-1}(t)=\mathbf{D}\tilde{\boldsymbol{\Lambda}}(t)\mathbf{D}\,.$$

(31)

The other important property of these invariant operators is that they correspond to the inverse evolution of the original operators, in other words,

$$\mathbf{R}_{j}=\hat{U}(t)\mathbf{r}_{j}\hat{U}^{\dagger}(t),$$

which in the most cases can be obtained from the Heisenberg picture operators by assuming a time reversal operation. This property implies that the entries of $\dot{\boldsymbol{\Lambda}}(t)$ in equation (29) satisfy the classical Hamilton equations after the time reversal operation, that is, after the change $p_{i}\rightarrow-p_{i}$ in the classical Hamilton equations.

By the use of these invariant operators, one can obtain the time dependence of the mean values of the operators in $\mathbf{r}$ ($\langle\mathbf{r}\rangle(t)$) and their covariances $\boldsymbol{\sigma}_{j,k}=\dfrac{1}{2}\langle\{\mathbf{r}_{j},\mathbf{r}_{k}\}\rangle-\langle\mathbf{r}_{j}\rangle\langle\mathbf{r}_{k}\rangle$. From the inverse of equation (28), one can demonstrate that

$$\langle\mathbf{r}\rangle(t)=\boldsymbol{\Lambda}^{-1}(t)(\langle\mathbf{R}\rangle-\boldsymbol{\Gamma}(t))\,,$$

as the invariant operators in $\mathbf{R}$ have a time derivative equal to zero, and they are equal to the standard operators $\mathbf{r}$ at zero time, then one can conclude that

$$\langle\mathbf{r}\rangle(t)=\boldsymbol{\Lambda}^{-1}(t)(\langle\mathbf{r}\rangle(0)-\boldsymbol{\Gamma}(t))\,.$$

(32)

From an analogous argument, one can see that the covariance matrix reads

$$\boldsymbol{\sigma}(t)=\boldsymbol{\Lambda}^{-1}(t)\boldsymbol{\sigma}(0)\tilde{\boldsymbol{\Lambda}}^{-1}(t)\,.$$

(33)

Then to obtain the expression for the time-derivative of the mean values $\langle\mathbf{r}\rangle(t)$ and the covariance matrix $\boldsymbol{\sigma}(t)$, we make use of equations (29), (31), (32), and (33) and arrive at the expressions

$$\frac{d}{dt}\langle\mathbf{r}\rangle(t)=-i\mathbf{D}(2\mathbf{B}(t)\langle\mathbf{r}\rangle(t)+\boldsymbol{\Delta}(t))\,,\qquad\frac{d}{dt}\boldsymbol{\sigma}(t)=2i(\boldsymbol{\sigma}(t)\mathbf{B}(t)\mathbf{D}-\mathbf{D}\mathbf{B}(t)\boldsymbol{\sigma}(t))\,.$$

(34)

These differential equations, being first obtained in [31, 32] are the generalization of the one-dimensional case discussed in the previous section. In our case, the nonlinear differential equations for the density matrix parameters can be obtained by explicit calculation of the covariances at time $t$. The resulting equations can then be solved by the substitution of the solution of equation (34) or by direct integration.

To make the relation easier to see, we point out that the $2N\times 2N$ symplectic matrix $\mathbf{D}$ contains in its diagonal, blocks made of the $2\times 2$ symplectic matrix $\boldsymbol{\Sigma}$, that is,

$$\mathbf{D}=-i\left(\begin{array}[]{ccccc}\boldsymbol{\Sigma}&0&0&\cdots&0\\ 0&\boldsymbol{\Sigma}&0&\cdots&0\\ \vdots&\ddots&\ddots&\ddots&\vdots\\ 0&\cdots&0&\boldsymbol{\Sigma}&0\\ 0&0&\cdots&0&\boldsymbol{\Sigma}\end{array}\right)\,.$$

One can also notice that the differential equations for the entries of the covariance matrix are linear and can be expressed in the following matrix form:

$$\frac{d}{dt}\mathbf{v}=\mathbf{M}\mathbf{v}\,,$$

(35)

where $\mathbf{v}$ is a $N(2N+1)$-dimensional vector, which contains all the independent covariances in the $N$-partite system, and the matrix $\mathbf{M}$ is a square matrix of the same dimension that contains the Hamiltonian coefficients.