Prerequisites for relevant spectral density and convergence of reduced density matrices at low temperatures

Hierarchical equations of motion approach with the Drude-Lorentz spectral density has been widely employed in investigating quantum dissipative phenomena. However, it is often computationally costly for low-temperature systems because a number of Matsubara frequencies are involved. In this note, we examine a prerequisite required for spectral density, and demonstrate that relevant spectral density may significantly reduce the number of Matsubara terms to obtain convergent results for low temperatures.


Akihito Ishizaki
Institute for Molecular Science, National Institutes of Natural Sciences and School of Physical Sciences, Graduate University for Advanced Studies, Okazaki, 444-8585,

Japan
Hierarchical equations of motion approach with the Drude-Lorentz spectral density has been widely employed in investigating quantum dissipative phenomena. However, it is often computationally costly for low-temperature systems because a number of Matsubara frequencies are involved. In this note, we examine a prerequisite required for spectral density, and demonstrate that relevant spectral density may significantly reduce the number of Matsubara terms to obtain convergent results for low temperatures.
In condensed phase molecular systems, chemical dynamics are significantly impacted by the dynamics of the surrounding environment. 1) A simple way to account for the environmental dynamics is to assume an exponential decay form for the correlation function of bath degrees of freedom, exp(−γt). 2) This corresponds to the Drude-Lorentz model and the spectral density is written as J DL (ω) = 2Λγω/(ω 2 + γ 2 ), where Λ quantifies the system-bath coupling strength. In quantum mechanical treatment, the correlation function is replaced with the symmetrized correlation function of the collective bath co- where ν k = 2πk/β is the bosonic Matsubara frequency andẼ β (ω) = Also, the response function, Φ(t) = (i/ ) [X(t),X(0)] = (2/π) ∞ 0 dω J(ω) sin ωt, is given by Note that transient behaviors of eqs. (1) and (2) are coarse-grained. 11) Particularly, the response function in eq. (2) gives a non-zero value at t = 0 although it should vanish by definition, Φ(0) = (i/ ) [X,X] = 0. Other problems arise owing to using the spectral density and the symmetrized correlation function. As shown in Fig. 1, the Drude-Lorentz spectral density exhibits a long tail in the high frequency region correspondingly to the coarse-grained nature. Consequently, a number of the Matsubara frequencies are required for the convergence of eq. (1). To compute quantum dynamics influenced by environments, the hierarchical equations of motion approach has been widely employed. [3][4][5][6] However, the long-tail and convergence problems make the hierarchical equations approach for low temperature computationally costly. 5) Hence, crude truncation was employed, 7,8) and a more sophisticated but applicable approach with Padé approximant was proposed. 9,10) To address the the long-tail and convergence problems, we start with the relaxation function defined with the canonical correlation function, 11) Ψ(t) = β X (t);X(0) . The response function and spectral density are derived as cos ωt, respectively. 11) Hence, a prerequisite required for the relaxation function is As an example to satisfy this, we examine the relaxation function of the form, In the limit of Γt ≪ 1, the function approaches the Gaussian form, providing the relevant short-time approximation. 11) The characteristic time constant is obtained as τ c = ∞ 0 dt Ψ(t)/Ψ(0) = 2/Γ, and thus, Γ/2 corresponds to γ in the Drude-Lorentz model. The response function and spectral density are obtained as Φ(t) = 2ΛΓ 2 te −Γt and J(ω) = 4ΛΓ 3 ω/(ω 2 + Γ 2 ) 2 , respectively. In the limit of ω ≪ Γ, the spectral density is written as the Ohmic form, J(ω) = 4ΛΓ −1 ω. As shown in Fig. 1, this spectral density does not exhibit a long tail, indicating that the tail of the Drude-Lorentz spectral density is governed by the condition in eq. (3). Therefore, it is expected that the symmetrized correlation function associated to eq. (4) converges with a smaller number of the Matsubara frequencies. form is a Brownian oscillator, 12,13) and the response function and spectral density are , and D (k) (t) is the k-th quantum correction term defined by withẼ (k) β (ω) = β −1 2ω 2 /(ω 2 − ν 2 k ). For numerical demonstration of the hierarchical equations for eq. (4), we con-temperature classical correlation 1st-order corr. 2nd-order corr.  The longitudinal relaxation can be modeled by a simple spin-boson Hamiltonian, The relaxation is described as the one-phonon process and hence the energy gap |E 1 − E 2 | needs to resonate with spectral density J(ω). Table I presents the density matrix element 1|ρ|1 at the equilibrium state at several temperatures for the fixed energy gap, The bath parameters are set to Λ = 50 cm −1 , Γ = 100 cm −1 and ǫ = 25 cm −1 , which correspond to the spectral density in Fig. 1. Hence, the energy gap resonates with the spectral density. At lower temperatures, the Matsubara frequency ν k is smaller, and thus the contribution of iJ(iν k ) in eq. The outer-sphere electron transfer reaction is described with the Hamiltonian, 14) H where E • 1/2 is the equilibrium state energy of the donor/accepter state. The system-bath coupling strength Λ is the reorganization energy for the electron transfer, and the collective bath coordinateX corresponds to the solvation coordinate. The electron transfer reaction is a thermally activated barrier crossing process; thus, resonance between the energy gap |E • 1 − E • 2 | and the spectral density J(ω) is unnecessary. Table II presents the density matrix element 1|ρ|1 at the equilibrium state at temperature 300 K for varying values of the energy gap, The inter-site coupling and the reorganization energy are J 12 = 50 cm −1 and Λ = 500 cm −1 , respectively. 14) The bath characteristic time (timescale of the solvation dynamics) is set to τ c = 100 fs with Γ = 100 cm −1 and ǫ = 25 cm −1 , correspondingly to the spectral density in Fig. 1. As presented in Table II, practically convergent equilibrium populations are obtained at low temperatures of β (E • 1 −E • 2 ) > 1 even in the absence of the low-temperature correction terms. This is explained as follows: At the fixed temperature 300 K, the Matsubara frequency ν k takes the value of 1310 cm −1 × k. Hence, iJ(iν k ) in eq. (5) is insignificant when the spectral density J(ω) presented in Fig. 1 is considered.
In summary, we discussed the the prerequisite required for relevant spectral density 5/7 or equivalently relaxation functions, eq. (3). The spectral density that satisfies the prerequisite does not exhibit a long tail in the high frequency region. Consequently, a smaller number of the Matsubara frequencies are required for the convergence in the symmetrized correlation function, potentially reducing the computational cost to obtain convergent results in the hierarchical equations approach for low temperature systems.

Acknowledgment
The author thanks Akihito Kato, Yuta Fujihashi, and Tatsushi