Variational Approach to Non‐Markovian Dynamics and Two‐Time Correlations Using Pseudomodes

Markovian pseudomode models are shown to accurately describe non‐Markovian and non‐perturbative reduced dynamics and multi‐time correlation functions under certain constraints. However, it is challenging to solve the Lindblad master equation for such models when a number of pseudomodes are involved. A variational approach based on multiple Davydov states is used to solve the master equation of pseudomode models for the multi‐mode Jaynes–Cummings (JC) model, and the spin‐boson system across four distinct configurations, enabling computation of the two‐time correlation function that incorporates bath phonon occupation dynamics. In parallel, calculations of the unitary evolution of full Hamiltonians are chosen as benchmark references. Benchmark comparisons demonstrate that the variational approach for Markovian pseudomode models successfully predicts both the single‐time and two‐time expectations within the multi‐mode JC model. In contrast, in the spin‐boson model, the variational simulations for the pseudomode system and the full Hamiltonian yield inconsistent bath occupation numbers in a non‐perturbative strong coupling regime despite that these two kinds of simulations predict almost identical reduced dynamics, resulting from a fitting artefact in the pseudomode model. Our findings offer insights into the applicability of the non‐unitary variational approach on the simulation of two‐time correlation functions within the pseudomode framework.