H-Theorem in an Isolated Quantum Harmonic Oscillator

Abstract

We consider the H-theorem in an isolated quantum harmonic oscillator through the time-dependent Schrödinger equation. The effect of potential in producing entropy is investigated in detail, and we found that including a barrier potential into a harmonic trap would lead to the thermalization of the system, while a harmonic trap alone would not thermalize the system. During thermalization, Shannon entropy increases, which shows that a microscopic quantum system still obeys the macroscopic thermodynamics law. Meanwhile, initial coherent mechanical energy transforms to incoherent thermal energy during thermalization, which exhibiting the decoherence of an oscillating wave packet featured by a large decreasing of autocorrelation length. When reaching thermal equilibrium, the wave packet comes to a halt, with the density distributions both in position and momentum spaces well-fitted by a microcanonical ensemble of statistical mechanics.

1. Introduction

H-theorem plays an important role in the time evolution of system in both the classical and quantum statistical mechanics. In particular, it identifies the entropy function appearing in the second law of thermodynamics. It is a longstanding question on the microscopic description of the second law of thermodynamics. One fundamental issue is whether an isolated system can reach thermal equilibrium, i.e., the state with maximum entropy [1,2,3,4,5,6,7,8]. In isolated quantum systems, it is especially important to see how the reversible microscopic quantum mechanics conceals the irreversible macroscopic phenomena of thermodynamics. Therefore, a directional time evolution of isolated quantum systems at nonequilibrium would often reveal the thermalization process concealed in quantum mechanics. Fortunately, ultracold quantum gases, which are pure and controllable, provide an excellent platform to study the nonequilibrium dynamics for isolated quantum systems. In atomic Bose–Einstein condensate (BEC) experiments, Kinoshita et al. [9] showed no evidence of thermalization by pairwise collision, from the Tonks–Girardeau limit to the intermediate coupling regime. Nevertheless, the dissipative motion of oscillating BEC in a disorder trap, first studied by Dries et al. [10] and further investigated by Hsueh et al. [11], does manifest the thermalization of an isolated quantum system. In particular, the noninteracting-limit case has been numerically studied in [12], which depicts Shannon entropy increasing to its maximum during thermalization, and at equilibrium, the initial mechanical energy is transformed to the thermal energy consisting of evenly distributed kinetic and potential energies.

In equilibrium statistical mechanics, one primary task is to determine the distribution function. Owing to the dependence of interaction, exact distribution function of a many-body system is usually not accessible. In nonequilibrium statistical mechanics, in contrast, the most difficult task is to determine the governing factors that describe the thermalization. It will be shown that the issue is not related to whether there is another new law describing the thermalization, but on how the governing microscopic laws can properly describe it and at the same time obey macroscopic thermodynamic laws. In this paper, we try to answer these two questions by studying an isolated quantum harmonic oscillator (QHO) that thermalizes through a single barrier. The disorder trap used for thermalization in [12] is constituted by a harmonic trap superposed with a Gaussian correlated disorder potential, which consists of many bumps in various sizes. Here, we demonstrate that a harmonic trap superposed with a single bump barrier alone can thermalize QHO as well through wave scattering. In addition, in [12], thermalization starts from a mobilized ground state of harmonic trap, while in this paper, besides a ground state of harmonic trap, we additionally employed a mixture of two pure states (such as ground and first excited states) of harmonic trap as the initial condition in order to exhibit different entropy evolution behaviors during thermalization. Moreover, here we presented a rigorous theoretical derivation of entropy production that shows that additional barrier potential can cause the increasing of entropy and leading to thermalization, while harmonic trap alone would not.

To be more specific, we studied this thermalization process by the time evolution of density distribution, for which we have numerically solved the time-dependent Schrödinger equation (TDSE), Equation (8), with initial condition, Equation (17), with $\hslash \omega$ and the harmonic trapping length $l=\sqrt{\hslash /m\omega }$ taken as the units of energy and length. To numerically integrate Equation (8), we use the method of lines with spatial discretization by a Fourier pseudospectral method and time integration by an RK45 scheme [13].

In this study, we will show that, after being released from the proposed initial coherent state with an excessive potential energy by placing the wave packet away from the center of trap initially, the wave packet would oscillate in this trap, wave-scattered by the single bump barrier located in the center of trap, and at the same time gradually transform the initial excessive potential energy into incoherent thermal energy. The whole phenomenon is analogous to a classical damped harmonic oscillator with the single bump barrier acting as the friction. When evolution time is long enough ($t\gg t_{\mathrm{th}}$ with $t_{\mathrm{th}}$ being the thermalization time), the system approaches the equilibrium. During the relaxation to its equilibrium, we found that an exchange between kinetic energy and potential energy continues to happen but with diminishing flux. In addition, we observed that the thermalization process of such a microscopic quantum system would still obey the second law of thermodynamics by the examination of increasing Shannon entropy. When equilibrium is reached, the oscillation stops, kinetic-potential energy exchange comes to a halt, Shannon entropy reaches its maximum, and the initial excessive potential energy is totally dissipated into incoherent thermal energy.

The paper is organized as follows: in Section 2, Shannon entropy in momentum space is adopted here to describe the thermalization of the current isolated quantum system, and an entropy production equation in momentum space is derived to analyze the effect of external potential. In Section 3, we found that harmonic trap alone will not generate any entropy, but with a barrier potential; in addition, entropy can then be produced, which means that the recruitment of a barrier potential into the harmonic trap would cause thermalization. In addition, the entropy production during thermalization is accompanied by a decoherence of wave packet featured by a large decreasing of autocorrelation length. In Section 4, when reaching thermal equilibrium eventually, we show that the density distributions both in position and momentum spaces are well fitted by a microcanonical ensemble of statistical mechanics. In Section 5, we study the time evolution of energy and investigate how different, from a pure QHO kinetic and potential, energies exchange with each other during thermalization, which describes how a coherent mechanical energy transforms to an incoherent thermal energy. Finally, we present the conclusions in Section 6.

2. Model and Methods

In quantum statistics, one often follows the von Neumann entropy:

$$S_v\left(t\right)=-k_{\mathrm{B}}\mathrm{Tr}\left[\widehat{\rho }\left(t\right)ln\widehat{\rho }\left(t\right)\right],$$

(1)

with $\widehat{\rho }$ being the density operator. The production rate of $S_v$ is then

$$\frac{d{S}_v\left(t\right)}{dt}=-k_{\mathrm{B}}\mathrm{Tr}\left[\frac{\partial \widehat{\rho }\left(t\right)}{\partial t}ln\widehat{\rho }\left(t\right)+\frac{\partial \widehat{\rho }\left(t\right)}{\partial t}\right].$$

(2)

The second term in Equation (2) vanishes due to the normalization of the total number of states. It is a general concept that $S_v\left(t\right)$ as well as $d{S}_v\left(t\right)/dt$ are basis independent. In terms of an arbitrary basis $|i\rangle$,

$$\frac{d{S}_v\left(t\right)}{dt}=-k_{\mathrm{B}}\sum _{i,j}\langle i|\frac{\partial \widehat{\rho }}{\partial t}|j\rangle \langle j|ln\widehat{\rho }|i\rangle .$$

(3)

However, if the basis $|i\rangle$ is chosen to be eigenkets of $\widehat{\rho }$ or energy eigenkets, then

$$\frac{d{S}_v\left(t\right)}{dt}=-k_{\mathrm{B}}\sum _{i,j}\langle i|\frac{\partial \widehat{\rho }}{\partial t}|j\rangle \langle j|ln\widehat{\rho }|i\rangle =\frac{i{k}_{\mathrm{B}}}{\hslash }\sum _{i,j}\langle i|[\widehat{H},\widehat{\rho }]|j\rangle \langle j|ln\widehat{\rho }|i\rangle =0.$$

(4)

Here, $\widehat{H}$ is the Hamiltonian. The vanishing of $d{S}_v\left(t\right)/dt$ in Equation (4) indicates that von Neumann entropy can not offer any entropy production for an isolated quantum system going through thermalization [14,15].

Alternatively, one can consider the Shannon entropy ($S_Q$) whose production rate is defined as [16,17,18]

$$\frac{d{S}_Q\left(t\right)}{dt}=-k_{\mathrm{B}}\sum _i\langle i|\frac{\partial \widehat{\rho }}{\partial t}|i\rangle ln\langle i|\widehat{\rho }|i\rangle =\frac{i{k}_{\mathrm{B}}}{\hslash }\sum _i\langle i|[\widehat{H},\widehat{\rho }]|i\rangle \langle i|ln\widehat{\rho }|i\rangle ,$$

(5)

where the basis $|i\rangle$ is the eigenkets of operator $\widehat{Q}$, which is arbitrary except for using the eigenkets of either $\widehat{\rho }$ or $\widehat{H}$. In addition, the basis $|i\rangle$ is not suitable to be chosen as the eigenkets of other complete sets of commuting observables (CSCO). In general, position $\widehat{x}$ and momentum $\widehat{p}$ are suitable choices of $\widehat{Q}$ since they are not both CSCO. When the basis $|i\rangle$ is chosen to be momentum eigenkets, the corresponding Shannon entropy and its production rate are

$$S_p\left(t\right)=-\frac{k_{\mathrm{B}}}{2\pi }{\int }_{-\infty }^{\infty }{\rho }_{kk}\left(t\right)ln\left[\frac{{\rho }_{kk}\left(t\right)}{2\pi l}\right]dk,$$

(6)

and

$$\frac{d{S}_p\left(t\right)}{dt}=-\frac{k_{\mathrm{B}}}{2\pi }{\int }_{-\infty }^{\infty }\frac{\partial {\rho }_{kk}\left(t\right)}{\partial t}ln\left[\frac{{\rho }_{kk}\left(t\right)}{2\pi l}\right]dk,$$

(7)

where ${\rho }_{kk}\equiv {\left|\tilde{\psi }(k,t)\right|}^2$ is the momentum distribution with $\tilde{\psi }(k,t)=\langle k|\psi \left(t\right)\rangle$.

In this paper, we focus on a one-dimensional quantum harmonic oscillator, initially released from an off-center position, that corresponds to an out-of-equilibrium state [10,19,20,21]. The system is well described by the time-dependent Schrödinger equation (TDSE)

$$i\hslash \frac{\partial \psi \left(x,t\right)}{\partial t}=\left[-\frac{{\hslash }^2}{2m}\frac{{\partial }^2}{\partial x^2}+V\left(x\right)\right]\psi \left(x,t\right),$$

(8)

where m is the particle mass, and the external potential is $V\left(x\right)=V_0\left(x\right)+V_1\left(x\right)$ consisting of a confining harmonic well $V_0\left(x\right)=m{\omega }^2{x}^2/2$ with $\omega$ the trapping frequency and a barrier potential $V_1\left(x\right)$ described later.

The real-space TDSE, Equation (8), can be Fourier transformed to momentum space,

$$i\hslash \frac{\partial \tilde{\psi }(k,t)}{\partial t}=\frac{{\hslash }^2{k}^2}{2m}\tilde{\psi }(k,t)+\frac{1}{2\pi }{\int }_{-\infty }^{\infty }\tilde{V}\left(k-k_1\right)\tilde{\psi }(k_1,t)d{k}_1,$$

(9)

where $\tilde{\psi }(k,t)$ and $\tilde{V}\left(k\right)$ are respectively Fourier transforms of wave function $\psi (x,t)$ and external potential $V\left(x\right)$. The Fourier transform of harmonic well $V_0\left(x\right)$ is

$$\begin{array}{ccc}\hfill {\tilde{V}}_0\left(k\right)& =& {\int }_{-\infty }^{\infty }{V}_0\left(x\right)exp(-ikx)dx=\frac{m{\omega }^2}{2}{\int }_{-\infty }^{\infty }{x}^{2}exp(-ikx)dx\hfill \\ & =& \frac{m{\omega }^2}{2{(-i)}^2}{\int }_{-\infty }^{\infty }\frac{{\partial }^{2}exp(-ikx)}{\partial k^2}dx=-\frac{m{\omega }^2}{2}\frac{{\partial }^2}{\partial k^2}{\int }_{-\infty }^{\infty }exp(-ikx)dx\hfill \\ & =& -\pi m{\omega }^2\frac{{\partial }^2\delta \left(k\right)}{\partial k^2},\hfill \end{array}$$

(10)

where $\delta \left(k\right)=1/\left(2\pi \right){\int }_{-\infty }^{\infty }exp(-ikx)dx$ is the Dirac delta function. The complex-conjugate of Equation (9) is

$$-i\hslash \frac{\partial {\tilde{\psi }}^{*}(k,t)}{\partial t}=\frac{{\hslash }^2{k}^2}{2m}{\tilde{\psi }}^{*}(k,t)+\frac{1}{2\pi }{\int }_{-\infty }^{\infty }{\tilde{V}}^{*}\left(k-k_1\right){\tilde{\psi }}^{*}(k_1,t)d{k}_1.$$

(11)

By ${\tilde{\psi }}^{*}(k,t)\times$ Equation (9)$-\tilde{\psi }(k,t)\times$ Equation (11), we obtain

$$\begin{array}{ccc}\hfill i\hslash \frac{\partial {\rho }_{kk}\left(t\right)}{\partial t}& =& i\hslash \left[{\tilde{\psi }}^{*}(k,t)\frac{\partial \tilde{\psi }(k,t)}{\partial t}+\tilde{\psi }(k,t)\frac{\partial {\tilde{\psi }}^{*}(k,t)}{\partial t}\right]\hfill \\ \hfill & =& \frac{1}{2\pi }{\int }_{-\infty }^{\infty }\left[\tilde{V}\left(k-k_1\right){\tilde{\psi }}^{*}(k,t)\tilde{\psi }(k_1,t)-{\tilde{V}}^{*}\left(k-k_1\right)\tilde{\psi }(k,t){\tilde{\psi }}^{*}(k_1,t)\right]d{k}_1,\hfill \end{array}$$

(12)

describing the evolution of momentum distribution ${\rho }_{kk}\left(t\right)$. By $\tilde{V}={\tilde{V}}_0+{\tilde{V}}_1$, Equation (12) can be reorganized as

$$\begin{array}{ccc}\hfill \frac{\partial {\rho }_{kk}\left(t\right)}{\partial t}& =& -\frac{1}{\pi \hslash }{\int }_{-\infty }^{\infty }\Im \left[{\tilde{V}}^{*}\left(k-k_1\right){\rho }_{k{k}_1}\right]d{k}_1\hfill \\ & =& -\frac{1}{\pi \hslash }{\int }_{-\infty }^{\infty }\Im \left[{\tilde{V}}_0^{*}\left(k-k_1\right){\rho }_{k{k}_1}\right]d{k}_1-\frac{1}{\pi \hslash }{\int }_{-\infty }^{\infty }\Im \left[{\tilde{V}}_1^{*}\left(k-k_1\right){\rho }_{k{k}_1}\right]d{k}_1\hfill \\ & =& \frac{\omega }{l^2}{\int }_{-\infty }^{\infty }\Im \left[\frac{{\partial }^2\delta \left(k-k_1\right)}{\partial k_1^2}\tilde{\psi }(k,t){\tilde{\psi }}^{*}(k_1,t)\right]d{k}_1-\frac{1}{\pi \hslash }{\int }_{-\infty }^{\infty }\Im \left[{\tilde{V}}_1^{*}\left(k-k_1\right){\rho }_{k{k}_1}\right]d{k}_1\hfill \\ & =& \frac{\omega }{l^2}{\int }_{-\infty }^{\infty }\Im \left[\delta \left(k-k_1\right)\tilde{\psi }(k,t)\frac{{\partial }^2{\tilde{\psi }}^{*}(k_1,t)}{\partial k_1^2}\right]d{k}_1-\frac{1}{\pi \hslash }{\int }_{-\infty }^{\infty }\Im \left[{\tilde{V}}_1^{*}\left(k-k_1\right){\rho }_{k{k}_1}\right]d{k}_1\hfill \\ & =& \frac{\omega }{l^2}\Im \left[\tilde{\psi }(k,t)\frac{{\partial }^2{\tilde{\psi }}^{*}(k,t)}{\partial k^2}\right]-\frac{1}{\pi \hslash }{\int }_{-\infty }^{\infty }\Im \left[{\tilde{V}}_1^{*}\left(k-k_1\right){\rho }_{k{k}_1}\right]d{k}_1,\hfill \end{array}$$

(13)

where $l=\sqrt{\hslash /m\omega }$ (harmonic trapping length), ${\rho }_{kq}\equiv \tilde{\psi }(k,t){\tilde{\psi }}^{*}(q,t)$, and $\Im [A]$ denotes the imaginary part of A. The first term of the right-hand side of Equation (13), which comes from the harmonic potential, indicates the contribution of diagonal elements of the density matrix, while the second term, from the barrier potential $V_1$, indicates the contribution of off-diagonal elements, which could lead to the decoherence of a density matrix.

Substitution of Equation (13) into Equation (7) yields

$$\begin{array}{ccc}\hfill \frac{d{S}_p\left(t\right)}{dt}& =& \frac{k_{\mathrm{B}}}{2{\pi }^2\hslash }{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }\Im \left[{\tilde{V}}^{*}\left(k-k_1\right){\rho }_{k{k}_1}\right]ln\left(\frac{{\rho }_{kk}}{2\pi l}\right)dkd{k}_1\hfill \\ & =& -\frac{k_{\mathrm{B}}\omega }{2\pi l^2}{\int }_{-\infty }^{\infty }\Im \left[\tilde{\psi }(k,t)\frac{{\partial }^2{\tilde{\psi }}^{*}(k,t)}{\partial k^2}\right]ln\left(\frac{{\rho }_{kk}}{2\pi l}\right)dk\hfill \\ & +& \frac{k_{\mathrm{B}}}{2{\pi }^2\hslash }{\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }\Im \left[{\tilde{V}}_1^{*}\left(k-k_1\right){\rho }_{k{k}_1}\right]ln\left(\frac{{\rho }_{kk}}{2\pi l}\right)dkd{k}_1.\hfill \end{array}$$

(14)

It will be shown later that only the second term of the right-hand side of the equation above, coming from barrier potential $V_1$, would contribute to the increasing of entropy.

For QHO, the generalized coherent states (GCSs) $|n,\alpha \rangle$, with n denoting the quantum number and $\alpha$ the complex quantum number, have the following forms in position and momentum spaces [22]:

$$\begin{array}{c}\hfill {\phi }_{\alpha ,n}(x,t)=\langle x|\alpha ,n\rangle =\frac{1}{\sqrt{l}}{\varphi }_n\left(\frac{x-\overline{x}}{l}\right)exp(-in\omega t)exp\left[-i\left(\frac{\omega t}{2}-x\overline{k}+\frac{\overline{x}\overline{k}}{2}\right)\right],\end{array}$$

(15)

$$\begin{array}{c}\hfill {\tilde{\phi }}_{\alpha ,n}(k,t)=\langle k|\alpha ,n\rangle =i^n\sqrt{2\pi l}{\varphi }_n\left[l(k-\overline{k})\right]exp(-in\omega t)exp\left[-i\left(\frac{\omega t}{2}-\overline{x}k+\frac{\overline{x}\overline{k}}{2}\right)\right]\end{array}$$

(16)

with $\overline{x}\left(t\right)/l=\sqrt{2}\alpha cos\left(\omega t\right)$, $l\overline{k}\left(t\right)=-\sqrt{2}\alpha sin\left(\omega t\right)$, and the eigenfunction of quantum harmonic oscillator ${\varphi }_n\left(x\right)=1/\left({\pi }^{1/4}\sqrt{2^{n}n!}\right)H_n\left(x\right)exp(-x^2/2)$, which carries energy ${\epsilon }_{\alpha ,n}=({\alpha }^2+n+1/2)\hslash \omega$. $H_n\left(x\right)$ is Hermite polynomial. Note that the corresponding wave packet oscillates in the harmonic potential (with time period $2\pi /\omega$) without changing the shape of density.

The initial conditions in x- and k-space corresponding to a mixture of GCSs would have the forms

$$\psi (x,0)=\sum _n{c}_n{\phi }_{\alpha ,n}(x,0)=\frac{1}{\sqrt{l}}\sum _n{c}_n{\varphi }_n\left[\frac{x-\overline{x}\left(0\right)}{l}\right]$$

(17)

and

$$\tilde{\psi }(k,0)=\sum _n{c}_n{\tilde{\phi }}_{\alpha ,n}(k,0)=\sqrt{2\pi l}\sum _n{i}^n{c}_n{\varphi }_n\left[l(k-\overline{k}\left(0\right))\right]exp\left[i\left(\overline{x}\left(0\right)k-\frac{\overline{x}\left(0\right)\overline{k}\left(0\right)}{2}\right)\right],$$

(18)

where $c_n$ is the coefficient of superposition. $\overline{x}\left(0\right)=\sqrt{2}\alpha l$ and $\overline{k}\left(0\right)=0$ are the initial centroid locations of wave packet in position and momentum spaces, respectively. Without $V_1$, for an initial condition with superposition of two CGSs at most, $\tilde{\psi }(k,t)=c_n{\tilde{\phi }}_{\alpha ,n}(k,t)+c_{n^{\prime }}{\tilde{\phi }}_{\alpha ,n^{\prime }}(k,t)$, ${\rho }_{kk}$ is read as

$$\begin{array}{ccc}\hfill \frac{{\rho }_{kk}\left(t\right)}{2\pi l}& =& |c_n{|}^2{\varphi }_n^2\left[l(k-\overline{k})\right]+{\left|c_{n^{\prime }}\right|}^2{\varphi }_{n^{\prime }}^2\left[l(k-\overline{k})\right]\hfill \\ & +& 2\Re \left[i^{(n-n^{\prime })}{c}_n{c}_{n^{\prime }}^{*}exp\left[-i\left(n-n^{\prime }\right)\omega t\right]\right]{\varphi }_n\left[l(k-\overline{k})\right]{\varphi }_{n^{\prime }}\left[l(k-\overline{k})\right],\hfill \end{array}$$

(19)

which basically still oscillates in the harmonic trap with period $2\pi /\omega$, since $\overline{k}$ is a periodic function of time, but with the shape modulated by a period $2\pi /[(n-n^{\prime })\omega ]$ due to the factor $exp[-i(n-n^{\prime })\omega t]$. Here, $\Re [A]$ denotes the real part of A. For the case of single CGS, distribution Equation (19) reduces to ${\rho }_{kk}\left(t\right)=2\pi l{\varphi }_n^2\left[l(k-\overline{k})\right]$, which oscillates with a period $2\pi /\omega$ and without changing its shape, so that the corresponding entropy $S_p\left(t\right)$ is just a constant according to Equation (6) shown as the blue horizontal line in Figure 1. This can also be verified by the vanishing of the first term on the right-hand side of Equation (14) (see Appendix A for details). For the case of a mixture of GCS, the wave packet again oscillates inside the harmonic trap but with the shape of density modulated periodically. Therefore, its corresponding entropy is periodic in time as shown with the red line in Figure 1.

3. Barrier Potential Produces Entropy and Decoherence

As postulated earlier, a bump (barrier) potential may cause the wave scattering of an oscillating wave packet, which leads to its decoherence featured by the increasing of entropy. Here, we consider a single central Gaussian barrier, $V_1\left(x\right)=V_{\mathrm{G}}exp(-x^2/l^2)$, of which the Fourier transform is ${\tilde{V}}_1\left(k\right)=V_{\mathrm{G}}l\sqrt{\pi }exp(-k^2{l}^2/4)$. One can picture that, with this ${\tilde{V}}_1\left(k\right)$ recruited to Equation (14), it can lead to the production of entropy. Since integrating the second term into the right-hand side of Equation (14) by hand is infeasible, we can only verify this entropy production by integrating TDSE, Equation (8), numerically.

Figure 1 shows the time development of Shannon entropies for a quantum oscillator in various initial conditions. For $V=V_0$ only, blue and red lines correspond to the initial conditions for an off-centered single GCS (${\phi }_{\alpha ,0}$) and a mixture of two GCS’s (${\phi }_{\alpha ,0}$ and ${\phi }_{\alpha ,2}$), respectively. Here, we took $\alpha =15/\sqrt{2}$, which means that the initial centroid position is $\overline{x}\left(0\right)=15l$, a position far away from the bump barrier such that the wave packet would not be immediately affected by the bump barrier right at the beginning of oscillation. The total energies of these two cases are $E=({15}^2+1)/2\hslash \omega$ and $E=({15}^2+3)/2\hslash \omega$, respectively. For $V=V_0+V_1$, with the recruitment of a Gaussian barrier $V_1$ having strength $V_{\mathrm{G}}=100\hslash \omega$, green and black lines correspond to the initial conditions for an off-centered single GCS and a mixture of two GCS’s as above. With the increasing trend of entropy, shown by green and black lines, it verifies the entropy production owing to the inclusion of barrier potential. Note that the fluctuations of green and black lines come from the modulation of the wave packet when oscillating within the trap. As time goes on, the system reaches thermal equilibrium and the entropies of green and black lines would grow to become saturated with diminishing fluctuations. The shape of wave packet, or density distribution, at equilibrium will be discussed in the next section.

Besides producing entropy, thermalization also causes the decoherence of a wave packet. Here, we define a spatial autocorrelation function (SACF) for wavefunction in order to measure the autocorrelation length, which can serve as a good indicator of coherence in a quantum system:

$$f(x,t)=\left|{\int }_{-\infty }^{\infty }{\psi }^{*}(x-x^{\prime },t)\psi (x^{\prime },t)d{x}^{\prime }\right|=\left|\frac{1}{2\pi }{\int }_{-\infty }^{\infty }{\rho }_{kk}\left(t\right)exp\left(ikx\right)dk\right|.$$

(20)

$f(x,t)$ depends intimately on the behaviors of $\psi (x,t)$ and measures the correlation between points separated by x at a given time t. Figure 2 shows $f(x,t)$ at an initial time and equilibrium with the blue curve corresponding to the initial one and the red curve corresponding to the equilibrium one.

At the beginning of oscillation, $f(x,t=0)$ exhibits a close-to-unity constant distribution in space, featuring a large autocorrelation length indicating that the system is highly coherent. In contrast, when the system comes to equilibrium, $f(x,t\gg t_{\mathrm{th}})$ is largely narrowed, featuring a small autocorrelation length indicating the strong incoherence of system. $t_{\mathrm{th}}$ is the so-called thermalization time, after which the system reaches equilibrium [12]. This significant decreasing of autocorrelation length reveals the decoherence of system during thermalization.

4. Density Distribution When Reaching Thermal Equilibrium

Figure 3 shows density distributions at initial time (green line) and equilibrium (blue line), with panels (a) and (c) in position space, and panels (b) and (d) in momentum space. One can see from panels (a) and (c), where the initial narrowly distributed $\rho \left(x\right)$, centered at $\overline{x}\left(0\right)=15l$, evolves to a centered-at-zero, symmetric, and much broader distribution at equilibrium. Similarity for the evolution of ${\rho }_{kk}$ is shown in panels (b) and (d), with the initial distribution centered at 0. These distributions are very close to the microcanonical distribution for a classical harmonic oscillator. In a microcanonical ensemble of a classical harmonic oscillator with a given energy E, the probability distribution in phase space is

$$P_{\mathrm{c}}(x,k)=\hslash \omega \delta \left(E-V\left(x\right)-\frac{{\hslash }^2{k}^2}{2m}\right)$$

(21)

with $1/\left(2\pi \right){\int }_{-\infty }^{\infty }{P}_{\mathrm{c}}(x,k)dxdk=1$, where $\delta (\cdots )$ is the Dirac delta function. The red lines in Figure 3 show the best fitting in position space, ${\rho }_{\mathrm{c}}\left(x\right)$, and momentum space, ${\tilde{\rho }}_{\mathrm{c}}\left(k\right)$:

$${\rho }_{\mathrm{c}}\left(x\right)=\frac{1}{2\pi }{\int }_{-\infty }^{\infty }{P}_{\mathrm{c}}(x,k)dk=\frac{m\omega }{\pi \sqrt{2m\left[E-V\left(x\right)\right]}}\mathrm{for}-x_{\mathrm{t}}<x<x_{\mathrm{t}},$$

(22)

and

$${\tilde{\rho }}_{\mathrm{c}}\left(k\right)={\int }_{-\infty }^{\infty }{P}_{\mathrm{c}}(x,k)dx\simeq \frac{2\hslash }{\sqrt{2mE-{\hslash }^2{k}^2}}\mathrm{for}-k_{\mathrm{t}}<k<k_{\mathrm{t}},$$

(23)

which are obtained from the microcanonical ensemble of statistical mechanics and are derived in Appendix B. Note that Equation (23) is approximated by neglecting the barrier potential $V_1\left(x\right)=V_{\mathrm{G}}exp(-x^2/l^2)$ and letting $V\left(x\right)=V_0\left(x\right)=m{\omega }^2{x}^2/2$. This is because, if the energy of system, E, is large enough, the effect of Gaussian barrier on phase-space distribution Equation (21) would be weak, and therefore the consequent distribution in momentum space can be approximated by ignoring it.

5. Energy Transformation

To monitor the evolution of kinetic and potential energies of the wave packet during its oscillation within the trap and the energy exchange between these two energy components, it is useful to consider the free energy of the system: $F=E-E_0$, where E is the total energy. Since total energy is conserved all the time, $E\left(t\right)=E\left(0\right)$, with $E\left(0\right)$ corresponding to the initial total energy of a wave packet, which is an GCS (or a mixture of GCSs) of the harmonic trap placed off-center initially. $E_0$ corresponds to the total energy of the initial wave packet located right at the center of trap, which would be stationary if the barrier potential $V_1\left(x\right)$ is neglected. Here, $E=K+U$ and similarly $E_0=K_0+U_0$ with K ($K_0$) and U ($U_0$) corresponding to kinetic and potential parts of energy. Basically, $K\left(0\right)=K_0$ since the shape of wave packet is the same, and $U\left(0\right)>U_0$ due to the initial off-center displacement of the wave packet. Therefore, $F\left(0\right)=E\left(0\right)-E_0=U\left(0\right)-U_0$, which means that the initial free energy is the initial excessive potential energy. When there is no barrier potential, the wave packet would keep its coherence with free energy $F\left(t\right)$ conserved and just re-distributed into kinetic and potential energies dynamically. It means that wave packet oscillates within the trap without thermalization as indicated by its non-increasing Shannon entropy depicted by the blue and red lines in Figure 1. Meanwhile, both $K\left(t\right)-K_0$ and $U\left(t\right)-U_0$ are periodic in time with their sum, $F\left(t\right)$, fixed in time, which means kinetic and potential energies fully exchange with each other during oscillation similar to a harmonic oscillator free of friction.

When barrier potential $V_1\left(x\right)$ exists, the oscillating wave packet is thermalized as indicated by its increasing Shannon entropy depicted by the green and black lines in Figure 1. Meanwhile, both $K\left(t\right)-K_0$ and $U\left(t\right)-U_0$ still oscillate with time but with diminishing amplitudes as shown in Figure 4. Note that their sum, $F\left(t\right)$, is still conserved in time. This means kinetic and potential energies exchange with each other less and less during thermalization, by which we can picture the wave packet oscillating less and less within the trap due to the wave scattering by the barrier. In this sense, the barrier acts like a friction to the oscillating wave packet that transforms its initial coherent mechanical energy into an incoherent thermal one. The decoherence of wave packet during thermalization can also be observed by its density distribution getting modulated by high-frequency fluctuations, besides a large decreasing of autocorrelation length as described in Figure 2. Note that the initial excessive coherent potential energy would not be equipartitioned into incoherent $K\left(t\right)-K_0$ and $U\left(t\right)-U_0$ as shown in Figure 4 when approaching thermal equilibrium. This is because the trap is no longer a harmonic one. We can picture that, when the height of barrier is smaller, energy would be more equipartitioned and sure it would also take a longer time to reach equilibrium due to a weaker barrier.

6. Conclusions

In conclusion, the thermalization in an isolated oscillating BEC in a harmonic trap superposed with a Gaussian barrier is studied. We solve the dynamic wavefunction $\psi (x,t)$ of the system and, from it, the relevant physics are obtained. The key results are: (i) The current microscopic quantum system is found to obey the macroscopic second law of thermodynamics by the increasing of Shannon entropy during thermalization, (ii) the production of entropy during thermalization is accompanied by the decoherence of system featured by the significant reduction of autocorrelation length, and (iii) the Gaussian barrier potential plays the role of dissipation with energy transformed from a coherent mechanical energy into an incoherent thermal one. (iv) At equilibrium, the density distributions both in position and momentum spaces are well fitted by a microcanonical ensemble of statistical mechanics.