From c32d75fed6594c4cf643488abf94eab425aa10ad Mon Sep 17 00:00:00 2001 From: Valentin Boettcher Date: Wed, 7 Sep 2022 14:23:47 +0200 Subject: [PATCH] update the analytical solution writeup --- src/analytical_solution.tex | 233 ++++++++++++++++---------------- src/flow.tex | 255 ++++++++++++++++++------------------ src/hops_tweak.tex | 67 ++++++++++ 3 files changed, 310 insertions(+), 245 deletions(-) diff --git a/src/analytical_solution.tex b/src/analytical_solution.tex index 55f9506..fd20c03 100644 --- a/src/analytical_solution.tex +++ b/src/analytical_solution.tex @@ -1,24 +1,27 @@ \chapter{An Analytical Solution for Quantum Brownian Motion Models} \label{chap:analytsol} -\begin{itemize} -\item \cref{chap:flow} is all nice and good, but need verification -\item many actual numerical approximations to control besides theory - verification -\end{itemize} + +The results of \cref{chap:flow} are promising from a numerical +perspective remain to be verified. Previous +work~\cite{Hartmann2017Dec,RichardDiss} has made it clear that the +reduced system dynamics, but it is an open question whether bath +related quantities can be calculated to a similar degree of accuracy. + +The best possible verification is the comparison with a soluble model, +ideally solved with a method completely different from the NMQSD. In +this chapter we will present the solution to the Heisenberg equations +of two quantum Brownian motion models for a single (\cref{sec:oneosc}) +and two baths (\cref{sec:twoosc}) and time independent Hamiltonians. + +These solutions will enable us to calculate the bath energy flow \(J\) +for one or two baths which will be compared with the numerical +solution in \cref{sec:hopsvsanalyt}. + \section{A Harmonic Oscillator coupled to a single Bath} \label{sec:oneosc} -\begin{itemize} -\item simplest model in the spirit of previous derivations -\item well known to be solvable, we will tease out the bath related - quantities in the heisenberg picture -\end{itemize} -The solution presented here is not entirely new, as the model is well -known. For reference see~\cite{Breuer2002Jun}. The trick with the -exponential expansion of the bath correlation function has been -arrived at independently, but may also be well known. - -The model is given by the quadratic hamiltonian +A simple quadratic model that is soluble~\cite{Breuer2002Jun} and of +the form \cref{eq:totalH} is given by \begin{equation} \label{eq:one_ho_hamiltonian} H = \frac{Ω}{4}\qty(p^2+q^2) + \frac{1}{2} q @@ -29,14 +32,13 @@ where \(a,a^†\) are the ladder operators of the harmonic oscillator, \(q=a+a^†\) and \(p=\frac{1}{\iu}\qty(a-a^†)\) so that \([q,p] = 2\iu\). -The Heisenberg equation yields +The Heisenberg equations for \cref{eq:one_ho_hamiltonian} \begin{align} \dot{q} &=Ω p \label{eq:qdot}\\ \dot{p} &= -Ω q - \int_0^t \Im[α_0(t-s)] q(s)\dd{s} + W(t) \label{eq:pdot} \\ \dot{b}_λ &= -\iu g_λ \frac{q}{2} - \iu\omega_λ b_λ \end{align} - with the operator noise \(W(t)=-\sum_λ \qty(g_λ^\ast b_λ(0) \eu^{-\iu\omega_λ t } + g_λ b_λ^†(0) @@ -79,9 +81,12 @@ leads to an algebraic formula \mathcal{L}\{G\}(z) = \qty(z-A + \mathcal{L}\{K\}(z))^{-1}. \end{equation} - \subsection{Solution} \label{sec:solution} +To solve \cref{eq:eqmotprop} and find the propagator \(G\), we have to +find an explicit expression for \cref{eq:galgebr}, a simple matrix +inversion, and then apply the inverse transformation. + We observe that \begin{equation} \label{eq:mdef} @@ -95,7 +100,7 @@ and therefore \mqty(z & Ω \\ -(Ω + \mathcal{L}\{\Im[α_0]\}(z)) & z). \end{equation} From this we can conclude that \(G_{11}=G_{22}\). -Because \(\ev{W(s)}=0\) for thermal initial states of the bath we have +Because \(\ev{W(s)}=0\) holds for thermal initial states of the bath we have \begin{equation} \label{eq:meanvals} \mqty(\ev{q(t)}\\ \ev{p(t)}) = G(t)\mqty(\ev{q(0)}\\ \ev{p(0)}). @@ -107,11 +112,16 @@ Knowing this, we can deduce from \(\ev{\dot{q}}= Ω \ev{p}\) that \end{align} These relations are true independent of the initial state of the system. It therefore suffices if we concern ourselves with -\(G_{12}\). We nevertheless continue in full generality. +\(G_{12}\). + +We nevertheless continue in full generality and approach the inverse +Laplace transformation by expanding the BCF in terms of functions that +have a simple Laplace transform. As we also use an exponential +expansion in HOPS and are only interested in finite times, we may +choose \(α_0(t)=\sum_{n=1}^N G_n \eu^{-W_n t - \i \varphi_n}\) with +\(W_n=\gamma_n + \i\delta_n\) and +\(G_n, \varphi_n, \gamma_n,\delta_n\in\RR\) for \(t\geq 0\). -Assume that \(α_0(t)=\sum_{n=1}^N G_n \eu^{-W_n t - \i - \varphi_n}\) with \(W_n=\gamma_n + \i\delta_n\) and \(G_n, -\varphi_n, \gamma_n,\delta_n\in\RR\) for \(t\geq 0\). This leads to a mathematically simple expression for the Laplace transform \begin{equation} @@ -166,10 +176,10 @@ of degree \(2(N+1)\) where we \textbf{assume that there are no roots with multiplicity greater than one}. With this we can now calculate the inverse laplace transform of expressions of the form \(\frac{f_0(z)g(z)}{p(z)}\) where \(g(z)\) is -any holonome function so that \(\frac{f_0(z)g(z)}{p(z)} \eu^{z\cdot +any holonomic function so that \(\frac{f_0(z)g(z)}{p(z)} \eu^{z\cdot t}\) falls off fast enough for \(t\geq 0\), \(\Re(z)>\max_l{\Re(\tilde{z}_l)}=\Delta\) and \(\Re(z) \rightarrow --\infty\). With this we can close the contour of the inverse Laplace +-\infty\). Now the contour of the inverse Laplace transform \begin{equation} \label{eq:invlap} @@ -177,7 +187,7 @@ transform \frac{1}{2\pi\i}\int_{\Delta - \i\infty}^{\Delta + \i\infty} \frac{f_0(z)g(z)}{p(z)} \eu^{z\cdot t}\dd{z} \end{equation} -to the left to obtain +can be closed to the left to obtain \begin{equation} \label{eq:simpleinvtrans} \mathcal{L}^{-1}\qty{\frac{f_0(z)g(z)}{p(z)}}(t) @@ -192,7 +202,8 @@ all our purposes. For completeness we give p'(z) = 2\mu\sum_{k=1}^{N+1}\qty[(z-\Re(\tilde{z}_k))\prod_{\substack{n=1\\ n\neq k}}^{N+1}(z-\tilde{z}_n)(z-\tilde{z}^\ast_n)]. \end{equation} -We can immediately conclude that all elements of \(G\) are sums of + +It can be immediately concluded that all elements of \(G\) are sums of exponentials, just like the BCF. In particular \begin{equation} \label{eq:gfinal} @@ -201,7 +212,17 @@ exponentials, just like the BCF. In particular \end{equation} with \(R_l={f_0(\tilde{z}_l)}/{p'(\tilde{z}_l)}\). +It may be noted that this solution does not contain any notion of +temperature, as we are working in the Heisenberg picture. + \subsubsection{Negative Times} +For completeness, it may be of interest to find a solution for +negative times. This solution is relatively unphysical, as the initial +condition of a product state plays a pivotal role in open system +dynamics~\cite{Rivas2012}. Therefore a system that starts out in some +entangled state just to reach the perfect product state at \(t=0\) is +not something that is likely to be applicable to physical questions. + The solution detailed above is only valid for positive times. Because we strive to employ the same formalism again for negative times, we will concern ourselves with the transformed quantities @@ -236,13 +257,20 @@ and \subsection{Applications} \label{sec:applications} -\begin{itemize} -\item the actual ``meat'' that is of interest to verify HOPS -\end{itemize} +Having found an expression for \(G\), we have in principle solved the +model. It remains however to apply that solution in a way that is +contributing towards our goal of validating the results of +\cref{chap:flow}. + + Knowing \(G\) and \(α\), we can calculate all observables of the -system. Simple closed form expressions of sums of exponentials can be -obtained by using an exponential expansions of \(α\). Throughout, -we assume a thermal bath initial state so that \(\ev{W(t)}=0\). +system with the ultimate goal of finding an expression for +\(J=-∂_{t}\ev{H_{B}}\). Simple closed form expressions of sums of +exponentials will be obtained by using an exponential expansion of +\(α\). + +Throughout, we assume a thermal bath initial state so that +\(\ev{W(t)}=0\). \subsubsection{Correlation Functions} \label{sec:correl} @@ -259,7 +287,7 @@ We proceed to calculate \(\ev{q(t)q(s)}\). For brevity we set \end{aligned} \end{equation} -For a pure harmonic oscillator initial state \(\ket{n}\) we have +For a Fock type initial state \(\ket{n}\) we have \begin{equation} \label{eq:hoexp} \begin{aligned} @@ -326,73 +354,17 @@ For further evaluation of \cref{eq:bathderiv_1} we have to calculate \biggr] \end{aligned} \\ - &= ∫_0^t\dd{r} g_1(t-r)\qty[g_2(t-r) + g_3(t,r)] = Λ_1(t) + Λ_2(t) + &= ∫_0^t\dd{r} g_1(t-r)\qty[g_2(t-r) + g_3(t,r)] = Λ_1(t) + Λ_2(t), \end{aligned} \end{equation} +where \(Λ(t,s)\) was defined in \cref{eq:qcorrel}. + This expression now only uses \(α(t)\) for \(t\geq 0\) so that we can -once again employ the exponential expansion for \(α\). In fact, all -other quantities in \cref{eq:lambdafold} have exponential expansion so -that we can now define\footnote{Note that this is inconsistent with - \cref{sec:solution}.} -\begin{equation} - \label{eq:expansions} - \begin{aligned} - α_0&=∑_k U_k\eu^{-Q_k t} & \dot{α}_0&=∑_k P_k\eu^{-L_k t} & α(t) - &= ∑_nG_n\eu^{-W_n t} \\ - A(t) &= ∑_l A_l\eu^{-C_l t} & B(t) &= ∑_l B_l\eu^{-C_l t}. - \end{aligned} -\end{equation} +once again employ the exponential expansion for \(α\). -With this we can calculate, -\begin{align} - \label{eq:lambdaintegrals} - ∫_r^t\dd{s}B(s-r)\dot{α}_0(t-s) - &=\sum_{m,k}\underbrace{\frac{B_mP_k}{L_k-C_m}}_{\equiv - Γ^1_{mk}}\qty[\eu^{-C_m(t-r)}-\eu^{-L_k(t-r)}]=g_1(t-r)\\ - ∫_0^{t-r}\dd{u}B(t-r-u)α(u) - &=\sum_{n,l}\underbrace{\frac{B_nG_l}{C_n-W_l}}_{\equiv - Γ^2_{nl}}\qty[\eu^{-W_l(t-r)}-\eu^{-C_n(t-r)}]=g_2(t-r)\\ - ∫_0^{r}\dd{u}B(t-r+u)α^\ast(u) - &=\sum_{n,l}\underbrace{\frac{B_nG_l^\ast}{C_n+W_l^\ast}}_{\equiv - Γ^3_{nl}}\qty[\eu^{-C_n(t-r)}-\eu^{-W_l^\ast r-C_n t}]=g_3(t,r) -\end{align} -and -\begin{align} - \label{eq:finalsummands} - Λ_1(t)&= ∑_{m,k,n,l}Γ^1_{mk}Γ^2_{nl}\qty[\frac{1-\eu^{-(C_m+W_l)t}}{C_m+W_l} - - - \frac{1-\eu^{-(C_m+C_n)t}}{C_m+C_n}- - \frac{1-\eu^{-(L_k+W_l)t}}{L_k+W_l} - + - \frac{1-\eu^{-(L_k+C_n)t}}{L_k+C_n}]\\ - Λ_2(t)&= - \begin{aligned}[t] - ∑_{m,k,n,l}Γ^1_{mk}Γ^3_{nl}\biggl[\frac{1-\eu^{-(C_m+C_n)t}}{C_m+C_n} - &-\frac{1-\eu^{-(L_k+C_n)t}}{L_k+C_n} - \\&-\frac{\eu^{-(C_n+W_l^\ast)t}-\eu^{-(C_m+C_n)t}}{C_m-W_l^\ast} - +\frac{\eu^{-(C_n+W_l^\ast)t}-\eu^{-(L_k+C_n)t}}{L_k-W^\ast_l}\biggr] - \end{aligned} -\end{align} - -Also required for \cref{eq:bathderiv_1} are -\begin{align} - \label{eq:ABconv} - ∫_0^t\dd{s}A(s)\dot{α}_0(t-s) &= ∑_{n,m}\underbrace{\frac{A_nP_m}{L_m-C_n}}_{\equiv - Γ^A_{nm}}\qty[\eu^{-C_n t}-\eu^{-L_m t}]\\ - ∫_0^t\dd{s}B(s)\dot{α}_0(t-s) &= ∑_{n,m}Γ^1_{nm}\qty[\eu^{-C_n t}-\eu^{-L_m t}] -\end{align} -and -\begin{multline} - \label{eq:nonzerotemplim} - ∫_0^t\dd{s}A(s)\qty(α(s)-α_0(s)) =\\ - ∑_{m,n}\frac{A_nG_m}{C_n+W_m}\qty(1-\eu^{-(C_n+W_m)t}) - ∑_{m,n}\frac{A_nU_m}{C_n+Q_m}\qty(1-\eu^{-(C_n+Q_m)t}). -\end{multline} - -This concludes the calculation. A possible measure of simplification -would be to write \cref{eq:bathderiv_1} as a sum of exponentials and -give explicit expressions for the coefficients and exponents. This is -not required for now. Code implementing this can be found under -\url{https://github.com/vale981/hopsflow}. +We will arrive at expressions that are weighted sums of exponentials +whose detailed calculation is quite tedious and can be found +in\cref{sec:explicit_flow}. \section{Two coupled Harmonic Oscillators coupled to two Baths}% \label{sec:twoosc} @@ -403,16 +375,19 @@ not required for now. Code implementing this can be found under not trivially additive and we can reuse the method from \cref{sec:oneosc} without alteration \end{itemize} +As we would like to verify our method also for more than one bath, a +model with two baths is required. The considerations of~\cref{sec:oneosc} can be straight forwardly generalized to the case of two coupled oscillators coupled in turn to -a bath each. +a bath each. This construction is chosen so that the previous results +can be reused and the coupling to the baths is not trivial. We will not give explicit formulas for the results in terms of sums of exponentials, as they are quite extensive and easily obtained via the use of a computer algebra system or the aforementioned code. -The model is again given by a quadratic hamiltonian +The model is again given by a quadratic Hamiltonian \begin{equation} \label{eq:hamiltonian_two_bath} \begin{aligned} @@ -427,28 +402,27 @@ The \(a_i^{(i)},a_i^{†}\) are the ladder operators of the harmonic oscillators and \(q_i=a_i+a_i^†\) and \(p=\frac{1}{\iu}\qty(a_i-a_i^†)\) so that \([q_i,p_j] = 2\iuδ_{ij}\) and \([q_i,q_j] = [p_i,p_j] = 0\). -\subsection{Equations of Motion} -\label{sec:eqmot_two} -The Heisenberg equation yields +The Heisenberg equations for \cref{eq:hamiltonian_two_bath} are \begin{align} \dot{q}_i &= Ω_i p_i \label{eq:qidot}\\ \dot{p}_i &= -(Ω_i+γ) q_i - \int_0^t \Im[α^{(i)}_0(t-s)] q_i(s)\dd{s} + W_i(t) \label{eq:pidot} \\ - \dot{b}^{(i)}_λ &= -\iu g^{(i)}_λ \frac{q_i}{2} - \iu\omega^{(i)}_λ b^{(i)}_λ + \dot{b}^{(i)}_λ &= -\iu g^{(i)}_λ \frac{q_i}{2} - \iu\omega^{(i)}_λ b^{(i)}_λ, \end{align} - with the operator noise \(W_i(t)=-\sum_λ \qty(g_λ^{(i),\ast} b^{(i)}_λ(0) \eu^{-\iu\omega^{(i)}_λ t } + g_λ^{(i)} b_λ^{(i),†}(0) \eu^{\iu\omega^{(i)}_λ t })\) satisfying \(\ev{W_i(s)}=0\) and -\(\ev{W_i(t)W_j(s)}=δ_{ij}α^{(i)}(t-s)\). We introduced \(α^{(i)}_0 -\equiv \eval{α^{(i)}}_{T=0}\). +\(\ev{W_i(t)W_j(s)}=δ_{ij}α^{(i)}(t-s)\). We introduced +\(α^{(i)}_0 \equiv \eval{α^{(i)}}_{T=0}\). We have given most quantities an extra index and accounted for the coupling between the two oscillators. Apart from this, the equations -of motion have the same structure as in \cref{seq:eqmot}. +of motion have the same structure as in \cref{sec:oneosc}. -Again, we obtain +\subsection{Solution} +\label{sec:eqmot_two} +With the same general program as before, we we first obtain \begin{equation} \label{eq:bsoltwo} b^{(i)}_λ(t) = b^{(i)}_λ(0) \eu^{-\iu ω^{(i)}_λ t} - \frac{\iu g^{(i)}_λ}{2}∫_0^t @@ -501,7 +475,7 @@ This can be done easily\footnote{We have use a computer algebra \begin{equation} \label{eq:minvtwo} M^{-1}(z) = \frac{1}{\det[M](z)} \tilde{M}(z) - \end{equation} +\end{equation} where \(\tilde{M}\) is a matrix containing only polynomials of \(z\) and of the Laplace transforms of the bath correlation functions. @@ -509,15 +483,16 @@ The numerator is \begin{equation} \label{eq:numerator} \begin{aligned} - \det[M](z)=a& b \Lambda \Omega +a \left(\gamma + \det[M](z)=a(z)& b(z) \Lambda \Omega +a(z) \left(\gamma \Lambda \Omega +\Lambda ^2 \Omega +\Omega z^2\right) - +b + +b(z) \left(\gamma \Lambda \Omega +\Lambda \Omega ^2+\Lambda z^2\right)\\ &+\gamma \Lambda ^2 \Omega +\gamma \Lambda \Omega ^2+\Lambda ^2 \Omega ^2+\gamma \Lambda z^2+\gamma \Omega z^2+\Lambda ^2 z^2+\Omega ^2 z^2+z^4, \end{aligned} \end{equation} -where \(a=\mathcal{L}\{\Im[α^{(1)}_0]\}\) and \(b=\mathcal{L}\{\Im[α^{(2)}_0]\}\). +where \(a(z)=\mathcal{L}\{\Im[α^{(1)}_0]\}(z)\) and +\(b(z)=\mathcal{L}\{\Im[α^{(2)}_0]\}(z)\). Using the same approach as in \cref{sec:solution}, we arrive at an expression similar to \cref{eq:prefactorrational} for @@ -534,17 +509,17 @@ exponentials G(t) = \sum_{l=1}^{2+N_1+N_2}\qty[R_l \tilde{M}(\tilde{z}_l)\eu^{\tilde{z}_l \cdot t} + \cc] \end{equation} -with \(R_l={f_0(\tilde{z}_l)}/{p'(\tilde{z}_l)}\). +with \(R_l={f_0(\tilde{z}_l)}/{p'(\tilde{z}_l)}\) as defined in +\cref{sec:solution}. \subsection{Applications} \subsubsection{Correlation Functions} \label{sec:correltwo} -\begin{itemize} -\item same game as in \cref{sec:applications} -\end{itemize} We can now proceed to calculate the correlation functions \(C(t,s) = \ev{x_i(t)x_j(s)}\) where the \(x_i\) are the phase space operators -of the two harmonic oscillators. +of the two harmonic oscillators. This will enable use to calculate the +system energies of the two oscillators (omitted here) and again the +bath energy flows of the two baths. We find \begin{equation} @@ -556,6 +531,18 @@ We find The matrix \(Θ_{ij}\) contains the bath-induced correlations and can be calculated as in the single-oscillator case. +For two oscillators that are initially in Fock states +\(\ket{ψ_{0}}=\ket{n}\otimes\ket{m}\) we have +\begin{equation} + \label{eq:initial_corr} + C(0,0) = + \begin{pmatrix} + 1 + 2 n & \i & 0 & 0 \\ + -\i & 1 + 2 n & 0 & 0 \\ + 0 & 0 & 1 + 2 m & \i \\ + 0 & 0 & -\iu & 1 + 2m + \end{pmatrix}. +\end{equation} \subsubsection{Bath Enery Derivative} \label{sec:bathflowtwo} @@ -568,6 +555,10 @@ Similar to the calculations in \cref{sec:bathflow} we find ∑_{k=1,2}G_{2n-1,2k}∂_s\qty(α^{(k)}(s)-α_0^{(k)}(s)). \end{equation} -This can be evaluated using the exponential expansions and yields -yet another sum of exponentials. The steady state flow can then be -found be setting all exponentials to zero. +This can be evaluated using the exponential expansions and yields yet +another sum of exponentials. The steady state flow can then be found +be setting all exponentials to zero, although care has to be taken, as +an exponential fit of the BCF may be only valid for finite times. + +This concludes the calculation. Python code implementing the solution +can be found under \url{https://github.com/vale981/hopsflow}. diff --git a/src/flow.tex b/src/flow.tex index f0d8e49..72da5a0 100644 --- a/src/flow.tex +++ b/src/flow.tex @@ -159,7 +159,7 @@ state. This can be generalized to any BCF that is a sum of exponentials. Interestingly one finds that \begin{equation} \label{eq:alternative} - \ev{L∂_t B^†(t)} = \i∫\frac{\dd[2]{z}}{\pi^N} + \ev{L∂_t B^†(t)} = \i \mathcal{M}_{η^\ast} \dot{η}_t^\ast \mel{\psi(η,t)}{L}{\psi(η^\ast,t)}. \end{equation} This expression is undesirable as it does not exist for all bath @@ -201,7 +201,8 @@ Following the usual derivation of the nonlinear NMQSD we write \end{aligned} \end{equation} where \(\tilde{z}_{\lambda}^{*}(t)=z_{\lambda}^{*}+\i g_{\lambda} ∫_{0}^{t} -\dd{s} \eu^{-\i ω_{\lambda} s}\ev{L^†}_{s}\). +\dd{s} \eu^{-\i ω_{\lambda} s}\ev{L^†}_{s}\) and \(\ev{L^\dag}_{t}=ψ(\tilde{η}_{t}^\ast)_{t}^\dag L^\dag +ψ(\tilde{η}_{t}^\ast)_{t}\) as in \cref{sec:nmqsd_basics}. It has to be shown now, that the term \({\braket{\psi}{\tilde{\vb{z}}(t)}\!\braket{\tilde{\vb{z}}(t)}{\psi}}\) @@ -438,26 +439,102 @@ This is an expression that we can easily evaluate with the HOPS method. We will however refrain from doing so, as it turns out in \cref{sec:hopsvsanalyt} that consistent results can be obtained using the derivative of the stochastic process \(ξ\). +\section{Generalization to Multiple Baths} +\label{sec:multibath} +Another requirement for thermodynamic application is the ability to +couple to multiple baths of possibly different structure and +temperature. + +Due to the structure of the NMQSD and HOPS the results above can be +generalized in straight-forward manner to models of the form +\begin{equation} + \label{eq:multimodel} + H = H_\sys + ∑_{n=1}^N \qty[H_\bath\nth + \qty(L_n^†B_n + \hc)], +\end{equation} +where \(N\) is the number of baths, \(H_\sys\) is the system +Hamiltonian, \(H_\bath\nth = ∑_λω_λ\nth a_λ^{(n),†}a_λ\nth\), +\(B_n=∑_{λ} g_λ\nth a_λ\nth\) and the \(L_n={(\vb{L})}_n\) are +arbitrary operators in the system Hilbert space. + +Note that this models a situation where each bath couples with the +system through exactly one spectral density and is therefore not fully +general. + +We refer to \cref{sec:hops_multibath} for an review of the NMQSD +theory and HOPS method for multiple baths. + +Because the bath energy change is being calculated directly and not +through energy conservation as in~\cite{Kato2016Dec}, we find +\begin{equation} + \label{eq:general_n_flow} + J_n=-\dv{\ev{H_\bath^{(n)}}}{t} = \iu\ev{[H_\bath^{(n)}, + H_\inter^{(n)}]} +\end{equation} +regardless of the (non-) commutativity\footnote{For example, the + three-level model used in \cite{Uzdin2015Sep,Klatzow2019Mar} has + non-commuting couplings.} of the interaction +Hamiltonians. Therefore, we can apply the formalism of the previous +sections almost unchanged, by taking care that all quantities involved +in the expression of \(J_n\) refer to the \(n\)th bath denoted by sub +and superscripts. + +This can be achieved by making the replacements +\begin{equation} + \label{eq:replacements} + \begin{aligned} + D_t \rightarrow D_t^{(n)} &\equiv + ∫_0^t\dd{s}α_n(t-s)\fdv{η^\ast_n(s)} \\ + ξ(t) \rightarrow ξ_n(t)&\equiv∑_{\lambda} g^{(n)}_{\lambda} + y_{\lambda} \eu^{-\mathrm{i} ω^{(n)}_{\lambda} t} + \end{aligned} +\end{equation} +in the previous sections, where the quantities involved are as in +\cref{sec:hops_multibath} and \cref{eq:xiproc}. + +Foreshadowing slightly it may be states that in the light of +\cref{sec:general_obs} it might be an interesting question what impact +mixed bath hierarchy states have. For a cyclic machine with long +strokes, where only one bath is coupled to the system at a time, it +might be efficient to truncate the hierarchy in a way that discards +mixed bath states more readily than single bath hierarchy states as +the correlations between the baths are expected to be small. + +\section{Generalization to Time Dependent Hamiltonians} +\label{sec:timedep} +To extract energy from a quantum thermal machine without an explicit +work reservoir, external modulation is required. + +The above discussion is based on the model \cref{eq:multimodel} which +did not include explicit time modulations of \(H_\sys\) or \(L\). As +we did not calculate any explicit time derivatives of those two +operators, the results of the previous sections remain valid when we +substitute \(H_\sys\rightarrow H_\sys(t)\) and \(L\rightarrow L(t)\). + +For the total power we find +\begin{equation} + \label{eq:power} + \dv{\ev{H}}{t} = \ev{\pdv{H_\inter}{t}} + \ev{\pdv{H_\sys}{t}}, +\end{equation} +which can be evaluated as we will find in \cref{sec:intener} by +replacing \(L(t)\) with \(\dot{L}(t)\). \section{General Collective Bath Observables} \label{sec:general_obs} -\begin{itemize} -\item now we have all the basic tools ready we can discuss what the - most general observables are that we can calculate -\item elucidates the meaning of the hierarchy states -\end{itemize} -The results obtained in \cref{sec:flow_lin,sec:nonlin_flow,sec:lin_finite} -can be generalized to calculate expectation values (and thus moments) -of arbitrary observables of the form +Now that we have introduced the formalism using the example of the +bath energy flow \(J\) in +\cref{sec:flow_lin,sec:nonlin_flow,sec:lin_finite,sec:multibath,sec:timedep}, +we may proceed to more general observables of the form can be +generalized to calculate expectation values (and thus moments) of +arbitrary observables of the form \begin{equation} \label{eq:collective_obs} O = f(B^†, B) = ∑_{α}F_α\qty(B^†)^{α_1}B^{α_2} \end{equation} -where \(α\) is a multi-index, \(B\) is as in~\cref{eq:totalH} and the -\(F_α\) are general observables acting on the system only. Note that -\(f\) is already normal-ordered. We will restrict the discussion to -the case of a single bath, as the generalization to multiple baths is -straight-forward. +where \(α\) is a two-dimensional multi-index, \(B\) is as +in~\cref{eq:totalH} and the \(F_α\) are general observables acting on +the system only. Note that \(f\) is already normal-ordered. We will +restrict the discussion to the case of a single bath, as the +generalization to multiple baths is straight-forward. To evaluate \(\ev{O}\) we have to find the value of \(\ev{\qty(B^†)^a B^b}\). This can be achieved by interjecting the @@ -475,12 +552,12 @@ For zero temperature, we find following the procedures of \sqrt{\frac{G^{\vb{k}}}{\vb{k}!}}ψ^{\vb{k}}\\ \label{eq:bdagmel}\mel{ψ}{\qty(B^†)^a}{z} &= \begin{aligned}[t] - \qty(\mel{z}{B^a}{ψ})^†&= \qty((-\iu D_t)^b\ket{ψ(η^\ast,t)})^\dag\\ + \qty(\mel{z}{B^a}{ψ})^†&= \qty((-\iu D_t)^a\ket{ψ(η^\ast,t)})^\dag\\ &= (\iu)^a∑_{\abs{\vb{k}}=a}\binom{a}{\vb{k}} (-\iu)^{\vb{k}} - \sqrt{\frac{\qty(G^{\vb{k}})^\ast}{\vb{k}!}}\qty(ψ^{\vb{k}})^\dag + \sqrt{\frac{\qty(G^{\vb{k}})^\ast}{\vb{k}!}}\qty(ψ^{\vb{k}})^\dag, \end{aligned} \end{align} -in ``fock-space'' normalization where \(\vb{k}! = k_1!k_2!\ldots\) and +where \(\vb{k}! = k_1!k_2!\ldots\) and \(G^{\vb{k}}=G_1^{k_1}G_2^{k_2}\ldots\) following the usual conventions of multi-indices. Thus, expressions involving the bath operator \(B\) to the \(b\)th power lead to expressions involving the @@ -513,7 +590,7 @@ involving the HOPS hierarchy states this reduces to dividing by the norm of the zeroth hierarchy state. The generalization to multiple baths may be performed in the same -manner as will be discussed in \cref{sec:multibath}. This allows to +manner as was discussed in \cref{sec:multibath}. This allows to calculate the expectation value involving multiple bath operators \(B^{(n)}\). Interestingly, the generalization of \cref{eq:bmel} to multiple baths immediately links hierarchy states of the form @@ -529,7 +606,7 @@ and inserting the coherent state resolution of unity we find terms of the form \begin{equation} \label{eq:with_process} - \mel{z}{\qty(B^\dag)^b}{ψ} \sim \qty(η^\ast(t))^b\ket{ψ(η^\ast,t)}. + \mel{z}{\qty(B^\dag)^b}{ψ} \sim \qty(η^\ast_{t})^b\ket{ψ(η^\ast,t)}. \end{equation} The corresponding version of~\cref{eq:f_ex_zero} would only depend on the zeroth order state and the stochastic processes. It has been @@ -541,6 +618,18 @@ their average dynamics whereas the stochastic process fluctuates around zero and does not contain much information about the actual dynamics. +% The process \(\pqty{η^\ast_{t}}^{b}\) has the autocorrelation function +% \(\ev{\pqty{η_{t}}^{b}\pqty{η^\ast_{s}}^{b}}=b! (α(t-s))^{b}\). Now +% for a decaying function +% \begin{equation} +% \label{eq:bcf_exponentiated} +% \pqty{\frac{α(τ)}{α(0)}}^{b}\xrightarrow{b\to ∞} +% \begin{cases} +% 1 & τ = 0 \\ +% 0 & τ \neq 0, +% \end{cases} +% \end{equation} +% so that we end up with a process that is some approximation of white noise. Also, this alternative method could be used convergence and consistency check, as expressions of the form~\cref{eq:with_process} only involve the hierarchy cutoff and the exponential expansion of the @@ -553,9 +642,10 @@ BCF in an indirect manner. \item can also be used to quantify how ``strong'' the coupling is \item simple application of the above formalism, even simpler than flow \end{itemize} -A simple application of the formalism discussed -in~\cref{sec:general_obs} is the expectation value of the interaction -Hamiltonian. +To access all contributions to the total energy \(\ev{H}\) a way to +calculate the expectation value of the interaction energy +\(\ev{H_{\inter}}\) is required. This is a application of the +formalism discussed in~\cref{sec:general_obs}. For zero temperature and the nonlinear method we arrive at \begin{equation} @@ -579,7 +669,7 @@ In HOPS terms \cref{eq:intexp} corresponds to For nonzero temperature an extra term \begin{equation} \label{eq:interexptherm} - \mathcal{M}_{\tilde{η}^\ast}\frac{\mel{\psi(\tilde{η},t)}{L^†ξ(t)}{\psi(\tilde{η}^\ast,t)}}{\braket{\psi(\tilde{η},t)}{\psi(\tilde{η}^\ast,t)}} + \mathcal{M}_{\tilde{η}^\ast,ξ}{\mel{\psi(\tilde{η},t)}{L^†ξ(t)}{\psi(\tilde{η}^\ast,t)}} + \cc \end{equation} has to be added to \cref{eq:intexp}, where \(ξ\) is the thermal @@ -587,10 +677,12 @@ stochastic process. \subsection{Higher Orders of the Coupling Hamiltonian} \label{sec:higher_order_coupling} -\begin{itemize} -\item a slight dive into the meaning and importance of the hierarchy - states for convergence -\end{itemize} +In this section, the question of how many hierarchy orders have to be +included in the simulation to consistently calculate the expectation +value of powers of the interaction Hamiltonian. Being nonessential for +the understanding of the rest of the work, this section may be +skipped. + For self adjoint coupling operators \(L=L^\dag\) we can use Wick's theorem to find a normally ordered expression for \(H_\inter^n=L^n(B^\dag + B)^n\). @@ -687,100 +779,15 @@ depends strongly on hierarchy states of order \(\sim \sqrt{n}\). \begin{figure}[h] \centering \plot{interaction_orders/k_weights} - \caption{\label{fig:kdist}The unnormalized distribution of the coefficients in - \cref{eq:interactionnormal} with respect to \(k\) for different - \(n\). As a particular \(k\) can appear multiple times in the sum, - only the maximal coefficient for each \(k\) is being considered in - the lines with the triangle markers. The maximum of this - distribution is given by \cref{eq:finalk}. The lines with the - circle markers show the full distribution. The dotted lines - correspond to binomial distributions centered at \(k_m\).} + \caption{\label{fig:kdist}The unnormalized distribution of the + coefficients in \cref{eq:interactionnormal} with respect to \(k\) + for different \(n\). As a particular \(k\) can appear multiple + times in the sum, only the maximal coefficient for each \(k\) is + being considered in the lines with the triangle markers. The + maximum of this distribution is given by \cref{eq:finalk}. The + lines with the circle markers show the full distribution and the + lines with the triangle markers show only the normalized + distribution of the maximum of \(l! 2^l k! (n-2l-k)!\) over \(l\). + The dotted lines correspond to binomial distributions centered at + \(k_m\).} \end{figure} - -\section{Generalization to Multiple Baths} -\label{sec:multibath} -\begin{itemize} -\item for thermodynamic machines: need multiple baths -\item due to the structure of NMQSD/HOPS generalization is - straightforward -\end{itemize} - -The results above can be generalized in straight-forward manner to -models of the form -\begin{equation} - \label{eq:multimodel} - H = H_\sys + ∑_{n=1}^N \qty[H_\bath\nth + \qty(L_n^†B_n + \hc)], -\end{equation} -where \(N\) is the number of baths, \(H_\sys\) is the (possibly time -dependent) system Hamiltonian, -\(H_\bath\nth = ∑_λω_λ\nth a_λ^{(n),†}a_λ\nth\), -\(B_n=∑_{λ} g_λ\nth a_λ\nth\) and the \(L_n={(\vb{L})}_n\) are -arbitrary operators in the system Hilbert space (again possibly time -dependent). This models a situation where each bath couples with the -system through exactly one spectral density and is therefore not fully -general. - -We refer to \cref{sec:hops_multibath} for an review of the NMQSD -theory and HOPS method for multiple baths. - -Because the bath energy change is being calculated directly and not -through energy conservation as in~\cite{Kato2016Dec}, we find -\begin{equation} - \label{eq:general_n_flow} - J_n=-\dv{\ev{H_\bath^{(n)}}}{t} = \iu\ev{[H_\bath^{(n)}, - H_\inter^{(n)}]} -\end{equation} -regardless of the (non-) commutativity\footnote{For example, the - three-level model used in \cite{Uzdin2015Sep,Klatzow2019Mar} has - non-commuting couplings.} of the interaction -Hamiltonians. Therefore, we can apply the formalism of the previous -sections almost unchanged, by taking care that all quantities involved -in the expression of \(J_n\) refer to the \(n\)th bath denoted by sub -and superscripts. - -This can be achieved by making the replacements -\begin{equation} - \label{eq:replacements} - \begin{aligned} - D_t \rightarrow D_t^{(n)} &\equiv - ∫_0^t\dd{s}α_n(t-s)\fdv{η^\ast_n(s)} \\ - ξ(t) \rightarrow ξ_n(t)&\equiv∑_{\lambda} g^{(n)}_{\lambda} - y_{\lambda} \eu^{-\mathrm{i} ω^{(n)}_{\lambda} t} - \end{aligned} -\end{equation} -in the previous sections, where the quantities involved are as in -\cref{sec:hops_multibath} and \cref{eq:xiproc}. - -In the light of \cref{sec:general_obs} it might be an interesting -question what impact mixed bath hierarchy states have. For a cyclic -machine with long strokes, where only one bath is coupled to the -system at a time, it might be efficient to truncate the hierarchy in a -way that discards mixed bath states more readily than single bath -hierarchy states as the correlations between the baths are expected to -be small. - -\section{Generalization to Time Dependent Hamiltonians} -\label{sec:timedep} -\begin{itemize} -\item second important ingredient: external modulation, at least for - our purposes here -\item again here -> no magic, straight forward -\end{itemize} - -The above discussion is based on the model \cref{eq:totalH} which did -not include explicit time modulations of \(H_\sys\) or \(L\). As we -did not calculate any explicit time derivatives of those two -operators, the results of the previous sections remain valid when we -substitute -\begin{align} - \label{eq:timedepsusbs} - H_\sys&\rightarrow H_\sys(t) & L\rightarrow L(t). -\end{align} - -For the total power we find -\begin{equation} - \label{eq:power} - \dv{\ev{H}}{t} = \ev{\pdv{H_\inter}{t}} + \ev{\pdv{H_\sys}{t}}, -\end{equation} -which can be evaluated as in \cref{sec:intener} by replacing \(L(t)\) -with \(\dot{L}(t)\). diff --git a/src/hops_tweak.tex b/src/hops_tweak.tex index 90df694..a0616a3 100644 --- a/src/hops_tweak.tex +++ b/src/hops_tweak.tex @@ -437,3 +437,70 @@ The polynomial expressions for the smoothstep functions are &{\text{if }}1\leq x\\ \end{cases}}. \end{equation} + +\subsection{Explicit Expressions for the Bath Energy Flow of the QBM + Model} +\label{sec:explicit_flow} +Here we detail the rest of the calculations omitted in +\cref{sec:bathflow}. + +All quantities in \cref{eq:lambdafold} have exponential expansion so +that we can now define\footnote{Note that this is inconsistent with + \cref{sec:solution}.} +\begin{equation} + \label{eq:expansions} + \begin{aligned} + α_0&=∑_k U_k\eu^{-Q_k t} & \dot{α}_0&=∑_k P_k\eu^{-L_k t} & α(t) + &= ∑_nG_n\eu^{-W_n t} \\ + A(t) &= ∑_l A_l\eu^{-C_l t} & B(t) &= ∑_l B_l\eu^{-C_l t}. + \end{aligned} +\end{equation} + +With this we can calculate, +\begin{align} + \label{eq:lambdaintegrals} + ∫_r^t\dd{s}B(s-r)\dot{α}_0(t-s) + &=\sum_{m,k}\underbrace{\frac{B_mP_k}{L_k-C_m}}_{\equiv + Γ^1_{mk}}\qty[\eu^{-C_m(t-r)}-\eu^{-L_k(t-r)}]=g_1(t-r)\\ + ∫_0^{t-r}\dd{u}B(t-r-u)α(u) + &=\sum_{n,l}\underbrace{\frac{B_nG_l}{C_n-W_l}}_{\equiv + Γ^2_{nl}}\qty[\eu^{-W_l(t-r)}-\eu^{-C_n(t-r)}]=g_2(t-r)\\ + ∫_0^{r}\dd{u}B(t-r+u)α^\ast(u) + &=\sum_{n,l}\underbrace{\frac{B_nG_l^\ast}{C_n+W_l^\ast}}_{\equiv + Γ^3_{nl}}\qty[\eu^{-C_n(t-r)}-\eu^{-W_l^\ast r-C_n t}]=g_3(t,r) +\end{align} +and +\begin{align} + \label{eq:finalsummands} + Λ_1(t)& = + \begin{aligned}[t] + ∑_{m,k,n,l}Γ^1_{mk}Γ^2_{nl}\biggl[\frac{1-\eu^{-(C_m+W_l)t}}{C_m+W_l} + &- + \frac{1-\eu^{-(C_m+C_n)t}}{C_m+C_n} + \\&- + \frac{1-\eu^{-(L_k+W_l)t}}{L_k+W_l} + + + \frac{1-\eu^{-(L_k+C_n)t}}{L_k+C_n}\biggr] + \end{aligned}\\ + Λ_2(t)&= + \begin{aligned}[t] + ∑_{m,k,n,l}Γ^1_{mk}&Γ^3_{nl}\biggl[\frac{1-\eu^{-(C_m+C_n)t}}{C_m+C_n} + -\frac{1-\eu^{-(L_k+C_n)t}}{L_k+C_n} + \\&-\frac{\eu^{-(C_n+W_l^\ast)t}-\eu^{-(C_m+C_n)t}}{C_m-W_l^\ast} + +\frac{\eu^{-(C_n+W_l^\ast)t}-\eu^{-(L_k+C_n)t}}{L_k-W^\ast_l}\biggr] + \end{aligned} +\end{align} + +Also required for \cref{eq:bathderiv_1} are +\begin{align} + \label{eq:ABconv} + ∫_0^t\dd{s}A(s)\dot{α}_0(t-s) &= ∑_{n,m}\underbrace{\frac{A_nP_m}{L_m-C_n}}_{\equiv + Γ^A_{nm}}\qty[\eu^{-C_n t}-\eu^{-L_m t}]\\ + ∫_0^t\dd{s}B(s)\dot{α}_0(t-s) &= ∑_{n,m}Γ^1_{nm}\qty[\eu^{-C_n t}-\eu^{-L_m t}] +\end{align} +and +\begin{multline} + \label{eq:nonzerotemplim} + ∫_0^t\dd{s}A(s)\qty(α(s)-α_0(s)) =\\ + ∑_{m,n}\frac{A_nG_m}{C_n+W_m}\qty(1-\eu^{-(C_n+W_m)t}) - ∑_{m,n}\frac{A_nU_m}{C_n+Q_m}\qty(1-\eu^{-(C_n+Q_m)t}). +\end{multline}