Fermionic Gaussian Operations

• Fermionic Gaussian operations are quadratic transformations on fermionic modes that preserve covariance matrix structure, enabling tractable simulation and analysis.
• They are implemented via linear Majorana evolution using matchgates, displaced unitaries, and Gaussian channels, providing efficient operator and state updates.
• These operations bridge free-fermion models with topological classifications, guiding efficient quantum simulations and innovative phase-space methods.

Searching arXiv for recent and foundational papers on fermionic Gaussian operations and closely related topics. Fermionic Gaussian operations are the free-fermionic transformations associated with quadratic structure in creation, annihilation, or Majorana operators. In the standard parity-preserving setting, they are generated by quadratic Hamiltonians and act linearly on Majorana modes; Gaussian states are then completely specified by their two-point correlators or covariance matrices, so their evolution closes on finite-dimensional matrix data (Cudby et al., 2024). The modern literature also distinguishes broader but related objects: Gaussian channels acting affinely on covariance matrices, displaced Gaussian operations with nonzero linear Majorana terms, Gaussian projectors and Gaussian operator kernels used in phase-space methods, and operations that sit just outside strict Gaussianity—most notably the fermionic partial transpose (Swingle et al., 2018, Lyu et al., 2024, Eisler et al., 2015).

1. Formal definition and covariance structure

A standard Majorana convention starts from fermionic mode operators and defines

$$\gamma_{2l-1}=a_l+a_l^\dagger,\qquad \gamma_{2l}=-i(a_l-a_l^\dagger),$$

or equivalently

$$a_{2m-1}=c_m+c_m^\dagger,\qquad a_{2m}=i(c_m-c_m^\dagger),$$

with canonical anticommutation relations $\{\gamma_j,\gamma_k\}=2\delta_{jk}$ and $\{a_k,a_l\}=2\delta_{kl}$ (Cudby et al., 2024, Eisler et al., 2015). In this language, fermionic Gaussian states are thermal states of quadratic Hamiltonians or, more generally, operator exponentials of quadratic forms. Representative forms are

$$\rho=\frac{1}{Z}\exp\!\left(\sum_{k,l}\frac{W_{kl}}{4}a_k a_l\right)$$

and

$\rho_h=\frac{e^{-H}}{\Tr(e^{-H})},$

with quadratic $H$ (Eisler et al., 2015, Cudby et al., 2024).

Their efficient descriptor is the covariance or correlation matrix. Two common conventions are

$$\Gamma_{kl}=\frac{1}{2}\langle [a_k,a_l]\rangle,\qquad \tanh\frac{W}{2}=\Gamma,$$

and

$\Gamma(\rho)_{jk}=\frac{i}{2}\Tr([\gamma_j,\gamma_k]\rho).$

In both cases, Wick’s theorem reconstructs higher correlators from the two-point data, so Gaussian states are completely determined by $\Gamma$ (Eisler et al., 2015, Cudby et al., 2024). For pure states, the exact algebraic constraint depends on convention: one formulation uses $$a_{2m-1}=c_m+c_m^\dagger,\qquad a_{2m}=i(c_m-c_m^\dagger),$$0, while another uses $$a_{2m-1}=c_m+c_m^\dagger,\qquad a_{2m}=i(c_m-c_m^\dagger),$$1 (Dias et al., 2023, Eisler et al., 2015).

This covariance-matrix closure is the basic reason fermionic Gaussian operations are tractable. It implies that state preparation, unitary evolution, reduction, overlap formulas, and many measurement updates can be reduced to orthogonal transformations, determinants, or Pfaffians rather than full Hilbert-space manipulations (Dias et al., 2023, Swingle et al., 2018).

2. Quadratic unitaries, matchgates, and displaced extensions

In the parity-preserving setting, a fermionic Gaussian unitary is generated by a quadratic Hamiltonian,

$$a_{2m-1}=c_m+c_m^\dagger,\qquad a_{2m}=i(c_m-c_m^\dagger),$$2

with $$a_{2m-1}=c_m+c_m^\dagger,\qquad a_{2m}=i(c_m-c_m^\dagger),$$3 real antisymmetric. Its defining property is linear Majorana evolution,

$$a_{2m-1}=c_m+c_m^\dagger,\qquad a_{2m}=i(c_m-c_m^\dagger),$$4

and, for Majorana monomials,

$$a_{2m-1}=c_m+c_m^\dagger,\qquad a_{2m}=i(c_m-c_m^\dagger),$$5

Thus Gaussian unitaries preserve monomial degree and are determined, up to phase, by their induced orthogonal matrix $$a_{2m-1}=c_m+c_m^\dagger,\qquad a_{2m}=i(c_m-c_m^\dagger),$$6 (Cudby et al., 2024).

Operationally, the same class is identified with fermionic linear optics and matchgates. The paper on learning Gaussian operations states that the Gaussian unitaries on $$a_{2m-1}=c_m+c_m^\dagger,\qquad a_{2m}=i(c_m-c_m^\dagger),$$7 modes form the matchgate group $$a_{2m-1}=c_m+c_m^\dagger,\qquad a_{2m}=i(c_m-c_m^\dagger),$$8, and that any Gaussian operation on $$a_{2m-1}=c_m+c_m^\dagger,\qquad a_{2m}=i(c_m-c_m^\dagger),$$9 qubits can be expressed as a circuit of at most $\{\gamma_j,\gamma_k\}=2\delta_{jk}$0 parity-preserving two-qubit matchgates $\{\gamma_j,\gamma_k\}=2\delta_{jk}$1 with $\{\gamma_j,\gamma_k\}=2\delta_{jk}$2 (Cudby et al., 2024). In number-conserving electronic-structure settings, the corresponding Gaussian basis changes are orbital rotations

$\{\gamma_j,\gamma_k\}=2\delta_{jk}$3

which act linearly on creation and annihilation operators and appear as the Gaussian layers in compressed many-body circuit constructions (Rubin et al., 2021).

A broader extension allows nonzero linear terms. Displaced Gaussian unitaries are defined by

$\{\gamma_j,\gamma_k\}=2\delta_{jk}$4

with $\{\gamma_j,\gamma_k\}=2\delta_{jk}$5 and $\{\gamma_j,\gamma_k\}=2\delta_{jk}$6. Their natural state descriptor is the extended covariance matrix

$\{\gamma_j,\gamma_k\}=2\delta_{jk}$7

and they act by orthogonal conjugation on this extended space: $\{\gamma_j,\gamma_k\}=2\delta_{jk}$8 (Lyu et al., 2024). The same work shows that displaced Gaussian circuits are computationally equivalent to nearest-neighbor matchgates augmented by single-qubit gates on the first line, and that every $\{\gamma_j,\gamma_k\}=2\delta_{jk}$9-qubit product state is a displaced Gaussian state (Lyu et al., 2024). This suggests that the standard even Gaussian formalism is not the maximal class of tractable free-fermion operations.

3. Gaussian channels, measurements, and local transformation classes

Fermionic Gaussian operations are not limited to unitaries. In Grassmann form, a Gaussian linear map is specified by a kernel with block data $\{a_k,a_l\}=2\delta_{kl}$0, and a Gaussian channel acts on covariance matrices by the affine rule

$\{a_k,a_l\}=2\delta_{kl}$1

Complete positivity is characterized by

$\{a_k,a_l\}=2\delta_{kl}$2

while trace preservation is equivalent to $\{a_k,a_l\}=2\delta_{kl}$3 and $\{a_k,a_l\}=2\delta_{kl}$4 (Swingle et al., 2018). Within this class, the Petz recovery map for a fermionic Gaussian channel and fermionic Gaussian reference state is itself Gaussian, with explicit channel matrices

$\{a_k,a_l\}=2\delta_{kl}$5

$\{a_k,a_l\}=2\delta_{kl}$6

(Swingle et al., 2018).

Projective measurements can also remain inside Gaussian formalism. For pure Gaussian states, occupation-number measurement has probability

$\{a_k,a_l\}=2\delta_{kl}$7

together with an explicit covariance-matrix update rule (Dias et al., 2023). In a more specialized doubled-system construction, post-selected Bell measurements on interlayer rungs are treated as Gaussian projectors, and the post-selected covariance matrix is updated through the Gaussian composition rule

$\{a_k,a_l\}=2\delta_{kl}$8

(Fang et al., 17 Dec 2025). In that half-filled number-conserving setting, the post-selected state factorizes into a product of Bell pairs on the unmeasured rungs (Fang et al., 17 Dec 2025).

Local Gaussian operations form a much narrower resource theory than generic LOCC. For $\{a_k,a_l\}=2\delta_{kl}$9-mode $$\rho=\frac{1}{Z}\exp\!\left(\sum_{k,l}\frac{W_{kl}}{4}a_k a_l\right)$$0-partite Gaussian fermionic states, Gaussian local unitaries (GLU) give a unique standard form, and there are no non-trivial Gaussian LOCC transformations: any deterministic GLOCC map between pure fully entangled Gaussian fermionic states is already a product of local Gaussian unitaries (Spee et al., 2017). By contrast, Gaussian stochastic LOCC and more general fermionic LOCC are richer and admit nontrivial entanglement-class structure (Spee et al., 2017).

4. The boundary of Gaussianity: partial transpose, negativity, and quasi-Gaussian structure

A central limitation of fermionic Gaussian operations is that some natural maps do not preserve Gaussianity. The most prominent example is the partial transpose. For fermions, even simple examples show that the partial transpose of a fermionic Gaussian state is typically not Gaussian; its correlations do not satisfy Wick factorization in any basis (Eisler et al., 2015). This failure is the basic obstruction to an exact covariance-matrix formula for logarithmic negativity in generic mixed fermionic Gaussian states (Kaliszewski et al., 2016).

The structure is nevertheless highly constrained. For a bipartition $$\rho=\frac{1}{Z}\exp\!\left(\sum_{k,l}\frac{W_{kl}}{4}a_k a_l\right)$$1, with covariance matrix block decomposition

$$\rho=\frac{1}{Z}\exp\!\left(\sum_{k,l}\frac{W_{kl}}{4}a_k a_l\right)$$2

the partial transpose can be written, in a suitable basis, as

$$\rho=\frac{1}{Z}\exp\!\left(\sum_{k,l}\frac{W_{kl}}{4}a_k a_l\right)$$3

where $$\rho=\frac{1}{Z}\exp\!\left(\sum_{k,l}\frac{W_{kl}}{4}a_k a_l\right)$$4 are Gaussian operators determined by

$$\rho=\frac{1}{Z}\exp\!\left(\sum_{k,l}\frac{W_{kl}}{4}a_k a_l\right)$$5

The paper characterizes this as a finite linear combination of exactly two Gaussian operators, and the summary identifies the resulting map as “quasi-Gaussian” rather than Gaussianity-preserving in the strict sense (Eisler et al., 2015).

This residual structure is sufficient for several entanglement calculations. Exact determinant formulas exist for traces of low moments such as $$\rho=\frac{1}{Z}\exp\!\left(\sum_{k,l}\frac{W_{kl}}{4}a_k a_l\right)$$6, and in reflection-symmetric geometries one obtains the lower bound

$$\rho=\frac{1}{Z}\exp\!\left(\sum_{k,l}\frac{W_{kl}}{4}a_k a_l\right)$$7

(Eisler et al., 2015). A complementary approach gives efficiently computable upper and lower negativity bounds via semidefinite programming, the Lagrangian formulation of fermionic linear optics, and products of Gaussian operators (Kaliszewski et al., 2016). The broader significance is precise: full spectral control is generally lost, but polynomially structured quantities remain accessible.

5. Simulation, learning, and matrix-element technology

Because fermionic Gaussian operations are classically simulable, they serve as the free layer for several algorithmic tasks. Unknown Gaussian unitaries can be learned efficiently from black-box access: for an unknown $$\rho=\frac{1}{Z}\exp\!\left(\sum_{k,l}\frac{W_{kl}}{4}a_k a_l\right)$$8, an estimate $$\rho=\frac{1}{Z}\exp\!\left(\sum_{k,l}\frac{W_{kl}}{4}a_k a_l\right)$$9 satisfying

$\rho_h=\frac{e^{-H}}{\Tr(e^{-H})},$0

can be determined with query complexity

$\rho_h=\frac{e^{-H}}{\Tr(e^{-H})},$1

for $\rho_h=\frac{e^{-H}}{\Tr(e^{-H})},$2, and the framework extends recursively to the Matchgate Hierarchy $\rho_h=\frac{e^{-H}}{\Tr(e^{-H})},$3, which contains the Clifford hierarchy with the shift $\rho_h=\frac{e^{-H}}{\Tr(e^{-H})},$4 (Cudby et al., 2024).

The same free-fermion structure underlies classical simulation beyond purely Gaussian inputs. A phase-sensitive Gaussian representation augments the covariance matrix by a reference amplitude,

$\rho_h=\frac{e^{-H}}{\Tr(e^{-H})},$5

allowing exact overlap evaluation in $\rho_h=\frac{e^{-H}}{\Tr(e^{-H})},$6 and simulation of fermionic linear optics on superpositions of Gaussian states. Exact runtimes scale polynomially in system size and Gaussian rank $\rho_h=\frac{e^{-H}}{\Tr(e^{-H})},$7, while approximate runtimes scale linearly in the fermionic Gaussian extent $\rho_h=\frac{e^{-H}}{\Tr(e^{-H})},$8 (Dias et al., 2023). Near-Gaussian circuit simulation pushes this further by decomposing non-Gaussian gates and channels into Gaussian ones. For example,

$\rho_h=\frac{e^{-H}}{\Tr(e^{-H})},$9

with optimal decompositions for diagonal two-qubit fermionic gates and gates acting on Jordan–Wigner-adjacent pairs (Dias et al., 19 Mar 2026). The same work derives channel decompositions showing that stochastic Pauli noise can reduce effective extent, while fermionic magic remains more robust to such noise than stabilizer magic (Dias et al., 19 Mar 2026).

A complementary operator-calculus advance gives exact matrix elements of fermionic Gaussian operators in arbitrary local Pauli product bases as a single Pfaffian. Starting from

$H$0

the matrix element between arbitrary Pauli product states is written as

$H$1

which turns arbitrary-basis amplitudes for Gaussian operators into polynomial-time Pfaffian evaluation (Rajabpour et al., 3 Jun 2025).

6. Operator bases, phase-space formulations, and topological viewpoints

Fermionic Gaussian operations also appear through Gaussian operator bases. Normalized Hermitian Gaussian operators

$H$2

admit a resolution of the identity

$H$3

for any even eigenvalue distribution $H$4 in the nonstandard Altland–Zirnbauer classes (Rosales-Zárate et al., 2012). This establishes Gaussian operators as a continuous positive operator basis for the full fermionic Hilbert space and supports phase-space representations analogous in spirit to bosonic coherent-state methods (Rosales-Zárate et al., 2012).

Phase-space dynamics exploits that basis to represent exact interacting evolution as stochastic motion over Gaussian kernels. In the Gaussian phase-space representation, the state is encoded by a distribution over variables such as normal occupations $H$5, anomalous correlators $H$6, and bosonic amplitudes $H$7, and the Liouville equation becomes a Fokker–Planck equation

$H$8

The deterministic drift alone reproduces pairing mean-field or Hartree–Fock dynamics, whereas the diffusion terms restore exact beyond-mean-field evolution (Ogren et al., 2010, Rousse et al., 2023). A recent development uses numerical Takagi factorizations $H$9 as diffusion gauges and reports that the resulting gauges can double the practical simulation time in benchmark fermionic problems (Rousse et al., 2023).

Finally, locality-preserving unitary Gaussian operations have a topological classification distinct from that of Gaussian states. In an operational framework based directly on covariance matrices and Majorana transformations, the classification of fermionic Gaussian states reproduces the periodic table, but the classification of fermionic Gaussian operations differs, and the relation between the two is mediated by a unitary-to-state homomorphism (Futami et al., 2021). Some topological Gaussian states can be disentangled by symmetry-preserving locality-preserving Gaussian operations, some cannot, and some topological Gaussian operations are “genuinely dynamical,” meaning they are not connected to the identity yet cannot generate any nontrivial topological Gaussian state from a trivial one (Futami et al., 2021). This suggests that the topology of fermionic Gaussian operations is not merely derivative of free-fermion state topology, but an autonomous dynamical structure.