Pseudomode Formalism in Open Quantum Systems

• Pseudomode formalism is a framework that represents non-Markovian open-system dynamics using a finite set of damped auxiliary modes.
• It reproduces bath correlation functions by matching exponential decay terms, transforming complex memory effects into a tractable Markovian embedding.
• Modern extensions employ coupled Lindblad, quasi-Lindblad, and purified constructions, enabling simulations of nonlinear and fermionic interactions with maintained accuracy.

Pseudomode formalism is a representation of non-Markovian open-system dynamics in which a structured environment is replaced by a finite set of auxiliary damped modes whose correlation functions reproduce the bath influence on the system. In its standard form, a system linearly coupled to a Gaussian bosonic bath is embedded into an enlarged Markovian dynamics for the system plus pseudomodes; tracing out the auxiliary modes then reproduces the original reduced dynamics. Modern versions extend this idea from independent damped modes to coupled Lindblad pseudomodes, quasi-Lindblad embeddings with relaxed complete positivity, purified and stochastic constructions, fermionic impurity models, nonlinear system-bath couplings, and input-output formulations for environmental observables (Huang et al., 12 Jun 2025, Park et al., 2024, Liang et al., 2024, Thoenniss et al., 2024, Zhang et al., 25 Feb 2025).

1. Foundational concept

The canonical setting is a system $S$ coupled linearly to a Gaussian environment $B$ with continuous degrees of freedom. A representative bosonic Hamiltonian is

$$\hat H=\hat H_S+\hat H_B+\hat H_{SB},\qquad \hat H_B=\int_0^\infty \omega\,\hat b_\omega^\dagger \hat b_\omega\,d\omega,\qquad \hat H_{SB}=\hat S\,\hat B,$$

with

$$\hat B=\int_0^\infty \sqrt{J(\omega)}\big(\hat b_\omega+\hat b_\omega^\dagger\big)\,d\omega.$$

The reduced system dynamics is determined by the bath correlation function

$$C(t)=\mathrm{Tr}\!\left(\hat B(t)\hat B(0)\,\hat\rho_B(0)\right).$$

Pseudomode theory replaces the continuous bath by finitely many auxiliary modes whose bath correlation function is designed to match $C(t)$ over the time interval of interest. If the pseudomode correlation matches $C(t)$ on $[0,T]$, then the reduced system dynamics is reproduced on $[0,T]$ for all bounded observables (Huang et al., 12 Jun 2025).

Historically, the formalism is especially transparent in single-excitation models. For a two-level system coupled to a structured reservoir within the rotating-wave approximation, the excited-state amplitude obeys

$$\frac{d}{dt}a_0(t)=-\int_0^t f(t-t')\,a_0(t')\,dt',$$

where the memory kernel is fixed by the spectral structure. When the reservoir structure is approximated by a sum of Lorentzians, each Lorentzian component yields one pseudomode, and the nonlocal equation is mapped to Markovian dynamics of the enlarged system consisting of the two-level system and the pseudomodes (Ohyama et al., 2017).

A central physical interpretation follows from this mapping: the pseudomode acts as the “memory part” of the reservoir. Information and excitation that have left the system can be stored in the pseudomode and later flow back, producing non-exponential decay, oscillations, revivals, and other signatures of non-Markovianity. In the damped Jaynes–Cummings case, the expectation value of the pseudomode occupation time provides a direct operational interpretation of reservoir memory time, $B$0 (Ohyama et al., 2017).

2. Correlation functions, poles, and exponential structure

The formalism is driven by correlation-function matching rather than by microscopic bath discretization. When the relevant bath response is meromorphic, the bath correlation function can be written exactly or approximately as a finite sum of damped oscillations. For a spectral function with poles $B$1 in the lower half-plane,

$B$2

so each pole corresponds to one damped auxiliary mode (Das et al., 24 Sep 2025).

This pole-based picture is the modern generalization of the Lorentzian construction. For a single Lorentzian spectral density,

$B$3

the correlation function is proportional to $B$4, which is exactly the correlation of one damped harmonic pseudomode of frequency $B$5 and damping rate $B$6. A sum of Lorentzians yields a sum of such exponential terms and therefore a multi-pseudomode embedding (Das et al., 24 Sep 2025).

More generally, pseudomode theory can be formulated as a realization problem for transfer functions. If $B$7, the Laplace transform of $B$8, admits a rational approximation

$B$9

then

$$\hat H=\hat H_S+\hat H_B+\hat H_{SB},\qquad \hat H_B=\int_0^\infty \omega\,\hat b_\omega^\dagger \hat b_\omega\,d\omega,\qquad \hat H_{SB}=\hat S\,\hat B,$$0

and the task is to realize this finite-dimensional representation as a physically meaningful open-system model (Huang et al., 12 Jun 2025).

A closely related viewpoint appears in many-body operator dynamics. In the Heisenberg recursion method, the autocorrelation of an observable is represented through the resolvent

$$\hat H=\hat H_S+\hat H_B+\hat H_{SB},\qquad \hat H_B=\int_0^\infty \omega\,\hat b_\omega^\dagger \hat b_\omega\,d\omega,\qquad \hat H_{SB}=\hat S\,\hat B,$$1

with $$\hat H=\hat H_S+\hat H_B+\hat H_{SB},\qquad \hat H_B=\int_0^\infty \omega\,\hat b_\omega^\dagger \hat b_\omega\,d\omega,\qquad \hat H_{SB}=\hat S\,\hat B,$$2. After tridiagonalization in Krylov space and the introduction of artificial dissipation, the poles of $$\hat H=\hat H_S+\hat H_B+\hat H_{SB},\qquad \hat H_B=\int_0^\infty \omega\,\hat b_\omega^\dagger \hat b_\omega\,d\omega,\qquad \hat H_{SB}=\hat S\,\hat B,$$3 yield a pseudomode expansion

$$\hat H=\hat H_S+\hat H_B+\hat H_{SB},\qquad \hat H_B=\int_0^\infty \omega\,\hat b_\omega^\dagger \hat b_\omega\,d\omega,\qquad \hat H_{SB}=\hat S\,\hat B,$$4

where the complex frequencies encode both oscillation and decay (Teretenkov et al., 2024).

3. Markovian embeddings and generator classes

The simplest realization is a Lindblad master equation for the system and independent damped pseudomodes. With bosonic pseudomodes $$\hat H=\hat H_S+\hat H_B+\hat H_{SB},\qquad \hat H_B=\int_0^\infty \omega\,\hat b_\omega^\dagger \hat b_\omega\,d\omega,\qquad \hat H_{SB}=\hat S\,\hat B,$$5,

$$\hat H=\hat H_S+\hat H_B+\hat H_{SB},\qquad \hat H_B=\int_0^\infty \omega\,\hat b_\omega^\dagger \hat b_\omega\,d\omega,\qquad \hat H_{SB}=\hat S\,\hat B,$$6

where

$$\hat H=\hat H_S+\hat H_B+\hat H_{SB},\qquad \hat H_B=\int_0^\infty \omega\,\hat b_\omega^\dagger \hat b_\omega\,d\omega,\qquad \hat H_{SB}=\hat S\,\hat B,$$7

and $$\hat H=\hat H_S+\hat H_B+\hat H_{SB},\qquad \hat H_B=\int_0^\infty \omega\,\hat b_\omega^\dagger \hat b_\omega\,d\omega,\qquad \hat H_{SB}=\hat S\,\hat B,$$8. In the stochastic pseudomode model, the same enlarged Lindbladian is supplemented by a classical stochastic drive $$\hat H=\hat H_S+\hat H_B+\hat H_{SB},\qquad \hat H_B=\int_0^\infty \omega\,\hat b_\omega^\dagger \hat b_\omega\,d\omega,\qquad \hat H_{SB}=\hat S\,\hat B,$$9 that reproduces the time-symmetric part of the bath correlation (Cirio et al., 2023, Luo et al., 2023).

Coupled Lindblad pseudomodes generalize this by allowing coherent and dissipative couplings among the auxiliary modes. The augmented system contains bosonic modes $$\hat B=\int_0^\infty \sqrt{J(\omega)}\big(\hat b_\omega+\hat b_\omega^\dagger\big)\,d\omega.$$0 with

$$\hat B=\int_0^\infty \sqrt{J(\omega)}\big(\hat b_\omega+\hat b_\omega^\dagger\big)\,d\omega.$$1

and a matrix-valued dissipator

$$\hat B=\int_0^\infty \sqrt{J(\omega)}\big(\hat b_\omega+\hat b_\omega^\dagger\big)\,d\omega.$$2

The resulting bath correlation function is

$$\hat B=\int_0^\infty \sqrt{J(\omega)}\big(\hat b_\omega+\hat b_\omega^\dagger\big)\,d\omega.$$3

The conditions $$\hat B=\int_0^\infty \sqrt{J(\omega)}\big(\hat b_\omega+\hat b_\omega^\dagger\big)\,d\omega.$$4 and $$\hat B=\int_0^\infty \sqrt{J(\omega)}\big(\hat b_\omega+\hat b_\omega^\dagger\big)\,d\omega.$$5 ensure completely positive and trace-preserving dynamics, so the enlarged evolution is a quantum channel (Huang et al., 12 Jun 2025).

Quasi-Lindblad pseudomode theory relaxes complete positivity in the auxiliary sector. In its bosonic form,

$$\hat B=\int_0^\infty \sqrt{J(\omega)}\big(\hat b_\omega+\hat b_\omega^\dagger\big)\,d\omega.$$6

which permits complex weights in the exponential representation and therefore faster convergence than positive-weight Lorentzian fits. Exact reduced dynamics is retained if the bath correlation function is matched, but stability depends on gauge choice and on the spectrum of the full auxiliary generator (Park et al., 2024).

Purified pseudomodes go further by splitting correlations into positive- and negative-time components and assigning them to auxiliary modes acting only on one side of the density operator in Liouville space. This yields a pure auxiliary state, direct access to environmental observables, and a natural treatment of non-Gaussian initial bath states and, more recently, nonlinear system-bath couplings (Liang et al., 2024, Zhang et al., 25 Feb 2025).

Variant	Defining feature	Characteristic property
Independent Lindblad pseudomodes	Each auxiliary mode is damped independently	Simple CPTP embedding; Lorentzian spectra
Coupled Lindblad pseudomodes	Dense $$\hat B=\int_0^\infty \sqrt{J(\omega)}\big(\hat b_\omega+\hat b_\omega^\dagger\big)\,d\omega.$$7 and positive semidefinite $$\hat B=\int_0^\infty \sqrt{J(\omega)}\big(\hat b_\omega+\hat b_\omega^\dagger\big)\,d\omega.$$8	CPTP by construction; non-Lorentzian effective response
Quasi-Lindblad pseudomodes	System-bath Lindblad-type couplings with relaxed CP	Complex weights and compact fits
Purified pseudomodes	Positive/negative-time Liouville-space split	Pure auxiliary states and environmental observables

4. Construction methods and scaling theory

For Lorentzian decompositions, the construction is direct: each Lorentzian peak defines one pseudomode, and the enlarged dynamics reproduces the target memory kernel exactly in the corresponding model class (Ohyama et al., 2017). More generally, modern pseudomode design is a system-realization problem.

A realization-based construction for coupled Lindblad pseudomodes begins with either the target time-domain correlation $$\hat B=\int_0^\infty \sqrt{J(\omega)}\big(\hat b_\omega+\hat b_\omega^\dagger\big)\,d\omega.$$9 or the frequency-domain response $$C(t)=\mathrm{Tr}\!\left(\hat B(t)\hat B(0)\,\hat\rho_B(0)\right).$$0. The frequency-domain fit is performed through a Loewner-type matrix built from sampled data and an SVD-based realization algorithm. Complete positivity is then enforced by a semidefinite program over a positive matrix $$C(t)=\mathrm{Tr}\!\left(\hat B(t)\hat B(0)\,\hat\rho_B(0)\right).$$1,

$$C(t)=\mathrm{Tr}\!\left(\hat B(t)\hat B(0)\,\hat\rho_B(0)\right).$$2

followed by

$$C(t)=\mathrm{Tr}\!\left(\hat B(t)\hat B(0)\,\hat\rho_B(0)\right).$$3

This avoids the non-convex optimization used in earlier constructions and yields pseudomode parameters $$C(t)=\mathrm{Tr}\!\left(\hat B(t)\hat B(0)\,\hat\rho_B(0)\right).$$4 directly (Huang et al., 12 Jun 2025).

A major theoretical result is that, under analyticity assumptions on the spectral density or correlation function, the number of coupled Lindblad pseudomodes needed to approximate the bath correlation function, and hence the reduced dynamics up to time $$C(t)=\mathrm{Tr}\!\left(\hat B(t)\hat B(0)\,\hat\rho_B(0)\right).$$5 with precision $$C(t)=\mathrm{Tr}\!\left(\hat B(t)\hat B(0)\,\hat\rho_B(0)\right).$$6, scales as

$$C(t)=\mathrm{Tr}\!\left(\hat B(t)\hat B(0)\,\hat\rho_B(0)\right).$$7

The same polylogarithmic structure appears in fermionic impurity models, where hybridization kernels are approximated by sums of complex exponentials and compressed further by interpolative decomposition. There, the analytic construction gives $$C(t)=\mathrm{Tr}\!\left(\hat B(t)\hat B(0)\,\hat\rho_B(0)\right).$$8, while the compressed count behaves empirically as

$$C(t)=\mathrm{Tr}\!\left(\hat B(t)\hat B(0)\,\hat\rho_B(0)\right).$$9

suggesting near-optimality (Huang et al., 12 Jun 2025, Thoenniss et al., 2024).

Alternative construction routes include AAA rational approximation, Prony-type exponential fitting, matrix-pencil and ESPRIT procedures for complex-frequency extraction, and Krylov-space pseudomode expansions for many-body correlation functions. This suggests that pseudomode theory is as much a realization theory for correlation kernels and self-energies as it is a specific master-equation ansatz (Thoenniss et al., 2024, Teretenkov et al., 2024).

5. Extensions, domains of application, and physical interpretation

The formalism now spans a broad range of models. In bosonic open systems, coupled Lindblad pseudomodes reproduce spin-boson population dynamics and absorption spectra with few auxiliary modes, including ultra-strong coupling and structured baths with sharp resonances (Huang et al., 12 Jun 2025). In cavity and nanophotonic settings, a Lorentzian near-field spectral density yields a single pseudomode description of a quantum emitter near a metal-dielectric interface, while a single auxiliary mode plus a constant background captures the Fano effect in dissipative cavity QED (Gonzalez-Tudela et al., 2014, Kobayashi et al., 15 Jan 2026).

In quantum metrology, a pseudomode-based Markovian embedding turns temporally and spatially correlated Gaussian bosonic noise into a finite-dimensional Lindbladian problem. This enables quantum-comb optimization and correlated-noise quantum Fisher information bounds for adaptive sensing protocols, with the pseudomode Hilbert space serving as inaccessible memory between control steps (Das et al., 24 Sep 2025).

Fermionic extensions map lesser and greater hybridization functions to auxiliary damped modes and corresponding Lindblad dynamics, providing compact embeddings for Anderson and related impurity models (Huang et al., 12 Jun 2025, Thoenniss et al., 2024). In strong-coupling thermodynamics, the pseudomode embedding permits heat, work, interaction energy, and entropy production to be written in terms of one-time expectation values of the system and pseudomodes rather than two-time bath observables (Albarelli et al., 2024).

Purified input-output pseudomodes broaden the scope from reduced states to environmental observables. They provide access to output fields, spectra, occupations, multiphoton transfer, and non-Gaussian initial bath states in delayed waveguide and cavity-array settings, while maintaining exact reduced dynamics when the fitted correlation structure is exact (Liang et al., 2024). Their nonlinear extension handles interactions of the form

$C(t)$0

provided the Gaussian bath correlations and equal-time contractions are reproduced; the method has been demonstrated for spontaneous emission in a lossy cavity and resonance fluorescence of a driven quantum dot in a phonon environment (Zhang et al., 25 Feb 2025).

6. Subtleties, misconceptions, and unresolved issues

A persistent misconception is that pseudomodes must correspond to physical oscillators in the ordinary sense. In fact, the formalism permits effective degrees of freedom with complex couplings, complex frequencies, non-Hermitian generators, or even non-positive auxiliary dynamics, provided the relevant correlation functions are reproduced. Non-Hermitian pseudomodes can reduce the number of modes needed, especially for finite-temperature underdamped and Drude baths, and exact reduced dynamics follows from matching the advanced and retarded auxiliary correlators rather than from positivity of the pseudomode state itself (Menczel et al., 2024).

A second subtlety concerns complete positivity. Coupled Lindblad pseudomodes are CP and TP by construction and therefore suitable for numerically stable simulation and quantum-channel implementation (Huang et al., 12 Jun 2025). Quasi-Lindblad and non-Hermitian pseudomodes may still reproduce the exact reduced system dynamics, but the global auxiliary evolution can become unstable, and gauge choice matters. This suggests a distinction between exactness of reduced dynamics and physical realizability of the embedding (Park et al., 2024, Menczel et al., 2024).

A third issue is basis efficiency. Independent damped modes generate sums of Lorentzians, but Lorentzians are not a complete basis for arbitrary structured spectral densities. In particular, the effective spectral density of many uncoupled pseudomodes does not necessarily converge in the limit of infinitely many modes. Couplings between pseudomodes, and even non-diagonalizable non-Hermitian single-particle Hamiltonians, generate higher-order pole structures such as derivative-of-Lorentzian and squared-Lorentzian terms that cannot be obtained from diagonalizable uncoupled models (Alford et al., 19 Sep 2025).

Finally, pseudomode formalism is no longer restricted to linear bosonic memory kernels. Current developments emphasize representability: if the eliminated sector produces a rational retarded self-energy, then a finite set of damped auxiliary modes can reproduce it, irrespective of nonlinearities in the retained subsystem. This suggests that the modern formalism is best understood as a general finite-dimensional realization theory for bath correlation functions, transfer functions, and self-energies, with Lindblad, quasi-Lindblad, purified, stochastic, fermionic, and nonlinear variants occupying different points on the spectrum between physical realizability, compactness, and numerical convenience (G. et al., 5 May 2026, Luo et al., 2023).