A method for the numerical solution of the dissipative Liouville–von Neumann equation is presented. The reduced density operator is expanded in a basis of time-dependent wave functions. This guarantees that the size of the basis required in the representation of the density operator is minimal. Equations of motion for the expansion wave functions and the density matrix elements are obtained from the Dirac–Frenkel variational principle. The numerical effort of the method scales proportional to N, where N is the dimension of the relevant Hilbert space. As a first example, the dynamics of a three-mode system with vibronic coupling (with N≈50 000) coupled to a dissipative bath has been studied, modeling the S1/S2 states of pyrazine. For the cw-absorption spectrum, fast convergence with respect to the number of expansion wave functions has been obtained.
A method for the numerical solution of the dissipative Liouville–von Neumann equation is presented. The reduced density operator is expanded in a basis of time-dependent wave functions. This guarantees that the size of the basis required in the representation of the density operator is minimal. Equations of motion for the expansion wave functions and the density matrix elements are obtained from the Dirac–Frenkel variational principle. The numerical effort of the method scales proportional to N, where N is the dimension of the relevant Hilbert space. As a first example, the dynamics of a three-mode system with vibronic coupling (with N≈50 000) coupled to a dissipative bath has been studied, modeling the S1/S2 states of pyrazine. For the cw-absorption spectrum, fast convergence with respect to the number of expansion wave functions has been obtained.
A mathematical and numerical framework has been worked out to represent the density operator in phase space and to propagate it in time under dissipative conditions. The representation of the density operator is based on the Fourier pseudospectral method which allows a description both in configuration as well as in momentum space. A new propagation scheme which treats the complex eigenvalue structure of the dissipative Liouville superoperator has been developed. The framework has been designed to incorporate modern computer architecture such as parallelism and vectorization. Comparing the results to closed-form solutions exponentially fast convergence characteristics in phase space as well as in the time propagation is demonstrated. As an example of its usefulness, the new method has been successfully applied to dissipation under the constraint of selection rules. More specifically, a harmonic oscillator which relaxes to equilibrium under the constraint of second-order coupling to the bath was studied.
Basing on our hierarchical equations of motion for time-dependent quantum transport [X. Zheng, G. H. Chen, Y. Mo, S. K. Koo, H. Tian, C. Y. Yam, and Y. J. Yan, J. Chem. Phys. 133, 114101 (2010)], we develop an efficient and accurate numerical algorithm to solve the Liouville-von-Neumann equation. We solve the real-time evolution of the reduced single-electron density matrix at the tight-binding level. Calculations are carried out to simulate the transient current through a linear chain of atoms, with each represented by a single orbital. The self-energy matrix is expanded in terms of multiple Lorentzian functions, and the Fermi distribution function is evaluated via the Padè spectrum decomposition. This Lorentzian-Padè decomposition scheme is employed to simulate the transient current. With sufficient Lorentzian functions used to fit the self-energy matrices, we show that the lead spectral function and the dynamics response can be treated accurately. Compared to the conventional master equation approaches, our method is much more efficient as the computational time scales cubically with the system size and linearly with the simulation time. As a result, the simulations of the transient currents through systems containing up to one hundred of atoms have been carried out. As density functional theory is also an effective one-particle theory, the Lorentzian-Padè decomposition scheme developed here can be generalized for first-principles simulation of realistic systems.
We introduce a non-Hermitian approximation of Bloch optical equations. This approximation provides a complete description of the excitation, relaxation, and decoherence dynamics of ensembles of coupled quantum systems in weak laser fields, taking into account collective effects and dephasing. In the proposed method, one propagates the wave function of the system instead of a complete density matrix. Relaxation and dephasing are taken into account via automatically adjusted time-dependent gain and decay rates. As an application, we compute the numerical wave packet solution of a time-dependent non-Hermitian Schrödinger equation describing the interaction of electromagnetic radiation with a quantum nano-structure, and compare the calculated transmission, reflection, and absorption spectra with those obtained from the numerical solution of the Liouville-von Neumann equation. It is shown that the proposed wave packet scheme is significantly faster than the propagation of the full density matrix while maintaining small error. We provide the key ingredients for easy-to-use implementation of the proposed scheme and identify the limits and error scaling of this approximation.
Methods of quantum nuclear wave-function dynamics have become very efficient in simulating large isolated systems using the time-dependent variational principle (TDVP). However, a straightforward extension of the TDVP to the density matrix framework gives rise to methods that do not conserve the energy in the isolated system limit and the total system population for open systems where only energy exchange with environment is allowed. These problems arise when the system density is in a mixed state and is simulated using an incomplete basis. Thus, the basis set incompleteness, which is inevitable in practical calculations, creates artificial channels for energy and population dissipation. To overcome this unphysical behavior, we have introduced a constrained Lagrangian formulation of TDVP applied to a non-stochastic open system Schrödinger equation [L. Joubert-Doriol, I. G. Ryabinkin, and A. F. Izmaylov, J. Chem. Phys. 141, 234112 (2014)]. While our formulation can be applied to any variational ansatz for the system density matrix, derivation of working equations and numerical assessment is done within the variational multiconfiguration Gaussian approach for a two-dimensional linear vibronic coupling model system interacting with a harmonic bath.
The feasibility of using a basis set approach to the study of quantum dissipative dynamics is investigated for the spin-boson model, a system of two discrete states linearly coupled to a harmonic bath. The infinite Hamiltonian is discretized to a finite number of degrees of freedom. Traditional basis set approach, in a multiconfiguration time-dependent Hartree context, is used to solve the time-dependent Schrödinger equations by explicitly including all the degrees of freedom (“system”+“bath”). Quantities such as the reduced density matrix are then evaluated via a quadrature summation/Monte Carlo procedure over a certain number of time-dependent wave functions. Numerically exact results are obtained by systematically increasing the number of bath modes used to represent the condensed phase environment, as well as other variational parameters (number of basis functions, configurations, etc.). The potential of the current method is briefly discussed.
We investigate the effect of inter-mode coupling on the vibrational relaxation dynamics of molecules in weak dissipative environments. The simulations are performed within the reduced density matrix formalism in the Markovian regime, assuming a Lindblad form for the system-bath interaction. The prototypical two-dimensional model system representing two CO molecules approaching a Cu(100) surface is adapted from an ab initio potential, while the diatom-diatom vibrational coupling strength is systematically varied. In the weak system-bath coupling limit and at low temperatures, only first order non-adiabatic uni-modal coupling terms contribute to surface-mediated vibrational relaxation. Since dissipative dynamics is non-unitary, the choice of representation will affect the evolution of the reduced density matrix. Two alternative representations for computing the relaxation rates and the associated operators are thus compared: the fully coupled spectral basis, and a factorizable ansatz. The former is well-established and serves as a benchmark for the solution of Liouville-von Neumann equation. In the latter, a contracted grid basis of potential-optimized discrete variable representation is tailored to incorporate most of the inter-mode coupling, while the Lindblad operators are represented as tensor products of one-dimensional operators, for consistency. This procedure results in a marked reduction of the grid size and in a much more advantageous scaling of the computational cost with respect to the increase of the dimensionality of the system. The factorizable method is found to provide an accurate description of the dissipative quantum dynamics of the model system, specifically of the time evolution of the state populations and of the probability density distribution of the molecular wave packet. The influence of intra-molecular vibrational energy redistribution appears to be properly taken into account by the new model on the whole range of coupling strengths. It demontrates that most of the mode mixing during relaxation is due to the potential part of the Hamiltonian and not to the coupling among relaxation operators.