Dynamics of a Quantum Dissipative System Coupled with an Oscillator

Emmanouil George Thrapsaniotis

Athens, Greece.

DOI: 10.4236/jamp.2024.124091                

Abstract

We study the dynamics of a quantum dissipative system. Besides its linear coupling with a harmonic bath modelling the dissipation, we suppose that it is coupled with an oscillator with an interaction of the form $s^2{x}^2$ . In our study, we integrate over the bath and the oscillator, extract the corresponding influence functionals and then solve the system’s sign problem. We apply the theory to the case of a double well and study the time evolution of the expectation value of the position.

Keywords

Double Well, Harmonic Bath, Oscillator, Influence Functional, Sign Solved Propagator, Position

1. Introduction

Path integral methods constitute an interesting and extended part of mathematical physics and there is considerable effort in their development [1] . The use of the central limit theorem in path integral methods is of interest as it leads to the solution of the sign problem [2] - [12] appearing in quantum physics. That solution is applicable to various systems even beyond quantum mechanics (see the conclusions in [6] ).

Here, we study the dissipative dynamics of a quantum mechanical system, which is coupled with an oscillator via an interaction term of the form $s^2{x}^2$ , where s is the coordinate of the system and x is one of the oscillators. The dissipation on the system is modelled via coupling the system’s particle with a harmonic heat bath of inverse temperature $\beta$ . Proceeding we first path integrate over the bath and the oscillator and obtain a path integral expression for the reduced system’s density matrix, which includes the bath and oscillator’s influence functionals. This kind of expression is well known to show a highly oscillatory phase leading to failure of numerical methods for their evaluation, namely the Monte Carlo method. That problem is known as the sign problem. We can solve that problem via extracting an alternative expression for the propagator called sign solved propagator where the oscillations are controlled. The whole method based on the use of the central limit theorem was developed by the author and was applied to other models in previous papers [2] - [12] . So, in that way, we can derive the time evolution of the system’s density matrix. In fact, here we study the time evolution of the position of a particle in a symmetric double well. We have chosen a double well potential as it incorporates tunnelling effects in the whole dynamics and behaves as a two-level system. In the final applications, we consider the particle interacting with only the harmonic bath, with only the oscillator or with both of them. Systems similar to the present one have been studied exhaustively. See for example [13] [14] [15] and references there. Other methods of study include the use of generalized Langevin equations or master equations [13] . However, the present path integral approach combined with the solution of the sign problem gives exact closed results from a fully quantum mechanical point of view.

The present paper proceeds as follows. In Section 2, we give the system and its Hamiltonian, consider the path integral that describes it and further path integrated over the bath and give the form of the corresponding influence functional. In Section 3, we derive the influence function of the interaction of the system with an oscillator. In Section 4, we solve the sign problem to study the time evolution of the system’s density matrix. In Section 5, we give results of the theory in the case of a double well coupled with a bath, or with an oscillator or with both of them and suppose initially a Gaussian wavefunction. In Section 6, we present our conclusions and finally, in Appendix, we solve the sign problem in the case of a density matrix.

2. Model Description

Since in the present paper, we consider the dynamics of a particle coupled on the one hand with a harmonic bath, which models a dissipative environment and on the other with a harmonic oscillator, the full Hamiltonian has the form:

$H_{tot}=H_s+H_{osc}+H_b$ (1)

$H_s$ is the system’s Hamiltonian given by:

$H_s=\frac{p_s^2}{2M}+V\left(s\right)$ (2)

The oscillator’s one has the form:

$H_{osc}=\frac{p^2}{2m}+\frac{1}{2}m{\omega }_0^2{x}^2+\frac{1}{2}\lambda s^2{x}^2$ (3)

supposing time independent parameters and coupling of the form $s^2{x}^2$ . That coupling can appear if we assume that the system deforms the oscillator’s potential. Then, we obtain the present effective interaction.

Finally, the harmonic bath’s Hamiltonian has the form:

$H_b={\displaystyle \underset{j}{\sum }\left[\frac{p_j^2}{2m_j}+\frac{1}{2}{m}_j{\omega }_j^2{\left(x_j-\frac{c_{j}s}{m_j{\omega }_j^2}\right)}^2\right]}$ (4)

which has a linear coupling. s may be interpreted as a reaction coordinate coupled to a large number of harmonic bath degrees of freedom. Equation (4) includes counterterms quadratic in s, which renormalize the system potential. That ensures that important potential features such as the barrier height do not change with the coupling strength.

Quantum mechanical observables of the system can be obtaining after tracing the full density matrix $W\left(t\right)$ over the bath and the oscillator. i.e.

$W_{red}\left(s^{+},s^{-};t\right)=T{r}_{bath}T{r}_x\langle s^{+}|W\left(t\right)|s^{-}\rangle =T{r}_{bath}T{r}_x\langle s^{+}|{\text{e}}^{-i{\displaystyle {\int }_0^t{H}_{tot}\text{d}\tau }}W\left(0\right){\text{e}}^{i{\displaystyle {\int }_0^t{H}_{tot}\text{d}\tau }}|s^{-}\rangle$ (5)

Moreover, we assume that the interaction of all the three subsystems is switched on at time $t=0$ . i.e. we assume that the density matrix at the initial time is:

$W\left(0\right)=W_s\left(0\right)W_{bath}\left(0\right)W_x\left(0\right)$ (6)

where $W_s\left(0\right)$ is the system’s initial density matrix, $W_{bath}\left(0\right)$ is the bath’s one and $W_x\left(0\right)$ the oscillator’s initial one.

The whole dynamics can be extracted via path integrating over the bath, the oscillator and the system. At first, we consider the integration over the bath. Its coordinates appear in the Hamiltonian (4). We assume it to be at inverse temperature $\beta$ . Then, according to standard methods [15] , we can obtain a corresponding influence functional in the form:

$\begin{array}{c}I=\exp \{-{\displaystyle \underset{0}{\overset{t}{\int }}\text{d}{t}^{\prime }}{\displaystyle \underset{0}{\overset{t^{\prime }}{\int }}\text{d}{t}^{″}}\left[s\left(t^{\prime }\right)-s^{\prime }\left(t^{\prime }\right)\right]\left[\alpha \left(t^{\prime }-t^{″}\right)s\left(t^{″}\right)-{\alpha }^{\ast }\left(t^{\prime }-t^{″}\right)s^{\prime }\left(t^{″}\right)\right]\\ -i{\displaystyle \underset{0}{\overset{t}{\int }}\text{d}{t}^{\prime }}{\displaystyle \underset{j}{\sum }\frac{c_j^2}{2m_j{\omega }_j^2}\left[s{\left(t^{\prime }\right)}^2-s^{\prime }{\left(t^{\prime }\right)}^2\right]}\}\end{array}$ (7)

where we set:

$\begin{array}{c}\alpha \left(t\right)={\displaystyle \underset{j}{\sum }\frac{c_j^2}{2m_j{\omega }_j}\left[\coth \left(\frac{{\omega }_j\beta }{2}\right)\cos \left({\omega }_{j}t\right)-i\sin \left({\omega }_{j}t\right)\right]}\\ =\frac{1}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}\text{d}\omega J\left(\omega \right)}\left[\coth \left(\frac{\omega \beta }{2}\right)\cos \left(\omega t\right)-i\sin \left(\omega t\right)\right]\end{array}$ (8)

where $J\left(\omega \right)$ is the spectral density. It incorporates the characteristics of the bath pertaining to the dynamics of the reaction coordinate corresponding to the system.

In the next section, we path integrate the oscillator Hamiltonian.

3. Integration over the Oscillator

We consider the Hamiltonian (3) of a harmonic oscillator of time independent frequencies and masses. We proceed to path integrate the Hamiltonian (3) via standard methods. We set:

$x=\frac{Q}{\sqrt{m}}$ (9)

$p=\sqrt{m}P$ (10)

to obtain the free of mass terms Hamiltonian:

$H_0\left(t\right)=\frac{P^2}{2}+\frac{1}{2}{\omega }_0^2{Q}^2+\frac{1}{2}\Gamma s^2{Q}^2$ (11)

where we set:

$\Gamma =\frac{\lambda }{m}$ (12)

If $K\left(x_f,x_i,t;s\right)$ is the propagator corresponding to the Hamiltonian (3) and $K_0\left(Q_f,Q_i,t;s\right)$ is the one corresponding to the Hamiltonian (11), then we can easily check that they obey the relation:

$K\left(x_f,x_i,t;s\right)=\sqrt{m}{K}_0\left(Q_f,Q_i,t;s\right)$ (13)

Therefore, we can concentrate our attention on the propagator of the Hamiltonian (11). It can be calculated via standard path integration. To proceed towards the integration, we perform the canonical transformation:

$Q=X\rho \left(s,t\right)$ (14)

$P=\frac{P_1}{\rho \left(s,t\right)}$ (15)

$\frac{\text{d}\tau }{\text{d}t}={\rho }^{-2}\left(s,t\right)$ (16)

The $\rho$ function depends on s through the differential Equation (32) (see below). i.e. s is supposed to be a function of time and to describe the coordinate of the double well. Moreover, during the present section’s evaluation, it is fixed. The above transformation is canonical since it preserves the Poisson brackets and therefore it preserves the volume element in phase space. As the present, transformation involves the generic time redefinition (16), we give more details. The $N+1$ time slices discrete form of $K_0\left(Q_f,Q_i,t;s\right)$ involves the times $t_n=n\epsilon$ $n=0,1,\cdots ,N+1$ , where the time step is $\epsilon =\frac{t}{N+1}$ . Now, on integrating the path integral expression on the momentums, it becomes (see below as well):

$K_0\left(Q_f,Q_i,t;s\right)={\left(\frac{1}{2\pi i\epsilon }\right)}^{\frac{N+1}{2}}{\displaystyle \int {\displaystyle \underset{n=1}{\overset{N}{\prod }}\left[\text{d}{Q}_n\right]}}\exp \left[i{S}^{\left(N\right)}\right]$ (17)

Then, under the transformations (14)-(16), the time step becomes ${\sigma }_n=\frac{\epsilon }{\rho \left(s,t_n\right)\rho \left(s,t_{n-1}\right)}$ , where we have symmetrised the expression in order to avoid any preference of the one time over the other. So, we conclude that the path differential measure takes the form:

$\begin{array}{c}{\left(\frac{1}{2\pi i\epsilon }\right)}^{\frac{N+1}{2}}{\displaystyle \underset{n=1}{\overset{N}{\prod }}\left[\text{d}{Q}_n\right]}={\displaystyle \underset{n=1}{\overset{N+1}{\prod }}{\left(\frac{1}{2\pi i{\sigma }_n\rho \left(s,t_n\right)\rho \left(s,t_{n-1}\right)}\right)}^{1/2}}{\displaystyle \underset{n=1}{\overset{N}{\prod }}\left[{\rho }_n\left(s,t_n\right)\text{d}{X}_n\right]}\\ =\frac{1}{{\left({\rho }_f{\rho }_i\right)}^{1/2}}{\displaystyle \underset{n=1}{\overset{N+1}{\prod }}{\left(\frac{1}{2\pi i{\sigma }_n}\right)}^{1/2}}{\displaystyle \underset{n=1}{\overset{N}{\prod }}\left[\text{d}{X}_n\right]}\end{array}$ (18)

and the discretized action appearing in Equation (17) is:

$\begin{array}{c}{S}^{\left(N\right)}={\displaystyle \underset{n=1}{\overset{N+1}{\sum }}\left[\frac{{\left(Q_n-Q_{n-1}\right)}^2}{2\epsilon }-\epsilon \frac{1}{2}{\omega }_0^2{Q}_n^2-\epsilon \frac{1}{2}\Gamma s_n^2{Q}_n^2\right]}\\ ={\displaystyle \underset{n=1}{\overset{N+1}{\sum }}[\frac{\bar{\rho }\left(s,{\tau }_n\right)X_n^2}{2{\sigma }_n\bar{\rho }\left(s,{\tau }_{n-1}\right)}+\frac{\bar{\rho }\left(s,{\tau }_{n-1}\right)X_{n-1}^2}{2{\sigma }_n\bar{\rho }\left(s,{\tau }_n\right)}-\frac{X_n{X}_{n-1}}{{\sigma }_n}}\\ -{\sigma }_n\frac{1}{2}{\omega }_0^2\bar{\rho }\left(s,{\tau }_{n-1}\right){\bar{\rho }}^3\left(s,{\tau }_n\right)X_n^2-{\sigma }_n\frac{1}{2}\Gamma \bar{\rho }\left(s,{\tau }_{n-1}\right){\bar{\rho }}^3\left(s,{\tau }_n\right)s_n^2{X}_n^2]\end{array}$ (19)

Therefore, on using the expansions $\frac{\bar{\rho }\left(s,\tau \right)}{\bar{\rho }\left(s,\tau \pm \sigma \right)}=1\mp \frac{\stackrel{\dot{}}{\bar{\rho }}\left(s,\tau \right)}{\bar{\rho }\left(s,\tau \right)}\sigma +\left(\frac{{\stackrel{\dot{}}{\bar{\rho }}}^2\left(s,\tau \right)}{{\bar{\rho }}^2\left(s,\tau \right)}-\frac{\ddot{\stackrel{¯}{\rho }}\left(s,\tau \right)}{2\bar{\rho }\left(s,\tau \right)}\right){\sigma }^2+O\left({\sigma }^3\right)$ , we find that the propagator $K_0\left(Q_f,Q_i,t;s\right)$ is related with the transformed one via:

$K_0\left(Q_f,Q_i,t;s\right)=\frac{1}{{\left({\rho }_f{\rho }_i\right)}^{1/2}}\exp \left\{\frac{i}{2}\left(\frac{{\stackrel{\dot{}}{\bar{\rho }}}_f}{{\rho }_f}{X}_f^2-\frac{{\stackrel{\dot{}}{\bar{\rho }}}_i}{{\rho }_i}{X}_i^2\right)\right\}{K}_0\left(X_f,X_i,\tau ;s\right)$ (20)

where on switching to a phase space path integral, we get:

$\begin{array}{l}{K}_0\left(X_f,X_i,\tau ;s\right)={\displaystyle \iint DX\frac{D{P}_1}{2\pi }}\\ \times \exp \left\{i{\displaystyle \underset{0}{\overset{\tau }{\int }}\text{d}\tau \left[P_1\stackrel{\dot{}}{X}-\left(\frac{P_1^2}{2}+\frac{1}{2}\left({\tilde{\omega }}^2\left(s,\tau \right)+{\omega }_0^2{\bar{\rho }}^4\left(s,\tau \right)\right)X^2+\frac{1}{2}\Gamma s^2\left(\tau \right){\bar{\rho }}^4\left(s,\tau \right)X^2\right)\right]}\right\}\end{array}$ (21)

where we set:

${\tilde{\omega }}^2\left(s,\tau \right)=\left[\frac{\ddot{\stackrel{¯}{\rho }}\left(s,\tau \right)}{\bar{\rho }\left(s,\tau \right)}-2{\left(\frac{\stackrel{\dot{}}{\bar{\rho }}\left(s,\tau \right)}{\bar{\rho }\left(s,\tau \right)}\right)}^2\right]={\rho }^3\left(s,\tau \right)\ddot{\rho }\left(s,\tau \right)$ (22)

and we have used the notation:

$\bar{\rho }\left(s,\tau \right)=\rho \left(s,t\right)$ (23)

$\stackrel{\dot{}}{\rho }\left(s,t\right)=\frac{\partial \rho \left(s,t\right)}{\partial t}$ (24)

$\stackrel{\dot{}}{\bar{\rho }}\left(s,\tau \right)=\frac{\partial \bar{\rho }\left(s,\tau \right)}{\partial \tau }$ (25)

Now, we impose constrain on $\rho$ by setting the global time-dependent term multiplying the $X^2$ terms in Equation (21) equal to a constant. i.e.

${\tilde{\omega }}^2\left(s,\tau \right)+Y^2\left(s\left(\tau \right),\tau \right){\bar{\rho }}^4\left(s,\tau \right)=1$ (26)

where we set:

$Y^2\left(s\left(\tau \right),\tau \right)={\omega }_0^2+\Gamma s^2\left(\tau \right)$ (27)

Further, the integration with respect the $\left(X,P_1\right)$ variables is performed and we find the propagator

$\begin{array}{c}{K}_0\left(Q_f,Q_i,t;s\right)=\sqrt{\frac{1}{2\pi i{\rho }_f\left(s\right){\rho }_i\left(s\right)\sin \phi \left(s,t\right)}}\exp \left\{\frac{i}{2}\left(\frac{{\stackrel{\dot{}}{\rho }}_f\left(s\right)}{{\rho }_f\left(s\right)}{Q}_f^2-\frac{{\stackrel{\dot{}}{\rho }}_i\left(s\right)}{{\rho }_i\left(s\right)}{Q}_i^2\right)\right\}\\ \times \exp \left\{\frac{i}{2\sin \phi \left(s,t\right)}\left[\left(\frac{Q_f^2}{{\rho }_f^2\left(s\right)}+\frac{Q_i^2}{{\rho }_i^2\left(s\right)}\right)\cos \phi \left(s,t\right)-\frac{2Q_f{Q}_i}{{\rho }_f\left(s\right){\rho }_i\left(s\right)}\right]\right\}\end{array}$ (28)

where $\rho \left(s,t\right)$ is the solution of the differential equation:

$\ddot{\rho }\left(s,t\right)+Y^2\left(s\left(t\right),t\right)\rho \left(s,t\right)=\frac{1}{{\rho }^3\left(s,t\right)}$ (29)

and

$\phi \left(s,t^{\prime }\right)={\displaystyle \underset{0}{\overset{t^{\prime }}{\int }}\text{d}\tau \frac{1}{{\rho }^2\left(s,\tau \right)}}$ (30)

In (28), we have set ${\rho }_i\left(s\right)=\rho \left(s,0\right)$ and ${\rho }_f\left(s\right)=\rho \left(s,t\right)$ .

So, finally, we can obtain the propagator $K\left(x_f,x_i,t;s\right)$ from Equations (9) and (10). It has the form:

$\begin{array}{c}K\left(x_f,x_i,t;s\right)=\sqrt{\frac{m}{2\pi i{\rho }_f\left(s\right){\rho }_i\left(s\right)\sin \phi \left(s,t\right)}}\exp \left\{\frac{im}{2}\left(\frac{{\stackrel{\dot{}}{\rho }}_f\left(s\right)}{{\rho }_f\left(s\right)}{x}_f^2-\frac{{\stackrel{\dot{}}{\rho }}_i\left(s\right)}{{\rho }_i\left(s\right)}{x}_i^2\right)\right\}\\ \times \exp \left\{\frac{im}{2\sin \phi \left(s,t\right)}\left[\left(\frac{x_f^2}{{\rho }_f^2\left(s\right)}+\frac{x_i^2}{{\rho }_i^2\left(s\right)}\right)\cos \phi \left(s,t\right)-\frac{2x_f{x}_i}{{\rho }_f\left(s\right){\rho }_i\left(s\right)}\right]\right\}\end{array}$ (31)

The differential Equation (29) takes the form:

$\ddot{\rho }\left(s,t\right)+\left({\omega }_0^2+\frac{\lambda }{m}{s}^2\right)\rho \left(s,t\right)=\frac{1}{{\rho }^3\left(s,t\right)}$ (32)

and

$\phi \left(s,t^{\prime }\right)={\displaystyle \underset{0}{\overset{t^{\prime }}{\int }}\frac{\text{d}\tau }{{\rho }^2\left(s,\tau \right)}}$ (33)

Therefore, the propagator can be derived from the system of Equations (31)-(33). We have to solve the differential Equation (32) with the variable s as parameter, evaluate Equation (33) and then apply (31). We notice that in Equations (3) and (31)-(33) instead of the square function $s^2\left(t\right)$ there could appear any function of s.

Now, we assume that the harmonic oscillator is initially at the state:

$\Phi \left(x\right)={\left(\frac{m{\omega }_0}{\pi }\right)}^{1/4}\exp \left[-\frac{m{\omega }_0}{2}{\left(x-\delta \right)}^2\right]$ (34)

Then, $W_x\left(0\right)$ in Equation (6) becomes:

$W_x\left(0\right)=|\Phi \rangle \langle \Phi |$ (35)

and the influence functional describing the effect of the oscillator on the system is:

$\begin{array}{l}R\left(s,s^{\prime }\right)={\displaystyle \int \text{d}{x}_f}\text{d}{x}_i\text{d}{x^{\prime }}_{i}K\left(x_f,x_i,t;s\right)K^{\ast }\left(x_f,x_i,t;s^{\prime }\right)\Phi \left(x_i\right){\Phi }^{\ast }\left({x^{\prime }}_i\right)\\ =\sqrt{i{\omega }_0}\exp \left[-m{\omega }_0{\delta }^2\right]\sqrt{\frac{2{\rho }_i\left(s\right){\rho }_i\left(s^{\prime }\right)}{{\rho }_f\left(s\right){\rho }_f\left(s^{\prime }\right)}}\sqrt{\sin \phi \left(s,t\right)\sin \phi \left(s^{\prime },t\right)}\\ \times \frac{1}{\sqrt{U\left(s,s^{\prime },t\right)T\left(s\right)T\left(s^{\prime }\right)+\frac{\sin \phi \left(s,t\right)}{{\rho }_f^2\left(s^{\prime }\right)}T\left(s\right)-\frac{\sin \phi \left(s^{\prime },t\right)}{{\rho }_f^2\left(s\right)}{T}^{\prime }\left(s^{\prime }\right)}}\\ \times \exp \{im\frac{{\delta }^2{\omega }_0^2}{2}[\frac{{\rho }_i^2\left(s\right)\sin \phi \left(s,t\right)}{T\left(s\right)}-\frac{{\rho }_i^2\left(s^{\prime }\right)\sin \phi \left(s^{\prime },t\right)}{T^{\prime }\left(s^{\prime }\right)}\begin{array}{c}\\ \\ \\ \end{array}\\ +\frac{\left(\frac{2{\rho }_i\left(s\right){\rho }_i\left(s^{\prime }\right)}{{\rho }_f\left(s\right){\rho }_f\left(s^{\prime }\right)}+\frac{{\rho }_i^2\left(s\right)T^{\prime }\left(s^{\prime }\right)}{{\rho }_f^2\left(s\right)T\left(s\right)}+\frac{{\rho }_i^2\left(s^{\prime }\right)T\left(s\right)}{{\rho }_f^2\left(s^{\prime }\right)T^{\prime }\left(s^{\prime }\right)}\right)\sin \phi \left(s,t\right)\sin \phi \left(s^{\prime },t\right)}{U\left(s,s^{\prime },t\right)T\left(s\right)T\left(s^{\prime }\right)+\frac{\sin \phi \left(s,t\right)}{{\rho }_f^2\left(s^{\prime }\right)}T\left(s\right)-\frac{\sin \phi \left(s^{\prime },t\right)}{{\rho }_f^2\left(s\right)}{T}^{\prime }\left(s^{\prime }\right)}]\}\end{array}$ (36)

where we set:

$\begin{array}{c}U\left(s,s^{\prime },t\right)=\left(\frac{{\stackrel{\dot{}}{\rho }}_f\left(s\right)}{{\rho }_f\left(s\right)}-\frac{{\stackrel{\dot{}}{\rho }}_f\left(s^{\prime }\right)}{{\rho }_f\left(s^{\prime }\right)}\right)\sin \phi \left(s,t\right)\sin \phi \left(s^{\prime },t\right)+\frac{\sin \phi \left(s^{\prime },t\right)\cos \phi \left(s,t\right)}{{\rho }_f^2\left(s\right)}\\ -\frac{\sin \phi \left(s,t\right)\cos \phi \left(s^{\prime },t\right)}{{\rho }_f^2\left(s^{\prime }\right)}\end{array}$ (37)

and

$T\left(s\right)=\cos \phi \left(s,t\right)+{\rho }_i\left(s\right)\left(-{\stackrel{\dot{}}{\rho }}_i\left(s\right)+i{\omega }_0{\rho }_i\left(s\right)\right)\sin \phi \left(s,t\right)$ (38)

$T^{\prime }\left(s^{\prime }\right)=\cos \phi \left(s^{\prime },t\right)-{\rho }_i\left(s^{\prime }\right)\left({\stackrel{\dot{}}{\rho }}_i\left(s^{\prime }\right)+i{\omega }_0{\rho }_i\left(s^{\prime }\right)\right)\sin \phi \left(s^{\prime },t\right)$ (39)

So, the effect of the bath on the system is described by the influence functional (7)-(8) while the effect of the oscillator by the influence functional (36). Now, we turn our attention to the system.

4. Solution of the System’s Sign Problem

According to standard path integral methods, as well as the discussion in the previous sections, the system’s density matrix at time t is going to have the path integral representation:

$\begin{array}{l}{W}_{red}\left(s^{+},s^{-};t\right)\\ ={\displaystyle \underset{n=0}{\overset{N}{\prod }}\left[{\displaystyle \int \text{d}{s}_n\text{d}{s^{\prime }}_n}\right]}{\displaystyle \underset{n=1}{\overset{N+1}{\prod }}\left[{\displaystyle \int \frac{\text{d}{p}_{sn}}{2\pi }\frac{\text{d}{p^{\prime }}_{sn}}{2\pi }}\right]}\langle s_0|W_s\left(0\right)|{s^{\prime }}_0\rangle h\left(s_0,s_1,s_2,\cdots ,s_N,s^{+},{s^{\prime }}_0,{s^{\prime }}_1,{s^{\prime }}_2,\cdots ,{s^{\prime }}_N,s^{-}\right)\\ \times \exp \left\{i{\displaystyle \underset{n=1}{\overset{N+1}{\sum }}\left[p_{sn}\left(s_n-s_{n-1}\right)-{p^{\prime }}_{sn}\left({s^{\prime }}_n-{s^{\prime }}_{n-1}\right)-\epsilon H_s\left(p_{sn},s_n\right)+\epsilon H_s\left({p^{\prime }}_{sn},{s^{\prime }}_n\right)\right]}\right\}\end{array}$ (40)

where we have set $s_{N+1}=s^{+}$ , ${s^{\prime }}_{N+1}=s^{-}$ and $\epsilon =\frac{t}{N+1}$ . $N+1$ is the number of time slices in the path integral and $h\left(s_0,s_1,s_2,\cdots ,s_N,s^{+},{s^{\prime }}_0,{s^{\prime }}_1,{s^{\prime }}_2,\cdots ,{s^{\prime }}_N,s^{-}\right)$ is the influence functional which describes the possible interactions of the system. Expression (40) is highly oscillatory and therefore standard Monte Carlo techniques fail to confront it. We can bypass the whole point if we interpret the Hamiltonians as random variables and apply the central limit theorem on the phase of Equation (40). In that way, we obtain an oscillation free expression called sign solved density matrix.

Now, we observe that the influence functional has the form:

$h\left(s,s^{\prime }\right)=I\left(s,s^{\prime }\right)R\left(s,s^{\prime }\right)$ (41)

(see Equations (7), (8) and (36)). It has a product form as the bath and the harmonic oscillator are not directly coupled.

The bath influence functional $I\left(s,s^{\prime }\right)$ in its discrete form is:

$I\left(s,s^{\prime }\right)=\exp \left[-{\displaystyle \underset{n=0}{\overset{N+1}{\sum }}{\displaystyle \underset{n^{\prime }=0}{\overset{n}{\sum }}\left(s_n-{s^{\prime }}_n\right)\left({\zeta }_{n,n^{\prime }}{s}_{n^{\prime }}-{\zeta }_{n,n^{\prime }}^{\ast }{s^{\prime }}_{n^{\prime }}\right)}}\right]$ (42)

The matrix elements ${\zeta }_{n,n^{\prime }}$ are given in [15] . As we can observe there, $N+2$ of them, corresponding to $n=n^{\prime }$ , are of order $\epsilon$ and the rest ones are of order ${\epsilon }^2$ . So, for N large enough, there are positive constants such that:

$\left|{\zeta }_{n,n^{\prime }}\right|<\frac{2C_1}{{\left(N+1\right)}^2}\begin{array}{cc}& n\ne n^{\prime }\end{array}$ (43)

$\left|{\zeta }_{n,n}\right|<\frac{C_2}{N+1}$ (44)

Further, expression (31) is bounded and therefore its matrix elements are bounded as well. So, $R\left(s,s^{\prime }\right)$ is bounded. As we prove in Appendix, a theory similar to the sign solved propagator theory of Ref. [7] and Ref. [12] applies. Eventually, we obtain the sign solved influence functional expression:

$\begin{array}{c}{W}_{red}\left(s^{+},s^{-};t\right)=\underset{N\to \infty }{\lim }h\left(\underset{N+2}{\underset{︸}{s^{+},s^{+},s^{+},\cdots ,s^{+},s^{+}}},\underset{N+2}{\underset{︸}{s^{-},s^{-},s^{-},\cdots ,s^{-},s^{-}}}\right)\\ \times \langle s^{+}|W_s\left(0\right)|s^{-}\rangle \exp \left\{-i\left[{\langle H_s\rangle }^{+}-{\langle H_s\rangle }^{-}\right]t\right\}\end{array}$ (45)

Now, we observe that the primed and unprimed variables of h in Equation (45) and therefore of I and R have a diagonal form and therefore $s\left(t\right)=s^{+}$ and $s^{\prime }\left(t\right)=s^{-}$ . So, eventually, the influence functional given by Equations (7) and (8) takes the diagonal form:

$I\left(s^{+},s^{-};t\right)=\exp \left[-\left(s^{+}-s^{-}\right)\left(\eta \left(t\right)s^{+}-{\eta }^{\ast }\left(t\right)s^{-}\right)\right]$ (46)

where we set:

$\eta \left(t\right)=\frac{1}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}\text{d}{\omega }_1\frac{J\left({\omega }_1\right)}{{\omega }_1^2}\left(2\coth \left(\frac{{\omega }_1\beta }{2}\right){\sin }^2\left(\frac{{\omega }_{1}t}{2}\right)+i\sin \left({\omega }_{1}t\right)\right)}$ (47)

We have evaluated the time integrals as the position variables are diagonal.

In the present paper, we use the ohmic spectral density:

$J\left(\omega \right)=\gamma \omega \exp \left[-\frac{\omega }{{\omega }_c}\right]$ (48)

where ${\omega }_c$ is a cutoff frequency.

Moreover, by setting $s_c$ to correspond either to $s^{+}$ or to $s^{-}$ , the propagator (31) becomes:

$\begin{array}{l}K\left(x_f,x_i,t;s_c\right)\\ =\sqrt{\frac{m\sqrt{{\omega }_0^2+\frac{\lambda s_c^2}{m}}}{2\pi i\sin \left(\sqrt{{\omega }_0^2+\frac{\lambda s_c^2}{m}}t\right)}}\exp \left\{\frac{im\sqrt{{\omega }_0^2+\frac{\lambda s_c^2}{m}}}{2\sin \left(\sqrt{{\omega }_0^2+\frac{\lambda s_c^2}{m}}t\right)}\left[\left(x_f^2+x_i^2\right)\cos \left(\sqrt{{\omega }_0^2+\frac{\lambda s_c^2}{m}}t\right)-2x_f{x}_i\right]\right\}\end{array}$ (49)

and according to Equation (36), the influence functional R has the diagonal form:

$\begin{array}{l}R\left(s^{+},s^{-}\right)=\sqrt{2i{\omega }_0}\exp \left[-m{\omega }_0{\delta }^2\right]\sqrt{\sin \phi \left(s^{+},t\right)\sin \phi \left(s^{-},t\right)}\\ \times \frac{1}{\sqrt{U_1\left(s^{+},s^{-},t\right)T_1\left(s^{+}\right){T^{\prime }}_1\left(s^{-}\right)+\frac{\sin \phi \left(s^{+},t\right)}{{\rho }_f^2\left(s^{-}\right)}{T}_1\left(s^{+}\right)-\frac{\sin \phi \left(s^{-},t\right)}{{\rho }_f^2\left(s^{+}\right)}{T^{\prime }}_1\left(s^{-}\right)}}\\ \times \exp [im\frac{{\delta }^2{\omega }_0^2}{2}\left(\frac{{\rho }_i^2\left(s^{+}\right)\sin \phi \left(s^{+},t\right)}{T_1\left(s^{+}\right)}-\frac{{\rho }_i^2\left(s^{-}\right)\sin \phi \left(s^{-},t\right)}{{T^{\prime }}_1\left(s^{-}\right)}\right)\begin{array}{c}\\ \\ \\ \end{array}\\ +\frac{im\frac{{\delta }^2{\omega }_0^2}{2}\left(2+\frac{{T^{\prime }}_1\left(s^{-}\right)}{T_1\left(s^{+}\right)}+\frac{T_1\left(s^{+}\right)}{{T^{\prime }}_1\left(s^{-}\right)}\right)\sin \phi \left(s^{+},t\right)\sin \phi \left(s^{-},t\right)}{U_1\left(s^{+},s^{-},t\right)T_1\left(s^{+}\right){T^{\prime }}_1\left(s^{-}\right)+\frac{\sin \phi \left(s^{+},t\right)}{{\rho }_f^2\left(s^{-}\right)}{T}_1\left(s^{+}\right)-\frac{\sin \phi \left(s^{-},t\right)}{{\rho }_f^2\left(s^{+}\right)}{T^{\prime }}_1\left(s^{-}\right)}]\end{array}$ (50)

where we set:

$U_1\left(s^{+},s^{-},t\right)=\frac{\sin \phi \left(s^{-},t\right)\cos \phi \left(s^{+},t\right)}{{\rho }_f^2\left(s^{+}\right)}-\frac{\sin \phi \left(s^{+},t\right)\cos \phi \left(s^{-},t\right)}{{\rho }_f^2\left(s^{-}\right)}$ (51)

and

$T_1\left(s^{+}\right)=\cos \phi \left(s^{+},t\right)+i{\omega }_0{\rho }_i^2\left(s^{+}\right)\sin \phi \left(s^{+},t\right)$ (52)

${T^{\prime }}_1\left(s^{-}\right)=\cos \phi \left(s^{-},t\right)-i{\omega }_0{\rho }_i^2\left(s^{-}\right)\sin \phi \left(s^{-},t\right)$ (53)

Here,

${\rho }_i\left(s_c\right)={\rho }_f\left(s_c\right)={\left({\omega }_0^2+\frac{\lambda s_c^2}{m}\right)}^{-\frac{1}{4}}$ (54)

and

$\phi \left(s_c,t\right)=\sqrt{{\omega }_0^2+\frac{\lambda s_c^2}{m}}t$ (55)

Finally, we obtain:

$W_{red}\left(s^{+},s^{-};t\right)=I\left(s^{+},s^{-};t\right)R\left(s^{+},s^{-};t\right)\langle s^{+}|W_s\left(0\right)|s^{-}\rangle \exp \left\{-i\left[{\langle H_s\rangle }^{+}-{\langle H_s\rangle }^{-}\right]t\right\}$ (56)

In the next section, we proceed to an application.

5. Application to a Symmetric Double Well

In the present section, we proceed to an application of the above theory. We consider the symmetric double well potential:

$V\left(s\right)=\frac{M{\omega }^2}{8a^2}{\left(s-a\right)}^2{\left(s+a\right)}^2$ (57)

and assume to have prepared a Gaussian wave packet centered on the right well with wavefunction:

$\Psi \left(s\right)={\left(\frac{M\omega }{\pi }\right)}^{1/4}\exp \left[-\frac{M\omega }{2}{\left(s-a\right)}^2\right]$ (58)

and energy $E_0$ . Due to tunneling, the level $E_0$ splits into the levels $E_1$ and $E_2$ with corresponding wavefunctions:

${\Psi }_1\left(s\right)=\frac{1}{\sqrt{2}}\left[\Psi \left(s\right)+\Psi \left(-s\right)\right]$ (59)

${\Psi }_2\left(s\right)=\frac{1}{\sqrt{2}}\left[\Psi \left(s\right)-\Psi \left(-s\right)\right]$ (60)

and energy differences:

$E_2-E_0=E_0-E_1=\frac{{\left(M\omega \right)}^{3/2}}{\sqrt{\pi }M}a{\text{e}}^{-M\omega a^2}$ (61)

To apply the theory of the previous section, we use the initial density matrix:

$W_s\left(0\right)=|\Psi \rangle \langle \Psi |$ (62)

and insert it in Equation (56).

Now, we are in position to generate the time evolution of the system’s density matrix. In that way, the evolution of the system’s observables under the preparation (58) can be studied.

We consider the matrix elements:

${\beta }_{ij}\left(t\right)=\langle {\Psi }_i|W_{red}\left(s^{+},s^{-};t\right)|{\Psi }_j\rangle \begin{array}{cc}& i,j=1,2\end{array}$ (63)

Then,

$W_{red}\left(t\right)={\displaystyle \underset{i,j=1}{\overset{2}{\sum }}|{\Psi }_i\rangle {\beta }_{ij}\left(t\right)\langle {\Psi }_j|}$ (64)

In the expectation values in Equation (56), we choose as sampling functions the expressions (58)-(60), so that ${\langle H_s\rangle }^{+}-{\langle H_s\rangle }^{-}=E_i-E_0$ $i=1,2$ where the i is the same as in Equations (63) and (64). Then, $E_2-E_0=-\left(E_1-E_0\right)=\hslash \Omega$ . $\Omega$ is the tunnelling frequency.

So, the expectation value of the position has the form:

$\langle s\left(t\right)\rangle =a\left({\beta }_{12}\left(t\right)+{\beta }_{21}\left(t\right)\right)$ (65)

We observe that the expectation value of the position is closely related with the inversion of a corresponding two-level system. In fact, as we can conclude from the above analysis the present double well potential can be interpreted as a two-level system.

Throughout the section, we consider the system initially in the state (62).

In Figure 1, we plot the expectation value $\langle s\left(t\right)\rangle$ in the case of the system interacting with only the harmonic bath. So, we use the influence functional $h=I$ . In that case, the matrix elements (63) are Gaussian and we obtain:

$\begin{array}{l}{\beta }_{12}\left(t\right)={\beta }_{21}^{\ast }\left(t\right)=\frac{1}{2}{\text{e}}^{i\Omega t}\frac{M\omega }{\sqrt{{\left|M\omega +\eta \left(t\right)\right|}^2-{\left(\mathrm{Re}\left[\eta \left(t\right)\right]\right)}^2}}\\ \times (-exp\left[-2M\omega a^2\right]+\exp \left[-2M\omega a^2\frac{{\left|\eta \left(t\right)\right|}^2-{\left(\mathrm{Re}\left[\eta \left(t\right)\right]\right)}^2}{{\left|M\omega +\eta \left(t\right)\right|}^2-{\left(\mathrm{Re}\left[\eta \left(t\right)\right]\right)}^2}\right]\\ -\text{exp}\left[-M\omega a^2\frac{{\left|M\omega +\eta \left(t\right)\right|}^2-2{\left(\mathrm{Re}\left[\eta \left(t\right)\right]\right)}^2+{\eta }^{\ast }\left(t\right)\left(M\omega +\eta \left(t\right)\right)}{{\left|M\omega +\eta \left(t\right)\right|}^2-{\left(\mathrm{Re}\left[\eta \left(t\right)\right]\right)}^2}\right]\\ +\exp \left[-M\omega a^2\frac{{\left|M\omega +\eta \left(t\right)\right|}^2-2{\left(\mathrm{Re}\left[\eta \left(t\right)\right]\right)}^2+\eta \left(t\right)\left(M\omega +{\eta }^{\ast }\left(t\right)\right)}{{\left|M\omega +\eta \left(t\right)\right|}^2-{\left(\mathrm{Re}\left[\eta \left(t\right)\right]\right)}^2}\right])\end{array}$ (66)

As the time increases, the bath causes decrease of the amplitude of the tunnelling oscillations.

In Figure 2(a) and Figure 2(b), we consider the system interacting with just the oscillator. So, we use the influence functional $h=R$ . At small times there appear extra oscillations, besides the tunnelling’s ones. Moreover, for fixed time the absolute value $\left|\langle s\left(t\right)\rangle \right|$ decreases as $\delta$ increases (see the initial wavefunction (34) of the oscillator and the influence functional (50)).

In Figure 3, we consider the full system. So, we use the influence functional $h=IR$ . There appears a combination of the effects described in the cases in Figure 1 and Figure 2. Oscillations at small times and decrease of the amplitude as time increases. We should expect such a result as the full influence functional is the product of the influence functionals corresponding to the system’s interaction with just the bath or the oscillator (see Equation (41)).

6. Conclusions

In the present paper, we study the dynamics of a quantum mechanical system interacting linearly with a bath and quadratically with an oscillator. In our study, we use influence functional methods derived previously and concerning the interaction of systems with baths as well as methods on the dynamics of oscillators and combine them with methods on the solution of the sign problem due to the author. As an application, we have considered a symmetric double well interacting with a bath or with an oscillator, or with both of them, and study the time evolution of the relevant density matrix. We focus on a double well potential as it incorporates tunnelling effects on the whole dynamics and behaves as a two-level system.

Figure 1. Mean position of the dissipative double well system. We use $\gamma =0.06$ (solid), $\gamma =0.4$ (dashed) and $\gamma =0.7$ (dotted). We have set: $\omega =1.0$ , $M=1.0$ , $a=2.0$ , ${\omega }_c=8.0\Omega$ , $\beta =\frac{1}{2.5\Omega }$ . Here, according to Equation (61), $\Omega =0.020667$ .

(a)(b)

Figure 2. Mean position of a double well system coupled with an oscillator. In (a), we consider the mean position as a function of $\Omega t$ ; while in (b), we consider the position as a function of $\delta$ for certain values of $\Omega t$ . We use: $\omega =1.0$ , $M=1.0$ , $a=2.0$ , ${\omega }_0=1.0$ , $m=1.0$ , $\lambda =0.1$ . Here, according to Equation (61), $\Omega =0.020667$ .

Figure 3. Mean position in the following cases. A double well: Dashed-dotted. A dissipative double well: Solid. A double well coupled with an oscillator: Dotted. A double well coupled with a harmonic bath and an oscillator: Dashed. The parameters have the values: $\omega =1.0$ , $M=1.0$ , $a=2.0$ , $\delta =0.0$ , ${\omega }_0=1.0$ , $m=1.0$ , $\lambda =0.1$ , $\gamma =0.5$ , ${\omega }_c=8.0\Omega$ , $\beta =\frac{1}{2.5\Omega }$ . Here, according to Equation (61), $\Omega =0.020667$ .

In conclusion, the present methods are capable for giving closed expressions on various systems’ density matrix time evolution and therefore interesting relevant dynamical information can be gained.

Appendix: Sign Solved Density Matrix

According to the solution of the sign problem for expression (40),

$\begin{array}{l}{W}_{red}\left(s^{+},s^{-};t\right)\\ ={\displaystyle \underset{n=0}{\overset{N}{\prod }}\left[{\displaystyle \int \text{d}{s}_n\text{d}{s^{\prime }}_n}\right]}{\displaystyle \underset{n=1}{\overset{N+1}{\prod }}\left[{\displaystyle \int \frac{\text{d}{p}_{sn}}{2\pi }\frac{\text{d}{p^{\prime }}_{sn}}{2\pi }}\right]}\langle s_0|W_s\left(0\right)|{s^{\prime }}_0\rangle h\left(s_0,s_1,s_2,\cdots ,s_N,s^{+},{s^{\prime }}_0,{s^{\prime }}_1,{s^{\prime }}_2,\cdots ,{s^{\prime }}_N,s^{-}\right)\\ \times \exp \left\{i{\displaystyle \underset{n=1}{\overset{N+1}{\sum }}\left[p_{sn}\left(s_n-s_{n-1}\right)-{p^{\prime }}_{sn}\left({s^{\prime }}_n-{s^{\prime }}_{n-1}\right)-\epsilon H_s\left(p_{sn},s_n\right)+\epsilon H_s\left({p^{\prime }}_{sn},{s^{\prime }}_n\right)\right]}\right\}\\ \cong h\left(\underset{N+2}{\underset{︸}{s^{+},s^{+},s^{+},\cdots ,s^{+},s^{+}}},\underset{N+2}{\underset{︸}{s^{-},s^{-},s^{-},\cdots ,s^{-},s^{-}}}\right)\langle s^{+}|W_s\left(0\right)|s^{-}\rangle \exp \left\{-i\left[{\langle H_s\rangle }^{+}-{\langle H_s\rangle }^{-}\right]t\right\}\\ +\frac{1}{\sqrt{N+1}}{\displaystyle \underset{n=0}{\overset{N}{\prod }}\left[{\displaystyle \int \text{d}{s}_n\text{d}{s^{\prime }}_n}\right]h}\left(s_0,s_1,s_2,\cdots ,s_N,s^{+},{s^{\prime }}_0,{s^{\prime }}_1,{s^{\prime }}_2,\cdots ,{s^{\prime }}_N,s^{-}\right)\\ \times \langle s_0|W_s\left(0\right)|{s^{\prime }}_0\rangle \left\{{\displaystyle \underset{n=1}{\overset{N+1}{\prod }}\left[f\left(s_n,s_{n-1}\right)f^{\prime }\left({s^{\prime }}_n,{s^{\prime }}_{n-1}\right)\right]}-i{\displaystyle \underset{n=1}{\overset{N+1}{\prod }}\left[g\left(s_n,s_{n-1}\right)g^{\prime }\left({s^{\prime }}_n,{s^{\prime }}_{n-1}\right)\right]}\right\}\end{array}$ (A1)

We have used the functions:

$\begin{array}{l}f\left(s_n,s_{n-1}\right)\\ =\frac{1}{2\pi \sqrt{2\pi }{\sigma }_V\left(t\right)t\sqrt{\left[{\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)\right]}}\\ \times \exp \{-\frac{1}{2}{\sigma }_m^2{t}^2\left[{\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)\right]{\left(s_n-s_{n-1}\right)}^2\begin{array}{c}\\ \\ \end{array}\\ -\frac{s_n^2}{2{\sigma }_V^2\left(t\right)t^2\left[{\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)\right]}\}\end{array}$ (A2)

$\begin{array}{l}g\left(s_n,s_{n-1}\right)\\ =\frac{1}{2\pi \sqrt{2\pi }{\sigma }_V\left(t\right)t\sqrt{\left[{\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)\right]}}\\ \times \exp \{-\frac{1}{2}{\sigma }_m^2{t}^2\left[{\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)\right]{\left(s_n-s_{n-1}\right)}^2\begin{array}{c}\\ \\ \end{array}\\ -\frac{s_n^2}{2{\sigma }_V^2\left(t\right)t^2\left[{\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)\right]}\}\end{array}$ (A3)

In (A2) and (A3), we use appropriate sampling functions. In the primed f and g in (A1), we use a primed sampling function for the variances.

In calculations, we are interested in integrals of the form ${\displaystyle \iint \text{d}{s}^{+}\text{d}{s}^{-}}{\Theta }_1^{\ast }\left(s^{+}\right)W_{red}\left(s^{+},s^{-};t\right){\Theta }_2\left(s^{-}\right)$ , where ${\Theta }_1$ , ${\Theta }_2$ are appropriate functions. So, we consider the expression $\wp$ corresponding to the term in Equation (A1) involving the $f,f^{\prime },g$ and $g^{\prime }$ functions after that integration. We intent to prove that only the diagonal term in the last expression in Equation (A1) can give the exact result as $N\to \infty$ because $\wp$ tends to zero. Then, we obtain Equation (45). To prove that, we diagonalize and integrate the Gaussian products ${\displaystyle \underset{n=1}{\overset{N+1}{\prod }}\left[f\left(s_n,s_{n-1}\right)f^{\prime }\left({s^{\prime }}_n,{s^{\prime }}_{n-1}\right)\right]I\left(s,s^{\prime }\right)}$ and ${\displaystyle \underset{n=1}{\overset{N+1}{\prod }}\left[g\left(s_n,s_{n-1}\right)g^{\prime }\left({s^{\prime }}_n,{s^{\prime }}_{n-1}\right)\right]I\left(s,s^{\prime }\right)}$ . We notice that according to Equation (50) the R is bounded (see Equation (41)).

According to Equations (42)-(44), the bath influence functional $I\left(s,s^{\prime }\right)$ has the form:

$I\left(s,s^{\prime }\right)=\exp \left[-{\displaystyle \underset{n=0}{\overset{N+1}{\sum }}{\displaystyle \underset{n^{\prime }=0}{\overset{n}{\sum }}\left(s_n-{s^{\prime }}_n\right)\left({\zeta }_{n,n^{\prime }}{s}_{n^{\prime }}-{\zeta }_{n,n^{\prime }}^{\ast }{s^{\prime }}_{n^{\prime }}\right)}}\right]$ (A4)

where for N large enough,

$\left|{\zeta }_{n,n^{\prime }}\right|<\frac{2C_1}{{\left(N+1\right)}^2}\begin{array}{cc}& n\ne n^{\prime }\end{array}$ (A5)

$\left|{\zeta }_{n,n}\right|<\frac{C_2}{N+1}$ (A6)

We want to bound the expression $\left|\wp \right|$ appropriately. We proceed via performing on the terms composed of the f functions, the change of variables:

$s_n\to \frac{s_n}{{\gamma }_1\left(t\right)}$ (A7)

and similarly, the change of variables:

$s_n\to \frac{s_n}{{\gamma }_2\left(t\right)}$ (A8)

on the terms composed of the g functions.

We have set:

${\gamma }_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\left(t\right)=\frac{{\sigma }_{m}t\left\{\begin{array}{l}\sqrt{\left[{\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)\right]}\\ \sqrt{\left[{\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)\right]}\end{array}\right\}}{\sqrt{2}}$ (A9)

Similar primed transformations apply to the case of primed variables.

Then, we obtain:

$\begin{array}{l}\left|\wp \right|\le \frac{2^{N+2}}{{\left[2\pi \right]}^{3N+3}{\left[{\sigma }_V\left(t\right)\right]}^{N+1}{\sigma }_m^{N+2}{\left[{\sigma }_{V^{\prime }}\left(t\right)\right]}^{N+1}{{\sigma }^{\prime }}_m^{N+2}}\\ \times \frac{1}{{\left[t^2\left({\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)\right)\right]}^{2N+3}}\\ \times \frac{1}{\sqrt{N+1}}{\displaystyle \underset{n=0}{\overset{N+1}{\prod }}\left[{\displaystyle \int \text{d}{s}_n\text{d}{s^{\prime }}_n}\right]{\Theta }_1^{\ast }}\left(\frac{s_{N+1}}{{\gamma }_1\left(t\right)}\right){\Theta }_2\left(\frac{{s^{\prime }}_{N+1}}{{{\gamma }^{\prime }}_1\left(t\right)}\right)\langle \frac{s_0}{{\gamma }_1\left(t\right)}|W_s\left(0\right)|\frac{{s^{\prime }}_0}{{{\gamma }^{\prime }}_1\left(t\right)}\rangle \\ \times \exp \left\{{\left(s^{+}\right)}^2+{\left(s^{-}\right)}^2+{\phi }_1{s}_0^2+{{\phi }^{\prime }}_1{s^{\prime }}_0^2+{\stackrel{\to }{\rho }}_1{\stackrel{\leftrightarrow }{M}}_{2N+4}\left({\beta }_1,{{\beta }^{\prime }}_1\right){\stackrel{\to }{\rho }}_1^{\text{T}}+\frac{1}{N+1}{\stackrel{\to }{\rho }}_{1ab}{\stackrel{\leftrightarrow }{M}}_{1,2N+4}^1{\stackrel{\to }{\rho }}_{1ab}^{\text{T}}\right\}\end{array}$

$\begin{array}{l}+\frac{2^{N+2}}{{\left[2\pi \right]}^{3N+3}{\left[{\sigma }_V\left(t\right)\right]}^{N+1}{\sigma }_m^{N+2}{\left[{\sigma }_{V^{\prime }}\left(t\right)\right]}^{N+1}{{\sigma }^{\prime }}_m^{N+2}}\\ \times \frac{1}{{\left[t^2\left({\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)\right)\right]}^{2N+3}}\\ \times \frac{1}{\sqrt{N+1}}{\displaystyle \underset{n=0}{\overset{N+1}{\prod }}\left[{\displaystyle \int \text{d}{s}_n\text{d}{s^{\prime }}_n}\right]{\Theta }_1^{\ast }}\left(\frac{s_{N+1}}{{\gamma }_2\left(t\right)}\right){\Theta }_2\left(\frac{{s^{\prime }}_{N+1}}{{{\gamma }^{\prime }}_2\left(t\right)}\right)\langle \frac{s_0}{{\gamma }_2\left(t\right)}|W_s\left(0\right)|\frac{{s^{\prime }}_0}{{{\gamma }^{\prime }}_2\left(t\right)}\rangle \\ \times \exp \left\{{\left(s^{+}\right)}^2+{\left(s^{-}\right)}^2+{\phi }_2{s}_0^2+{{\phi }^{\prime }}_2{s^{\prime }}_0^2+{\stackrel{\to }{\rho }}_1{\stackrel{\leftrightarrow }{M}}_{2N+4}\left({\beta }_2,{{\beta }^{\prime }}_2\right){\stackrel{\to }{\rho }}_1^{\text{T}}+\frac{1}{N+1}{\stackrel{\to }{\rho }}_{1ab}{\stackrel{\leftrightarrow }{M}}_{2,2N+4}^1{\stackrel{\to }{\rho }}_{1ab}^{\text{T}}\right\}\end{array}$ (A10)

${\phi }_{\left\{\begin{array}{c}1\\ 2\end{array}\right\}}$ and ${{\phi }^{\prime }}_{\left\{\begin{array}{c}1\\ 2\end{array}\right\}}$ are appropriate time dependent functions. We have set:

${\stackrel{\to }{\rho }}_1=\left(s_0,s_1,s_2,\cdots ,s_N,s^{+},{s^{\prime }}_0,{s^{\prime }}_1,{s^{\prime }}_2,\cdots ,{s^{\prime }}_N,s^{-}\right)$ (A11)

and

${\stackrel{\to }{\rho }}_{1ab}=\left(\left|s_0\right|,\left|s_1\right|,\left|s_2\right|,\cdots ,\left|s_N\right|,\left|s^{+}\right|,\left|{s^{\prime }}_0\right|,\left|{s^{\prime }}_1\right|,\left|{s^{\prime }}_2\right|,\cdots ,\left|{s^{\prime }}_N\right|,\left|s^{-}\right|\right)$ (A12)

where $s_{N+1}=s^{+}$ , ${s^{\prime }}_{N+1}=s^{-}$ . The matrices in Equation (A10) correspond to the symmetric matrices:

${\stackrel{\leftrightarrow }{M}}_{2N+4}\left({\beta }_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}},{{\beta }^{\prime }}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\right)=\left[\begin{array}{cc}{\stackrel{\leftrightarrow }{M}}_{N+2}\left({\beta }_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\right)& \stackrel{\leftrightarrow }{0}\\ \stackrel{\leftrightarrow }{0}& {\stackrel{\leftrightarrow }{M}}_{N+2}\left({{\beta }^{\prime }}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\right)\end{array}\right]$ (A13)

where we set:

${\left({\stackrel{\leftrightarrow }{M}}_{N+2}\left({\beta }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}\right)\right)}_{ij}=\{\begin{array}{ll}1\hfill & \text{if}i=j\pm 1\hfill \\ {\beta }_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\hfill & \text{if}i=j\hfill \\ 0\hfill & \text{otherwise}\hfill \end{array}$ (A14)

and similarly for the primed variables. Moreover,

${\left({\stackrel{\leftrightarrow }{M}}_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}2N+4}^1\right)}_{ij}=\{\begin{array}{l}\frac{C_2}{{\gamma }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}{\gamma }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}}\text{if}i=j\text{or}i=j\pm \left(N+2\right)\\ \frac{C_1}{N+1}\text{otherwise}\end{array}$ (A15)

${\beta }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}\left(t\right)$ has the form:

${\beta }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}\left(t\right)=-2-\frac{1}{{\sigma }_m^2{\sigma }_V^2\left(t\right)t^4\left\{\begin{array}{l}{\left({\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)\right)}^2\\ {\left({\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)\right)}^2\end{array}\right\}}$ (A16)

Depending on the s being primed, unprimed or mixed we use primed, unprimed or mixed variables in the coefficients (A15) of the quadratic forms in Equation (A10). For example, for the term corresponding to $s_i{s^{\prime }}_j$ in Equation (A10), we use in the expression (A15) the product ${\gamma }_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}{{\gamma }^{\prime }}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}$ .

Now, we study the matrices (A13, A14). We observe that on the one hand,

$\det \left({\stackrel{\leftrightarrow }{M}}_{2N+4}\left({\beta }_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}},{{\beta }^{\prime }}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\right)\right)=\det \left({\stackrel{\leftrightarrow }{M}}_{N+2}\left({\beta }_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\right)\right)\det \left({\stackrel{\leftrightarrow }{M}}_{N+2}\left({{\beta }^{\prime }}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\right)\right)$ (A17)

and on the other for each determinant,

$\det \left({\stackrel{\leftrightarrow }{M}}_{N+2}\left({\beta }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}\right)\right)=U_{N+2}\left(\frac{{\beta }_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}{2}\right)$ (A18)

$U_{N+2}\left(\frac{{\beta }_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}{2}\right)$ is a Chebyshev polynomial of the second kind of order $N+2$ . More particularly let the numbers ${\xi }_{N+2}^{\ast \left(n\right)}$ , $n=0,\cdots ,N+1$ , be the roots of the equation $U_{N+2}\left(x\right)=0$ . They are simple, real roots and ${\xi }_{N+2}^{\ast \left(n\right)}\in \left(-1,1\right)$ , $n=0,\cdots ,N+1$ . Then, the eigenvalues of the matrices ${\stackrel{\leftrightarrow }{M}}_{2N+4}\left({\beta }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}},{{\beta }^{\prime }}_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}\right)$ are going to be given by the expressions:

${\lambda }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}n}^{\ast \left(2N+4\right)}=-2{\xi }_{N+2}^{\ast \left(n\right)}+{\beta }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}\left(t\right)\begin{array}{cc}& n=0,\cdots ,N+1\end{array}$ (A19)

${\lambda }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}n}^{\ast \left(2N+4\right)}=-2{\xi }_{N+2}^{\ast \left(n-N-2\right)}+{{\beta }^{\prime }}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\left(t\right)\begin{array}{cc}& n=N+2,\cdots ,2N+3\end{array}$ (A20)

Further, the diagonal quadratic forms ${\left(s^{+}\right)}^2+{\left(s^{-}\right)}^2+{\phi }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}{s}_0^2+{{\phi }^{\prime }}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}{s^{\prime }}_0^2$ can be diagonalized simultaneously with the quadratic forms corresponding to the matrices ${\stackrel{\leftrightarrow }{M}}_{2N+4}\left({\beta }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}},{{\beta }^{\prime }}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\right)$ . So, we conclude that the eigenvalues of the quadratic forms:

${\left(s^{+}\right)}^2+{\left(s^{-}\right)}^2+{\phi }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}{s}_0^2+{{\phi }^{\prime }}_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}{s^{\prime }}_0^2+{\stackrel{\to }{\rho }}_1{\stackrel{\leftrightarrow }{M}}_{2N+4}\left({\beta }_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}},{{\beta }^{\prime }}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\right){\stackrel{\to }{\rho }}_1^{\text{T}}$ (A21)

are going to have the form:

${\lambda }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}n}^{\left(2N+4\right)}=-2{\xi }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}2N+4}^{\left(n\right)}+{\beta }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}\left(t\right)\begin{array}{cc}& n=0,\cdots ,N+1\end{array}$ (A22)

${\lambda }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}n}^{\left(2N+4\right)}=-2{\xi }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}2N+4}^{\left(n\right)}+{{\beta }^{\prime }}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\left(t\right)\begin{array}{cc}& n=N+2,\cdots ,2N+3\end{array}$ (A23)

where ${\xi }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}2N+4}^{\left(n\right)}$ are appropriate real numbers with ${\xi }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}2N+4}^{\left(n\right)}>-1$ .

Further, for N large enough, the expressions corresponding to the matrices $\frac{1}{N+1}{\stackrel{\leftrightarrow }{M}}_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}2N+4}^1$ in (A10) are perturbation terms. So, eventually the full matrices on the two exponentials will have eigenvalues:

${\mu }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}n}^{\left(2N+4\right)}=-2{\xi }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}2N+4}^{\left(n\right)}+{\beta }_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}\left(t\right)+O\left(\frac{1}{N+1}\right)\begin{array}{cc}& n=0,\cdots ,N+1\end{array}$ (A24)

${\mu }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}n}^{\left(2N+4\right)}=-2{\xi }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}2N+4}^{\left(n\right)}+{{\beta }^{\prime }}_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}}\left(t\right)+O\left(\frac{1}{N+1}\right)\begin{array}{cc}& n=N+2,\cdots ,2N+3\end{array}$ (A25)

According to the whole above discussion after a Gaussian integration, we obtain:

$\begin{array}{l}\left|\wp \right|\le b_1{\left(2\pi \right)}^3{\sigma }_V\left(t\right){\sigma }_{V^{\prime }}\left(t\right)t^2\left({\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)\right)\\ \left\{\frac{1}{\sqrt{N+1}}{\displaystyle \underset{n=0}{\overset{N+1}{\prod }}\left[\frac{1}{2\pi \sqrt{{\Lambda }_{1,n}^{\left(2N+4\right)}}}\right]}{\displaystyle \underset{n=0}{\overset{N+1}{\prod }}\left[\frac{1}{2\pi \sqrt{{{\Lambda }^{\prime }}_{1,n+N+2}^{\left(2N+4\right)}}}\right]}\right\}\\ +b_2{\left(2\pi \right)}^3{\sigma }_V\left(t\right){\sigma }_{V^{\prime }}\left(t\right)t^2\left({\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)\right)\\ \left\{\frac{1}{\sqrt{N+1}}{\displaystyle \underset{n=0}{\overset{N+1}{\prod }}\left[\frac{1}{2\pi \sqrt{{\Lambda }_{2,n}^{\left(2N+4\right)}}}\right]}{\displaystyle \underset{n=0}{\overset{N+1}{\prod }}\left[\frac{1}{2\pi \sqrt{{{\Lambda }^{\prime }}_{2,n+N+2}^{\left(2N+4\right)}}}\right]}\right\}\end{array}$ (A26)

where we set:

$\begin{array}{l}{\Lambda }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}n}^{\left(2N+4\right)}=1+2\left(1+{\xi }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}2N+4}^{\left(n\right)}\right){\sigma }_m^2{\sigma }_V^2\left(t\right)t^4\\ \times \left\{\begin{array}{l}{\left({\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)\right)}^2\\ {\left({\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)\right)}^2\end{array}\right\}\\ +O\left(\frac{1}{N+1}\right),n=0,\cdots ,N+1\end{array}$ (A27)

$\begin{array}{l}{{\Lambda }^{\prime }}_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}n}^{\left(2N+4\right)}=1+2\left(1+{\xi }_{{}_{\left\{\begin{array}{l}1\\ 2\end{array}\right\}}2N+4}^{\left(n\right)}\right){{\sigma }^{\prime }}_m^2{\sigma }_{V^{\prime }}^2\left(t\right)t^4\\ \times \left\{\begin{array}{l}{\left({\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)\right)}^2\\ {\left({\sin }^2\left(t{\langle H_s\rangle }^{+}\right){\sin }^2\left(t{\langle H_s\rangle }^{-}\right)+{\cos }^2\left(t{\langle H_s\rangle }^{+}\right){\cos }^2\left(t{\langle H_s\rangle }^{-}\right)\right)}^2\end{array}\right\}\\ +O\left(\frac{1}{N+1}\right),n=N+2,\cdots ,2N+3\end{array}$ (A28)

The constants $b_1$ , $b_2$ depend on the form of the functions ${\Theta }_1\left(s\right)$ and ${\Theta }_2\left(s\right)$ as well as $\langle s_0|W_s\left(0\right)|{s^{\prime }}_0\rangle$ .

Finally, we infer that since the terms in the curly brackets in Equation (A26) tend to zero as $N\to \infty$ , the first term in eq. (A1) is exact as $N\to \infty$ and corresponds to the sign solved time evolution of the density matrix.

Conflicts of Interest

The author declares no conflicts of interest regarding the publication of this paper.

References

[1]	Johnson, G.W. and Lapidus, M.L. (2000) The Feynman Integral and Feynman’s Operational Calculus. Oxford Mathematical Monographs, New York.
[2]	Thrapsaniotis, E.G. (2003) Monte Carlo Methods via the Use of the Central-Limit Theorem. Europhysics Letters, 63, 479-484. https://doi.org/10.1209/epl/i2003-00491-5
[3]	Thrapsaniotis, E.G. (2004) Path Integral Approach of Ionization by Ultrashort Laser Pulses. Physical Review A, 70, Article ID: 033410. https://doi.org/10.1103/PhysRevA.70.033410
[4]	Thrapsaniotis, E.G. (2006) Atom in an Overdamped Cavity: Influence Functional Approach. Journal of Modern Optics, 53, 1501-1522. https://doi.org/10.1080/09500340600572701
[5]	Thrapsaniotis, E.G. (2007) Coherent State Path Integral Methods and Central Limit Theorem. Physics Letters A, 365, 191-199. https://doi.org/10.1016/j.physleta.2007.01.011
[6]	Thrapsaniotis, E.G. (2008) Sign Solved Propagator of Two Coupled and Forced Oscillators with Coupling of the Form . Far East Journal of Dynamical Systems, 10, 23-46.
[7]	Thrapsaniotis, E.G. (2008) Path Integral Methods via the Use of the Central Limit Theorem and Application. Journal of Physics A: Mathematical and Theoretical, 41, Article ID: 205202. https://doi.org/10.1088/1751-8113/41/20/205202
[8]	Thrapsaniotis, E.G. (2008) Chapter 12. Path Integral Approach to the Interaction of One Active Electron Atoms with Ultrashort Squeezed Pulses. In: Nilsson, V.P., Ed., Photonic Research Developments, Nova Science Publishers, New York, 351-383.
[9]	Thrapsaniotis, E.G. (2011) Chapter 19. Path Integral Approach to the Interaction of One Active Electron Atoms with Ultrashort Squeezed Pulses. In: George, T.F., Letfullin, R.R. and Zhang, G., Eds., Perspectives in Theoretical Physics, Nova Science Publishers, New York, 249-272.
[10]	Thrapsaniotis, E.G. (2010) Sign Solved Propagator of a Two Electrons 2D Atom and Application. Far East Journal of Dynamical Systems, 14, 171-182.
[11]	Thrapsaniotis, E.G. (2010) Path Integral Approach to the Interaction of One Active Electron Atoms with Ultrashort Squeezed Pulses. International Journal of Theoretical Physics, Group Theory, and Nonlinear Optics, 13, 221-244.
[12]	Thrapsaniotis, E.G. (2022) Angular Decomposition and Sign Solved Propagator of the 2D and 3D Helium Atom Path Integrals. Journal of Applied Mathematics and Physics, 10, 538-557. https://doi.org/10.4236/jamp.2022.102040
[13]	Weiss, U. (2012) Quantum Dissipative Systems. 4th Edition, World Scientific, Singapore City. https://doi.org/10.1142/9789814374927
[14]	Muller, K., Luoma, K. and Strunz, W.T. (2020) Geometric Phase Gates in Dissipative Quantum Dynamics. Physical Review A, 102, Article ID: 032611. https://doi.org/10.1103/PhysRevA.102.032611
[15]	Makri, N. (1995) Numerical Path Integral Techniques for Long Time Dynamics of Quantum Dissipative Systems. Journal of Mathematical Physics, 36, 2430-2457. https://doi.org/10.1063/1.531046