arXiv:2602.05058

Learning fermionic linear optics with Heisenberg scaling and physical operations

Aria Christensen, Andrew Zhao

quant-phcs.DScs.LGSubmitted Feb 4, 2026• Announced Feb 6, 2026• 0 citations

Abstract

We revisit the problem of learning fermionic linear optics (FLO), also known as fermionic Gaussian unitaries. Given black-box query access to an unknown FLO, previous proposals required $\widetilde{\mathcal{O}}(n^5 / \varepsilon^2)$ queries, where $n$ is the system size and $\varepsilon$ is the error in diamond distance. These algorithms also use unphysical operations (i.e., violating fermionic superselection rules) and/or $n$ auxiliary modes to prepare Choi states of the FLO. In this work, we establish efficient and experimentally friendly protocols that obey superselection, use minimal ancilla (at most $1$ extra mode), and exhibit improved dependence on both parameters $n$ and $\varepsilon$. For arbitrary (active) FLOs this algorithm makes at most $\widetilde{\mathcal{O}}(n^4 / \varepsilon)$ queries, while for number-conserving (passive) FLOs we show that $\mathcal{O}(n^3 / \varepsilon)$ queries suffice. The complexity of the active case can be further reduced to $\widetilde{\mathcal{O}}(n^3 / \varepsilon)$ at the cost of using $n$ ancilla. This marks the first FLO learning algorithm that attains Heisenberg scaling in precision. As a side result, we also demonstrate an improved copy complexity of $\widetilde{\mathcal{O}}(n η^2 / \varepsilon^2)$ for time-efficient state tomography of $η$-particle Slater determinants in $\varepsilon$ trace distance, which may be of independent interest.

Comments: 56 pages

Related papers

Learning interacting fermionic Hamiltonians at the Heisenberg limit

Arjun Mirani, Patrick M. Hayden

arXiv:2403.00069Physical Review ADec 23, 20244 citations

Query-optimal unitary channel tomography in diamond distance with parallel access

Entong He, Zihao Li, Yuxiang Yang

quant-pharXiv:2609.38145Sep 30, 20260 citations

Efficient learning of bosonic Gaussian unitaries

Marco Fanizza, Iyer, Vishnu, Junseo Lee +2 more

arXiv:2510.05531Oct 7, 20250 citations

Particle-preserving fermionic shadows with mode-independent sample complexity

Maxwell West, M. Cerezo, Martin Larocca

quant-pharXiv:2606.27254Jun 26, 20260 citations

Mildly-Interacting Fermionic Unitaries are Efficiently Learnable

Vishnu Iyer

arXiv:2504.11318Apr 15, 20250 citations

Learning the closest Slater determinant

Nisarga Paul, Haimeng Zhao, David D. Dai

quant-pharXiv:2607.20623Jul 24, 20260 citations

arXiv:2602.05058

Learning fermionic linear optics with Heisenberg scaling and physical operations

Aria Christensen, Andrew Zhao

quant-phcs.DScs.LGSubmitted Feb 4, 2026• Announced Feb 6, 2026• 0 citations

arXiv HTML PDF

Abstract

We revisit the problem of learning fermionic linear optics (FLO), also known as fermionic Gaussian unitaries. Given black-box query access to an unknown FLO, previous proposals required

\widetilde{\mathcal{O}}(n^5 / \varepsilon^2)

queries, where

n

is the system size and

\varepsilon

is the error in diamond distance. These algorithms also use unphysical operations (i.e., violating fermionic superselection rules) and/or

n

auxiliary modes to prepare Choi states of the FLO. In this work, we establish efficient and experimentally friendly protocols that obey superselection, use minimal ancilla (at most

1

extra mode), and exhibit improved dependence on both parameters

n

and

\varepsilon

. For arbitrary (active) FLOs this algorithm makes at most

\widetilde{\mathcal{O}}(n^4 / \varepsilon)

queries, while for number-conserving (passive) FLOs we show that

\mathcal{O}(n^3 / \varepsilon)

queries suffice. The complexity of the active case can be further reduced to

\widetilde{\mathcal{O}}(n^3 / \varepsilon)

at the cost of using

n

ancilla. This marks the first FLO learning algorithm that attains Heisenberg scaling in precision. As a side result, we also demonstrate an improved copy complexity of

\widetilde{\mathcal{O}}(n η^2 / \varepsilon^2)

for time-efficient state tomography of

η

-particle Slater determinants in

\varepsilon

trace distance, which may be of independent interest.

Comments: 56 pages

Related papers

Learning interacting fermionic Hamiltonians at the Heisenberg limit

Arjun Mirani, Patrick M. Hayden

Efficiently learning an unknown Hamiltonian given access to its dynamics is a problem of interest for quantum metrology, many-body physics, and machine learning. A fundamental question is whether learning can be performed at the Heisenberg limit, where the Hamiltonian evolution time scales inversely with the error $\ensuremath{\varepsilon}$ in the reconstructed parameters. The Heisenberg limit has previously been shown to be achievable for certain classes of qubit and bosonic Hamiltonians. Most recently, a Heisenberg-limited learning algorithm was proposed for a simplified class of fermionic Hubbard Hamiltonians restricted to real hopping amplitudes and zero chemical potential at all sites, along with on-site interactions. In this work, we provide an algorithm to learn a more general class of fermionic Hubbard Hamiltonians at the Heisenberg limit, allowing complex hopping amplitudes and nonzero chemical potentials in addition to the on-site interactions, thereby including several models of physical interest. The required evolution time across all experiments in our protocol is $O(1/\ensuremath{\varepsilon})$ and the number of experiments required to learn all the Hamiltonian parameters is $O[\text{polylog}(1/\ensuremath{\varepsilon})]$, independent of system size as long as each fermionic mode interacts with $O(1)$ other modes. Unlike prior algorithms for bosonic and fermionic Hamiltonians, to obey fermionic parity superselection constraints in our more general setting, our protocol utilizes $O(N)$ ancillary fermionic modes, where $N$ is the system size. Each experiment involves preparing fermionic Gaussian states, interleaving time evolution with fermionic linear optics unitaries, and performing local occupation number measurements on the fermionic modes. The protocol is robust to a constant amount of state preparation and measurement error.Efficiently learning an unknown Hamiltonian given access to its dynamics is a problem of interest for quantum metrology, many-body physics, and machine learning. A fundamental question is whether learning can be performed at the Heisenberg limit, where the Hamiltonian evolution time scales inversely with the error $\ensuremath{\varepsilon}$ in the reconstructed parameters. The Heisenberg limit has previously been shown to be achievable for certain classes of qubit and bosonic Hamiltonians. Most recently, a Heisenberg-limited learning algorithm was proposed for a simplified class of fermionic Hubbard Hamiltonians restricted to real hopping amplitudes and zero chemical potential at all sites, along with on-site interactions. In this work, we provide an algorithm to learn a more general class of fermionic Hubbard Hamiltonians at the Heisenberg limit, allowing complex hopping amplitudes and nonzero chemical potentials in addition to the on-site interactions, thereby including several models of physical interest. The required evolution time across all experiments in our protocol is $O(1/\ensuremath{\varepsilon})$ and the number of experiments required to learn all the Hamiltonian parameters is $O[\text{polylog}(1/\ensuremath{\varepsilon})]$, independent of system size as long as each fermionic mode interacts with $O(1)$ other modes. Unlike prior algorithms for bosonic and fermionic Hamiltonians, to obey fermionic parity superselection constraints in our more general setting, our protocol utilizes $O(N)$ ancillary fermionic modes, where $N$ is the system size. Each experiment involves preparing fermionic Gaussian states, interleaving time evolution with fermionic linear optics unitaries, and performing local occupation number measurements on the fermionic modes. The protocol is robust to a constant amount of state preparation and measurement error.

arXiv:2403.00069Physical Review ADec 23, 20244 citations

Query-optimal unitary channel tomography in diamond distance with parallel access

Entong He, Zihao Li, Yuxiang Yang

In the study of quantum process tomography, it has remained open how to learn a unitary channel in diamond distance with strictly parallel queries as efficiently as with sequential queries. Sequential queries allow one to exploit adaptivity to refine the estimate and adjust the learning strategy accordingly, whereas such adjustments are impossible for parallel queries. From the perspective of higher-order quantum operations, any parallel learning strategy can be simulated by a sequential one, whereas the converse does not hold in general, suggesting that sequential learning strategies are potentially more powerful. Contrary to this intuition, we present a new protocol for learning an unknown $d$ -dimensional unitary channel to within $\varepsilon$ in diamond distance using $\mathcal{O}(d^2/\varepsilon)$ parallel queries, matching the lower bound presented in [Haah, Kothari, O'Donnell, and Tang, FOCS '23]. Our protocol thus closes the gap between parallel and sequential strategies for unitary channel tomography. As applications, it achieves optimal query complexity for boundary-regime quantum channel tomography and improves the query efficiency of the best-known protocol for tomography of fermionic linear optics.In the study of quantum process tomography, it has remained open how to learn a unitary channel in diamond distance with strictly parallel queries as efficiently as with sequential queries. Sequential queries allow one to exploit adaptivity to refine the estimate and adjust the learning strategy accordingly, whereas such adjustments are impossible for parallel queries. From the perspective of higher-order quantum operations, any parallel learning strategy can be simulated by a sequential one, whereas the converse does not hold in general, suggesting that sequential learning strategies are potentially more powerful. Contrary to this intuition, we present a new protocol for learning an unknown $d$-dimensional unitary channel to within $\varepsilon$ in diamond distance using $\mathcal{O}(d^2/\varepsilon)$ parallel queries, matching the lower bound presented in [Haah, Kothari, O'Donnell, and Tang, FOCS '23]. Our protocol thus closes the gap between parallel and sequential strategies for unitary channel tomography. As applications, it achieves optimal query complexity for boundary-regime quantum channel tomography and improves the query efficiency of the best-known protocol for tomography of fermionic linear optics.

quant-pharXiv:2609.38145Sep 30, 20260 citations

Efficient learning of bosonic Gaussian unitaries

Marco Fanizza, Iyer, Vishnu, Junseo Lee +2 more

Bosonic Gaussian unitaries are fundamental building blocks of central continuous-variable quantum technologies such as quantum-optic interferometry and bosonic error-correction schemes. In this work, we present the first time-efficient algorithm for learning bosonic Gaussian unitaries with a rigorous analysis. Our algorithm produces an estimate of the unknown unitary that is accurate to small worst-case error, measured by the physically motivated energy-constrained diamond distance. Its runtime and query complexity scale polynomially with the number of modes, the inverse target accuracy, and natural energy parameters quantifying the allowed input energy and the unitary's output-energy growth. The protocol uses only experimentally friendly photonic resources: coherent and squeezed probes, passive linear optics, and heterodyne/homodyne detection. We then employ an efficient classical post-processing routine that leverages a symplectic regularization step to project matrix estimates onto the symplectic group. In the limit of unbounded input energy, our procedure attains arbitrarily high precision using only $2m+2$ queries, where $m$ is the number of modes. To our knowledge, this is the first provably efficient learning algorithm for a multiparameter family of continuous-variable unitaries.Bosonic Gaussian unitaries are fundamental building blocks of central continuous-variable quantum technologies such as quantum-optic interferometry and bosonic error-correction schemes. In this work, we present the first time-efficient algorithm for learning bosonic Gaussian unitaries with a rigorous analysis. Our algorithm produces an estimate of the unknown unitary that is accurate to small worst-case error, measured by the physically motivated energy-constrained diamond distance. Its runtime and query complexity scale polynomially with the number of modes, the inverse target accuracy, and natural energy parameters quantifying the allowed input energy and the unitary's output-energy growth. The protocol uses only experimentally friendly photonic resources: coherent and squeezed probes, passive linear optics, and heterodyne/homodyne detection. We then employ an efficient classical post-processing routine that leverages a symplectic regularization step to project matrix estimates onto the symplectic group. In the limit of unbounded input energy, our procedure attains arbitrarily high precision using only $2m+2$ queries, where $m$ is the number of modes. To our knowledge, this is the first provably efficient learning algorithm for a multiparameter family of continuous-variable unitaries.

arXiv:2510.05531Oct 7, 20250 citations

Particle-preserving fermionic shadows with mode-independent sample complexity

Maxwell West, M. Cerezo, Martin Larocca

We consider the problem of learning expectation values of particle-preserving operators with respect to an unknown $η$ -particle $n$ -mode fermionic state via classical shadows. Our main application is to estimating overlaps with arbitrary Slater determinant states: While it is known that such overlaps can, in the average case, be learnt to a fixed additive precision with a constant number of samples, the best-known worst case bound is $\mathcal{O}(\sqrt n \log n)$ ; here we improve this to $\mathcal{O}(η\logη)$ , achieving a mode-independent sample cost. Our procedure is also computationally efficient, requiring only classical post-processing which for a generic dense orbital runs in time $\mathcal{O}(nη^2)$ . For the task of estimating the expectation value of a general particle-preserving quadratic fermionic observable $h$ , we prove a sample complexity bound of $\mathcal{O}(η\|h_0\|_2^2)$ , where $h_0$ is the traceless component of $h$ ; the associated classical post-processing scales as $\mathcal{O}(n^2η)$ . Finally, we discuss implementation of the required randomization: in a first-quantized encoding, approximate unitary designs give circuit depths polylogarithmic in the number of modes, contrasting with linear-depth requirements for nearest-neighbor second-quantized matchgate implementations. On the technical side, our proof reduces the extremal shadow variance to harmonic analysis on the AIII symmetric space $U(n)/(U(η)\times U(n-η))$ and evaluates the resulting integral using techniques from the theories of Jacobi ensembles and orthogonal polynomials, in a calculation which may be of independent interest.We consider the problem of learning expectation values of particle-preserving operators with respect to an unknown $η$-particle $n$-mode fermionic state via classical shadows. Our main application is to estimating overlaps with arbitrary Slater determinant states: While it is known that such overlaps can, in the average case, be learnt to a fixed additive precision with a constant number of samples, the best-known worst case bound is $\mathcal{O}(\sqrt n \log n)$; here we improve this to $\mathcal{O}(η\logη)$, achieving a mode-independent sample cost. Our procedure is also computationally efficient, requiring only classical post-processing which for a generic dense orbital runs in time $\mathcal{O}(nη^2)$. For the task of estimating the expectation value of a general particle-preserving quadratic fermionic observable $h$, we prove a sample complexity bound of $\mathcal{O}(η\|h_0\|_2^2)$, where $h_0$ is the traceless component of $h$; the associated classical post-processing scales as $\mathcal{O}(n^2η)$. Finally, we discuss implementation of the required randomization: in a first-quantized encoding, approximate unitary designs give circuit depths polylogarithmic in the number of modes, contrasting with linear-depth requirements for nearest-neighbor second-quantized matchgate implementations. On the technical side, our proof reduces the extremal shadow variance to harmonic analysis on the AIII symmetric space $U(n)/(U(η)\times U(n-η))$ and evaluates the resulting integral using techniques from the theories of Jacobi ensembles and orthogonal polynomials, in a calculation which may be of independent interest.

quant-pharXiv:2606.27254Jun 26, 20260 citations

Mildly-Interacting Fermionic Unitaries are Efficiently Learnable

Vishnu Iyer

Recent work has shown that one can efficiently learn fermionic Gaussian unitaries, also commonly known as nearest-neighbor matchcircuits or non-interacting fermionic unitaries. However, one could ask a similar question about unitaries that are near Gaussian: for example, unitaries prepared with a small number of non-Gaussian circuit elements. These operators find significance in quantum chemistry and many-body physics, yet no algorithm exists to learn them. We give the first such result by devising an algorithm which makes queries to an $n$ -mode fermionic unitary $U$ prepared by at most $O(t)$ non-Gaussian gates and returns a circuit approximating $U$ to diamond distance $\varepsilon$ in time $\textrm{poly}(n,2^t,1/\varepsilon)$ . This resolves a central open question of Mele and Herasymenko under the strongest distance metric. In fact, our algorithm is much more general: we define a property of unitary Gaussianity known as unitary Gaussian dimension and show that our algorithm can learn $n$ -mode unitaries of Gaussian dimension at least $2n - O(t)$ in time $\textrm{poly}(n,2^t,1/\varepsilon)$ . Indeed, this class subsumes unitaries prepared by at most $O(t)$ non-Gaussian gates but also includes several unitaries that require up to $2^{O(t)}$ non-Gaussian gates to construct. In addition, we give a $\textrm{poly}(n,1/\varepsilon)$ -time algorithm to distinguish whether an $n$ -mode unitary is of Gaussian dimension at least $k$ or $\varepsilon$ -far from all such unitaries in Frobenius distance, promised that one is the case. Along the way, we prove structural results about near-Gaussian fermionic unitaries that are likely to be of independent interest.Recent work has shown that one can efficiently learn fermionic Gaussian unitaries, also commonly known as nearest-neighbor matchcircuits or non-interacting fermionic unitaries. However, one could ask a similar question about unitaries that are near Gaussian: for example, unitaries prepared with a small number of non-Gaussian circuit elements. These operators find significance in quantum chemistry and many-body physics, yet no algorithm exists to learn them. We give the first such result by devising an algorithm which makes queries to an $n$-mode fermionic unitary $U$ prepared by at most $O(t)$ non-Gaussian gates and returns a circuit approximating $U$ to diamond distance $\varepsilon$ in time $\textrm{poly}(n,2^t,1/\varepsilon)$. This resolves a central open question of Mele and Herasymenko under the strongest distance metric. In fact, our algorithm is much more general: we define a property of unitary Gaussianity known as unitary Gaussian dimension and show that our algorithm can learn $n$-mode unitaries of Gaussian dimension at least $2n - O(t)$ in time $\textrm{poly}(n,2^t,1/\varepsilon)$. Indeed, this class subsumes unitaries prepared by at most $O(t)$ non-Gaussian gates but also includes several unitaries that require up to $2^{O(t)}$ non-Gaussian gates to construct. In addition, we give a $\textrm{poly}(n,1/\varepsilon)$-time algorithm to distinguish whether an $n$-mode unitary is of Gaussian dimension at least $k$ or $\varepsilon$-far from all such unitaries in Frobenius distance, promised that one is the case. Along the way, we prove structural results about near-Gaussian fermionic unitaries that are likely to be of independent interest.

arXiv:2504.11318Apr 15, 20250 citations

Learning the closest Slater determinant

Nisarga Paul, Haimeng Zhao, David D. Dai

Learning compact, interpretable descriptions of quantum many-body states is an important task in quantum science. We study the task of learning the Slater determinant with maximum fidelity to an arbitrary fermionic many-body state, with motivation from both Hartree-Fock methods and agnostic tomography. Given an $n$ -fermion wavefunction built from $m$ fermionic modes, we provide classical and quantum algorithms returning a Slater determinant with fidelity within $\varepsilon$ of maximal in time $m^{\text{poly}(n,1/\varepsilon)}$ . We prove matching hardness lower bounds, assuming standard complexity conjectures, along some parameter axes. Given access to quantum copies, we prove this can be accomplished with $\text{poly}(m,n,1/\varepsilon)$ copies of $ρ$ . We also show that above a fidelity of $2/3$ any stationary point is the unique global maximum while below $2/3$ the optimization landscape can have spurious stationary points, and hence $2/3$ marks a transition point in the optimization landscape for this problem. We apply the algorithm to the Fermi-Hubbard model, extracting the closest Slater determinant from neural quantum state solutions. Together, our results provide algorithmic tools with provable guarantees in understanding fermionic many-body systems with classical or quantum simulation.Learning compact, interpretable descriptions of quantum many-body states is an important task in quantum science. We study the task of learning the Slater determinant with maximum fidelity to an arbitrary fermionic many-body state, with motivation from both Hartree-Fock methods and agnostic tomography. Given an $n$-fermion wavefunction built from $m$ fermionic modes, we provide classical and quantum algorithms returning a Slater determinant with fidelity within $\varepsilon$ of maximal in time $m^{\text{poly}(n,1/\varepsilon)}$. We prove matching hardness lower bounds, assuming standard complexity conjectures, along some parameter axes. Given access to quantum copies, we prove this can be accomplished with $\text{poly}(m,n,1/\varepsilon)$ copies of $ρ$. We also show that above a fidelity of $2/3$ any stationary point is the unique global maximum while below $2/3$ the optimization landscape can have spurious stationary points, and hence $2/3$ marks a transition point in the optimization landscape for this problem. We apply the algorithm to the Fermi-Hubbard model, extracting the closest Slater determinant from neural quantum state solutions. Together, our results provide algorithmic tools with provable guarantees in understanding fermionic many-body systems with classical or quantum simulation.

quant-pharXiv:2607.20623Jul 24, 20260 citations

We revisit the problem of learning fermionic linear optics (FLO), also known as fermionic Gaussian unitaries. Given black-box query access to an unknown FLO, previous proposals required

\widetilde{\mathcal{O}}(n^5 / \varepsilon^2)

queries, where

n

is the system size and

\varepsilon

is the error in diamond distance. These algorithms also use unphysical operations (i.e., violating fermionic superselection rules) and/or

n

auxiliary modes to prepare Choi states of the FLO. In this work, we establish efficient and experimentally friendly protocols that obey superselection, use minimal ancilla (at most

1

extra mode), and exhibit improved dependence on both parameters

n

and

\varepsilon

. For arbitrary (active) FLOs this algorithm makes at most

\widetilde{\mathcal{O}}(n^4 / \varepsilon)

queries, while for number-conserving (passive) FLOs we show that

\mathcal{O}(n^3 / \varepsilon)

queries suffice. The complexity of the active case can be further reduced to

\widetilde{\mathcal{O}}(n^3 / \varepsilon)

at the cost of using

n

ancilla. This marks the first FLO learning algorithm that attains Heisenberg scaling in precision. As a side result, we also demonstrate an improved copy complexity of

\widetilde{\mathcal{O}}(n η^2 / \varepsilon^2)

for time-efficient state tomography of

η

-particle Slater determinants in

\varepsilon

trace distance, which may be of independent interest.

O(1)

other modes. Unlike prior algorithms for bosonic and fermionic Hamiltonians, to obey fermionic parity superselection constraints in our more general setting, our protocol utilizes

O(N)

ancillary fermionic modes, where

N

is the system size. Each experiment involves preparing fermionic Gaussian states, interleaving time evolution with fermionic linear optics unitaries, and performing local occupation number measurements on the fermionic modes. The protocol is robust to a constant amount of state preparation and measurement error.

d

-dimensional unitary channel to within

\varepsilon

in diamond distance using

\mathcal{O}(d^2/\varepsilon)

parallel queries, matching the lower bound presented in [Haah, Kothari, O'Donnell, and Tang, FOCS '23]. Our protocol thus closes the gap between parallel and sequential strategies for unitary channel tomography. As applications, it achieves optimal query complexity for boundary-regime quantum channel tomography and improves the query efficiency of the best-known protocol for tomography of fermionic linear optics.

2m+2

queries, where

m

is the number of modes. To our knowledge, this is the first provably efficient learning algorithm for a multiparameter family of continuous-variable unitaries.

We consider the problem of learning expectation values of particle-preserving operators with respect to an unknown

η

-particle

n

-mode fermionic state via classical shadows. Our main application is to estimating overlaps with arbitrary Slater determinant states: While it is known that such overlaps can, in the average case, be learnt to a fixed additive precision with a constant number of samples, the best-known worst case bound is

\mathcal{O}(\sqrt n \log n)

; here we improve this to

\mathcal{O}(η\logη)

, achieving a mode-independent sample cost. Our procedure is also computationally efficient, requiring only classical post-processing which for a generic dense orbital runs in time

\mathcal{O}(nη^2)

. For the task of estimating the expectation value of a general particle-preserving quadratic fermionic observable

h

, we prove a sample complexity bound of

\mathcal{O}(η\|h_0\|_2^2)

, where

h_0

is the traceless component of

h

; the associated classical post-processing scales as

\mathcal{O}(n^2η)

. Finally, we discuss implementation of the required randomization: in a first-quantized encoding, approximate unitary designs give circuit depths polylogarithmic in the number of modes, contrasting with linear-depth requirements for nearest-neighbor second-quantized matchgate implementations. On the technical side, our proof reduces the extremal shadow variance to harmonic analysis on the AIII symmetric space

U(n)/(U(η)\times U(n-η))

and evaluates the resulting integral using techniques from the theories of Jacobi ensembles and orthogonal polynomials, in a calculation which may be of independent interest.

n

-mode fermionic unitary

U

prepared by at most

O(t)

non-Gaussian gates and returns a circuit approximating

U

to diamond distance

\varepsilon

in time

\textrm{poly}(n,2^t,1/\varepsilon)

. This resolves a central open question of Mele and Herasymenko under the strongest distance metric. In fact, our algorithm is much more general: we define a property of unitary Gaussianity known as unitary Gaussian dimension and show that our algorithm can learn

n

-mode unitaries of Gaussian dimension at least

2n - O(t)

in time

\textrm{poly}(n,2^t,1/\varepsilon)

. Indeed, this class subsumes unitaries prepared by at most

O(t)

non-Gaussian gates but also includes several unitaries that require up to

2^{O(t)}

non-Gaussian gates to construct. In addition, we give a

\textrm{poly}(n,1/\varepsilon)

-time algorithm to distinguish whether an

n

-mode unitary is of Gaussian dimension at least

k

\varepsilon

-far from all such unitaries in Frobenius distance, promised that one is the case. Along the way, we prove structural results about near-Gaussian fermionic unitaries that are likely to be of independent interest.

n

-fermion wavefunction built from

m

fermionic modes, we provide classical and quantum algorithms returning a Slater determinant with fidelity within

\varepsilon

of maximal in time

m^{\text{poly}(n,1/\varepsilon)}

. We prove matching hardness lower bounds, assuming standard complexity conjectures, along some parameter axes. Given access to quantum copies, we prove this can be accomplished with

\text{poly}(m,n,1/\varepsilon)

copies of

ρ

. We also show that above a fidelity of

2/3

any stationary point is the unique global maximum while below

2/3

the optimization landscape can have spurious stationary points, and hence

2/3

marks a transition point in the optimization landscape for this problem. We apply the algorithm to the Fermi-Hubbard model, extracting the closest Slater determinant from neural quantum state solutions. Together, our results provide algorithmic tools with provable guarantees in understanding fermionic many-body systems with classical or quantum simulation.