[challenge]: How cold can a purified tensor-network Anderson impurity solver go?

Released by

Weiyi Guo, Univerisity of Amsterdam

Contact email

w.guo@uva.nl

Method

MPS Based Algorithm

Update at 4th Aug 2026

Author's note: erratum on the research target, and where implicit log‑grid stepping actually belongs

Now that the submission window has closed, I want to correct something in how I framed the research target, because I think the framing sent people toward the wrong bottleneck.

The original wording was:

Can implicit imaginary‑time evolution on a logarithmic grid, combined with adaptive bond expansion, stabilize a fully purified Anderson impurity solver down to $\beta = 100$ at controlled error?

That question packs together a conflation I should have caught.

$\tau$ and $\beta$ are not the same resource

The motivating result—Zima, Stoudenmire, White, Parcollet, and Kaye, arXiv:2606.02930—reaches imaginary time $\tau = 1024$ on a single‑site Anderson impurity model. That is a $\tau$, not a $\beta$. It is pure‑state imaginary‑time propagation,

$$
e^{-\tau H}|\psi_0\rangle,
$$

and large $\tau$ is being used to obtain energy resolution: it resolves the exponentially small Kondo scale in a ground‑state problem.

In a purification, $\beta$ is the physical temperature, and it is paid for in entanglement. The symbols look similar, but they represent different computational currencies. Writing "down to $\beta = 100$" next to a citation whose headline number is $\tau = 1024$ made the target look much closer than it is, and I do not think that was obvious from the issue text.

While I am at it, let me state the number I should have stated in the first place. For a representative parameter set, for example

$$
D = 1,\qquad U = 0.8,\qquad \Gamma = 0.1,\qquad \epsilon_d = -\frac{U}{2},\qquad \mu = 0,
$$

Haldane's formula gives

T_K =\sqrt{\frac{U\Gamma}{2}}\exp\left[-\frac{\pi U}{8\Gamma}+\frac{\pi\Gamma}{2U}\right]\approx 0.0105,

and therefore

$$
\beta_K = T_K^{-1} \approx 95.
$$

So $\beta = 100$ was never an arbitrary round number: to within five per cent, it is the inverse Kondo temperature $\beta_K$. The stretch goal was to reach the scale at which the impurity moment becomes screened. Correspondingly, $\beta = 16$ sits at $T \approx 6,T_K$ and $\beta = 32$ at $T \approx 3,T_K$—equivalently, at $\beta_K/6$ and $\beta_K/3$. The four‑day acceptance line does not yet reach Kondo physics, which is by design, but nobody could tell that from the issue.

Why the log grid is the wrong first lever for purification

The log‑grid plus A‑stable implicit method reduces the number of steps,

$$
O(\tau_{\max}) \longrightarrow O(\log \tau_{\max}).
$$

It says nothing about the cost per step. That distinction is harmless when the cost per step is roughly constant. It is not constant for purification.

Cost structure. For ground‑state projection, $D(\tau)$ rises and then saturates at the ground‑state entanglement, so the total cost is approximately

$$
N_{\text{steps}} \times \mathrm{cost}(D_\infty),
$$

and reducing $N_{\text{steps}}$ reduces the total cost directly. For purification, $D(\beta)$ grows along the trajectory, so the total cost is dominated by the final few steps. Removing $10^4$ cheap early steps while retaining the expensive late steps does not deliver the advertised orders‑of‑magnitude improvement.

Error control. A logarithmic grid deliberately takes $\Delta\tau \sim \tau$: each step is intended to resolve a new decade of energy scales, so the state may change substantially within one step. A‑stability guarantees that the ODE integrator does not blow up. It guarantees nothing about projection onto the fixed‑bond‑dimension MPS manifold, which is the nonlinear component and the part that actually limits the calculation.

This collides directly with deliverable No. 4 of the challenge, which asks for a four‑axis error budget with the time‑step axis separately resolved and reported. A method whose selling point is taking single steps of unbounded size is structurally hostile to that deliverable.

Conditioning. The implicit solve has the form

$$
\left(1 - c,\Delta\tau,H\right)|\psi\rangle = |b\rangle.
$$

Its condition number grows roughly with $\Delta\tau$ times the spectral width of $H$. For an impurity model, the spectral width is extensive in the bath size. The inner‑iteration count can therefore consume part of the outer step‑count saving.

The steelman. None of this makes implicit log‑grid stepping wrong for purification; it makes it premature. If the ancilla gauge is fixed first, the bath's thermal contribution can collapse and the problem becomes ground‑state‑like again. The bond dimension may then saturate, and the log‑grid argument recovers much of its force.

The correct claim is about ordering:

$$
\text{gauge first} ;\longrightarrow; \text{integrator second}.
$$

My phrasing invited people to do this in the opposite order.

What I left out of the problem statement entirely

The purification is not unique. Since $H$ acts only on the physical sites, any unitary acting exclusively on the ancilla space commutes with the Hamiltonian and leaves all physical observables invariant. The entanglement of $|\Psi_\beta\rangle$ is therefore partly a property of the gauge choice, not solely of $\rho_\beta$. It is also the single highest‑leverage knob in the problem.

Three important possibilities were absent from the issue:

• Thermofield / Bogoliubov transformation. For a non‑interacting bath, each bath mode and its ancilla can be rotated by a $\beta$‑dependent angle so that the bath's thermal state becomes the exact vacuum of new modes—a product state at every $\beta$. In a single‑impurity model, the interaction $U$ acts on exactly one site while the rest is free, so this can remove most of the cost. See de Vega and Bañuls, Phys. Rev. A 92, 052116 (2015).

• Disentanglers. Ancilla unitaries can be applied during the evolution to minimize entanglement. See Hauschild et al. Phys. Rev. B 98, 235163 (2018).

• Geometry. An impurity is a hub, not a link in a chain. Fork tensor‑product states place the impurity at a junction with one branch per spin or orbital bath. See Bauernfeind, et al. , Phys. Rev. X 7, 031013 (2017).

An MPS also has no natural room for a dangling ancilla leaf attached to every site: that structure is a comb, not a chain. One must either fuse the ancilla into the site tensor, raising the local dimension from $d$ to $d^2$; flatten the ancillas into the chain, roughly doubling its length and, with fermionic ancillas, dragging Jordan–Wigner strings through them; or leave the MPS geometry in favour of a tree. (There is also a framing that dissolves the question rather than answering it; see the XTRG/tanTRG section below.)

The third option is the one I would most like to see pursued, because for an impurity problem a tree is not a workaround but the natural topology — the impurity is a hub, and the impurity–bath entanglement is intrinsically multiscale. Two recent anchors. Grundner et al. , PRB 109, 155124 (2024), use a T3NS as the working substrate for multi-orbital impurity Green's functions and tame entanglement growth by deforming the evolution onto complex-time contours. Zhan et al. , PRB 113, 195144 (2026), decompose the bath itself onto a Cayley tree precisely because that geometry natively captures the multiscale entanglement structure of impurity–bath systems. Grundner's thesis, Tensor Network Impurity Solvers: Simulating Quantum Materials, is the natural review for complex time solver on tree topology and imaginary time evolution with purified thermal state.

I want to be explicit that neither of these is yet a finite-temperature purification on a tree: the first is complex-time, the second real-time, and neither carries ancillas. Attaching the ancillas as leaves of such a tree — so that the purification's gauge freedom and the bath's multiscale geometry are handled by one and the same network — is, as far as I am aware, open. That is precisely why I now think it is a better target for this challenge than the integrator question I originally posed.

Where implicit log‑grid stepping does belong: METTS

I still think the method is excellent. I simply pointed it at the wrong algorithm. Its natural home in this problem space is METTS, and the fit is close to perfect:

1. Each METTS sample is exactly the paper's setting. One evolves

$$ e^{-\beta H/2}|i\rangle$$

from a classical product state. This is pure‑state imaginary‑time propagation of a sum of decaying exponentials, which is the setting on which the log‑grid argument is based.

2. The cost per step is roughly flat. The point of "minimally entangled" thermal states is that a typical sample carries much less entanglement than a purification, so reducing the number of steps translates more directly into wall‑clock savings.

3. Only the endpoint matters. Each sample is a one‑shot calculation: obtain the state at $\beta/2$. There is no requirement that every intermediate state be accurate, which is precisely what a log grid cannot guarantee.

4. A‑stability matters more, not less. A product initial state has support across the entire spectrum, so the stability limit of an explicit scheme is set by $E_{\max}$. This is exactly the regime in which A‑stable implicit stepping earns its keep.

5. There is a statistical error floor. METTS already carries Monte Carlo noise. The systematic time‑stepping error only needs to remain below that floor. Trading fine‑grained per‑step control for speed is therefore reasonable, whereas it is much less appropriate for deterministic purification with a mandated error budget.

6. It is embarrassingly parallel. The integrator can be validated by rerunning a subset of samples on a fine uniform grid. Two METTS references already in circulation for this challenge are Bauernfeind et al., Phys. Rev. B 105, 195107 (2022) on Matsubara Green functions, and Cao et al., Phys. Rev. B 109, 245113 (2024) on a finite‑temperature METTS impurity solver on ForkTPS topology.

The field has already run this experiment: XTRG versus tanTRG

There is a framing of the ancilla‑geometry question that dissolves it rather than answering it, and it carries a lesson bearing directly on the erratum above.

Write the purification out against the maximally entangled reference state $|I\rangle = \sum_i |i\rangle|i\rangle$:

$$
|\Psi_\beta\rangle = \sum_{ij}\left(e^{-\beta K/2}\right)_{ij},|i\rangle_p,|j\rangle_a .
$$

The coefficient tensor is the matrix $e^{-\beta K/2}$. An MPS carrying one ancilla leg per site and an MPO representing $\rho^{1/2}$ are the same object written in two notations, related by Choi–Jamiołkowski. In operator space the question "where do I attach the leaf?" cannot even be posed, because the leg is already an index of the tensor. This is not a fourth option—it is the $d \to d^2$ option in its natural coordinates. But working there is what makes the following visible.

XTRG has provided a logarithmic grid for thermal states since 2018. Chen et al. Phys. Rev. X 8, 031082 (2018), evolve $\rho$ as an MPO by

$$
\rho(2\beta) = \rho(\beta)^2,
$$

doubling $\beta$ at every step and therefore reaching $\beta = 2^n$ in $n$ steps. The reason this is available for thermal states but not for the pure‑state problem is structural: $e^{-\beta H}$ is an operator and obeys the semigroup property, so it can be squared. A state $e^{-\tau H}|\psi_0\rangle$ cannot be squared, which is precisely why Zima et al. had to reach for A‑stable implicit machinery in order to obtain the same logarithmic saving. For the thermal branch of this challenge, the logarithmic grid was never the missing piece.

And then the same group went the other way. Five years later, Li et al., Phys. Rev. Lett. 130, 226502 (2023), introduced tanTRG, which achieves "an optimal evolution of the density operator" with "a mild $O(D^3)$ complexity" by tangent‑space projection—remaining on the $D$‑dimensional manifold throughout and never inflating, at the price of returning to fine $\beta$ steps. Its targets are width‑8 cylinders and $10 \times 10$ lattices of the 2D Hubbard model, down to $T \approx 1/24$.

Why abandon an exponential saving in step count? Because squaring inflates the bond dimension, $D \to D^2$, before compression. When $D$ is modest, the step saving dominates and XTRG wins. When $D$ must be large—2D, fermions, low $T$—the scaling in $D$ dominates, the intermediate $D^2$ object becomes the memory wall, and a factor of $\log_2 \beta \approx 7$ in step count cannot compensate a factor of $D$ in cost per step.

That is the same trade‑off as this erratum, decided by the same community, on thermal states, before the paper I cited was written:

	Regime	Binding constraint	Choice
Zima et al. (2026)	pure state, $D$ saturates	step count	fewer steps ✓
XTRG (2018)	thermal, 1D, small $D$	step count	fewer steps ✓
tanTRG (2023)	thermal, 2D, large $D$	$D$	cheaper steps ✓
This challenge, as posed	thermal impurity, $D(\beta)$ growing	$D$	fewer steps ✗

Two caveats I would rather record than paper over. First, the abstracts do not state XTRG's cost exponent explicitly; the inflate‑then‑compress mechanism is standard MPO arithmetic, but anyone quoting a specific power of $D$ should check the body text. Second, squaring yields $\rho$ only at dyadic $\beta$, whereas this challenge requires $G(\tau)$ on a $\tau$ grid inside $[0,\beta]$; assembling those from dyadic pieces is a genuine engineering problem, not a drop‑in substitution. XTRG is not a ready‑made impurity solver, and I am pointing at the lesson rather than at the code.

And the caveat that matters most: which axis binds is itself a consequence of the gauge. In a naive interleaved purification, $D(\beta)$ grows and tanTRG‑style reasoning applies. After a thermofield transformation the bath's contribution collapses, $D$ may saturate, and the step‑count argument recovers. This is the same "gauge first, integrator second" ordering as above, now with independent support from the thermal tensor‑network literature.

For what it is worth, the XTRG paper is already present in this repository's own knowledge base, under .knowledge/literature/ltrg/.

Revised guidance for follow‑up work

If anyone picks this up:

• For the purification track, stay with TEBD, TDVP, MPO-$W^{\mathrm{II}}$, or Krylov time evolution, or move into operator space and use a tanTRG‑style tangent‑space evolution of $\rho$. Use uniform or mildly adaptive time steps, with per‑step truncation errors monitored and reported. That is what makes the four‑axis error budget meaningful. Direct-star geometry produces long-range impurity–bath couplings, so plain nearest-neighbour TEBD cannot be applied directly. One may use a chain mapping, swap-gate Trotter/TEBD as in Bauernfeind et al. Phys. Rev. X 7, 031013 (2017), or a full-MPO method such as $W^{\mathrm{II}}$/TDVP.

• Spend the first effort on the gauge, not the integrator. As we discussed a lot above, either into Ancilla leaf site of tree tensor network or otherwise like thermofield transformation. I would now consider a thermofield/rotated-gauge/tree-tensor ‑based $\beta = 16$ result with a clean error budget a stronger submission than a naive‑gauge $\beta = 16$ result.

• Treat implicit log‑grid stepping as a METTS project. It is a good project, and I would be glad to see it attempted in that form.

• The $\beta = 100$ frontier still stands. Read it as "reach $T_K$", and expect the gauge choice, but not the integrator, to be what gets you there. Both routes are now explicitly in scope: a deterministic purification and a METTS sampling approach are different answers to the same question, and I would be glad to see either, or a comparison of the two.

Post‑mortem

For what it is worth, the submission that got furthest hit exactly the wall described above: a textbook interleaved‑identity purification with fermionic Ancillas, no thermofield transformation, no disentangler, and two‑site TDVP. Its $\beta = 16$ runs saturated the $D = 128$ and $D = 256$ cells, reaching approximately $D \approx 254\text{–}382$ against a $D = 512$ cap before they were stopped. A full‑acceptance job ran out of memory at 31 GB against a 32 GB allocation.

The engineering was, if anything, more careful than the challenge required: atomic hash‑bound checkpoints, fail‑closed capability gates, and a time-step sweep that honestly reported non‑monotonic convergence—with $\Delta t = 0.01$ performing worse than $\Delta t = 0.02$—rather than quietly claiming the passing point. That non‑monotonicity is itself informative: at fixed cutoff and bond dimension, halving $\Delta t$ doubles the number of truncation events, so the accumulated projection error can outgrow the reduced integration error. The time‑step and bond‑dimension axes of the error budget are not independent, and a smaller step buys nothing unless $D$ is raised with it.

That is not a competence failure. It is what happens when the problem statement names the wrong bottleneck and the folklore fix is nowhere in the text. That one is on me, and this note is the correction.

Original Challenge issue on 15th July 2026

The question

Build a trustworthy finite-temperature tensor-network solver for the single-orbital Anderson impurity model with a continuous fermionic bath, using deterministic purification — then find out how cold it can go.

Central research question: can implicit imaginary-time evolution on a logarithmic grid, combined with adaptive bond expansion, stabilize a fully purified Anderson impurity solver down to $\beta = 100$ at controlled error?

The four-day core is deliberately smaller than the research question: one solver, one bath, verified against independent references, with an explicit error budget. The colder regimes and the bosonic bath are layered on top.

Model

Spinful single-impurity Anderson model at particle-hole symmetry, $\epsilon_d = -U/2$, $\mu = 0$:

K = \sum_\sigma (\epsilon_d-\mu)\,n_{d\sigma}
  + U n_{d\uparrow}n_{d\downarrow}
  + \sum_{k\sigma}\epsilon_k c^\dagger_{k\sigma}c_{k\sigma}
  + \sum_{k\sigma}\left(V_k d^\dagger_\sigma c_{k\sigma}+\mathrm{h.c.}\right),

with the semicircular hybridization

\Gamma_f(\omega)=\Gamma\sqrt{1-(\omega/D)^2}\,\Theta(D-|\omega|),
\qquad D = 1,\; U = 0.8,\; \Gamma = 0.1 .

Core temperatures $\beta \in {1, 4, 16, 32}$; $\beta = 100$ is the stretch frontier. A submission must serialize the hybridization it actually uses on a common frequency grid, so that different bath fits remain comparable.

Thermal-state convention. The baseline is deterministic purification of the complete interacting finite Hamiltonian obtained after bath discretization,

|\Psi_\beta\rangle=\left(e^{-\beta K/2}\otimes I_a\right)|I\rangle,

with impurity, bath sites and thermal ancillas in one purified state. METTS is not required. A thermofield rotation may be used as a quadratic bath-basis transformation, but not as a replacement for preparing the interacting Gibbs state.

Core deliverable (four-day acceptance line): a trusted fermionic solver

1. Bath fitting. Discretize/fit $\Gamma_f$ with a documented convergence study (number of bath levels, low-frequency resolution).
2. Finite-bath validation. For a small number of bath levels per spin, match exact diagonalization on $n_d$, double occupancy $D = \langle n_{d\uparrow} n_{d\downarrow}\rangle$ and $G(\tau)$ to $10^{-6}$ maximum error, after separately converging time step and bond truncation.
3. Continuous-bath run. One calculation of the reference model at $\beta = 16$ or $\beta = 32$, cross-checked against CT-HYB (e.g. TRIQS) or a fermionic GTEMPO reference.
4. Error budget. Report $n_d$, double occupancy and $G(\tau)$ on a shared tau grid, with the error split into bath discretization, chain length, bond truncation, and time-step/residual contributions — plus wall time, peak memory and per-bond dimensions.

An honest convergence-or-failure report at $\beta = 16 ~ 32$, regenerated by an automated script, is a valid school deliverable.

Research target: implicit logarithmic evolution with adaptive bonds

Uniform-step two-site TDVP is the reference, but the low-temperature target motivates A-stable implicit stepping on a logarithmic imaginary-time grid: each step solves $A x = b$ with $A = I + (h/2)K$, $b = (I - (h/2)K)\psi_n$.

• Bootstrap. The purified $\beta = 0$ state has minimal bond dimensions, so at initialization only, use a few Krylov-informed global subspace expansion passes to open the bonds.
• Main loop. Solve each implicit step inside the current manifold. When the tangential residual has converged but the full residual $r = b - A x$ has not, expand only the responsible bonds via local residual-driven controlled bond expansion. This residual-driven scheme is an open algorithmic problem, not a solved recipe.
• Frontier. $\beta = 100$ is the stretch target, not a gate. Bracketing the coldest reachable $\beta$ at $10^{-3}$ common observable error is itself the research result.
• Comparison. Plot $G(\tau)$ error versus wall time, peak memory and maximum bond dimension for uniform TDVP, the implicit-log scheme and, where available, GTEMPO, on identical physical inputs.

Bosonic extension (stretch / future work)

The first layer is a single local Einstein phonon coupled to the impurity charge (Anderson-Holstein with one mode), validated against finite-mode exact diagonalization including $\chi_{nn}(\tau)$. The continuous phonon bath $J_b(\omega)=2\alpha\omega e^{-\omega/\omega_c}$, projected purification, and a mixed electron-phonon GTEMPO comparison are explicitly future work: they introduce boson cutoff, bosonic bath discretization, projected-purification and mixed-bath benchmark errors all at once, which would obscure the fermionic solver milestone within four days.

Out of scope

DMFT/EDMFT self-consistency, multiorbital or spin-orbit impurities, real-time dynamics, analytic continuation, and implementing METTS. However, we note that a tree tensor network is a natural network geometry to adapt those novel applications in terms of DMFT with mulitiorbital spatial disentanglement initialised by Bauernfeind et al. and Cao et al..

Key references

• Bauernfeind et al., Minimally Entangled Typical Thermal States Algorithms for Finite Temperature Matsubara Green Functions (2022).
• Zima et al., Fast Tensor Network Imaginary Time Evolution by Implicit Stepping on Logarithmic Grids (2026).
• Li, Gleis and von Delft, Time-dependent variational principle with controlled bond expansion for matrix product
  states (PRL 2024).
• Yang and White, Time Dependent Variational Principle with Ancillary Krylov Subspace (PRB 2020).
• Kohn and Santoro, Efficient mapping for Anderson impurity problems with matrix product states (2020) — used here only as a quadratic bath representation.
• Chen, Gu and Guo, Tensor Network Algorithm to Solve Polaron Impurity Problems (2025) — the GTEMPO comparator; its mixed-bath setting is the future extension.
• Grundner, Tensor Network Impurity Solvers: Simulating Quantum Materials, PhD thesis, LMU Munich (2025) — global-then-local subspace expansion.
• TRIQS/ctseg one-boson retarded-interaction tutorial -- compact $G(\tau)$/ $\chi_{nn}(\tau)$ cross-check for the finite-mode cell.
• Harnessing Quantum 2026 challenge-idea page and participant guide.

No repository files, tests, or entry points are named. Start with the corrected research target and its METTS, purification, and tensor-network geometry sections, then review the cited methods before choosing a concrete implementation path. Done would require a validated solver study against the challenge's stated error and temperature goals.