bond_update_bug: symmetry-native Basis-Update & Galerkin time integrator (+ TDVP2) - #3
Open
madhavmen wants to merge 16 commits into
Open
Conversation
Add `alice.algorithm.two_site_bug`, a gate-based two-site BUG (Basis-Update & Galerkin) integrator for real- and imaginary-time evolution of an MPS under a nearest-neighbour Hamiltonian. It is built entirely on the existing Alice/Nicole stack (MPS, the AutoMPO interaction list, decomp, the PyTorch backend) and adds no new tensor infrastructure. - run(mps, interactions, opts): even/odd Trotter sweeps of two-site bond gates with a truncated-SVD split that adapts the bond dimension; Lie (first-order) and Strang (second-order) steps, plus imaginary-time cooling. - Bond gates reuse the AutoMPO interaction list: each nearest-neighbour Interaction2Site's leading/terminal tensors are contracted over their operator channel and exponentiated block-wise on the PyTorch backend, preserving the symmetry block structure exactly. - Options (TOML-loadable) and Summary mirror the DMRG interface; Summary records per-step kept and proposed (augmented) bond dimensions. - Tests validate against exact diagonalization: state fidelity, exact norm conservation, U(1) charge conservation, and second-order Trotter scaling. - Docs: API reference pages, nav entry, and changelog. ruff check passes; full suite: 835 passed.
Replace the gate-based two-site update with the faithful rank-adaptive Basis-Update & Galerkin (Ceruti–Kusch–Lubich) K/L/S local bond update. - Vendor the Nicole-native, symmetry-aware faithful-KLS kernel under `_kernel/` (K-augment, L-augment, Galerkin S-step with Krylov expv). - Alice bridge: bond Hamiltonians from the AutoMPO interaction list (bond.py), canonical MPS snapshot + odd/even Strang sweep (scheme.py) driving the kernel. - Remove the gate-based propagator path (gate.py and its tests). - `trunc_thresh` (S-step SVD threshold) controls how far each bond grows; the per-step discarded weight is tracked in `Summary.disc_weights`. - Tests: full-phys domain wall vs exact diagonalization (1 - F ~ 1.9e-8 at dt = 0.05, clean 2nd-order Strang), U(1) Sz and norm conservation, and imaginary-time cooling.
Adds alice.algorithm.tdvp2, a rank-adaptive two-site TDVP integrator (Haegeman et al., arXiv:1408.5056) for real- and imaginary-time evolution of an MPS under a Hamiltonian MPO. Symmetric Strang sweeps of effective-Hamiltonian exponentials with an inverse-free one-site backward correction; the per-bond SVD truncation adapts the bond dimension. Reuses the DMRG environment machinery and 1-/2-site effective-Hamiltonian contractions; the local Krylov expv and the evolution-prefactor handling are self-contained in the package, so it depends only on alice.network and alice.algorithm.dmrg. Validated against exact diagonalization on the Heisenberg chain: state fidelity, exact norm conservation, U(1) total-Sz conservation, imaginary-time cooling, and bond-dimension growth as a domain wall melts. Adds API docs pages and a changelog entry.
A rank-adaptive BUG variant derived from the faithful Ceruti-Kusch-Lubich scheme, differing only in the local bond update: the discarded (orthogonal-complement) projector is applied to the K/L generator before the exponential (project-before), and the basis is grown by a plain direct sum [U0 | Qk] / [V0 ; Ql] with no augmented overlap matrices. The project-before generator is non-Hermitian, so the K/L substep uses a symmetry-preserving tensor Arnoldi exponential; the Hermitian S-step reuses the faithful kernel's tensor Lanczos. The odd/even Trotter sweep, AutoMPO bond Hamiltonians, and Options/Summary are reused unchanged from two_site_bug, so only the per-bond candidate is new. Validated on the symmetric Heisenberg chain against exact diagonalization (state fidelity, O(dt^2) Strang convergence, norm and total-Sz conservation, rank growth, imaginary-time cooling) — 22 tests. Adds the discarded_bug package and tests, docs/api/discarded-bug pages with mkdocs nav and api/index entries (and the previously missing two-site-bug overview section), a changelog entry, and the module exports.
…ture/tdvp2-and-discarded-bug # Conflicts: # docs/api/index.md # docs/getting-started/changelog.md # mkdocs.yml # src/alice/__init__.py
Reimplement the discarded-projector BUG as the MPS specialisation of the Lubich tree-tensor-network BUG. Like two-site TDVP and DMRG it takes a Hamiltonian MPO and exponentiates the two-site effective Hamiltonian with the left/right MPO environments (no Trotter splitting), reusing the DMRG environment machinery; it is inverse-free. A step recursively bisects the chain — the MPS realisation of the reference's balanced-binary-tree Step (the tree is built by recursive bisection of the 1D modes). At each bisection bond the two-site block is evolved once, the K-step and L-step grow the left/right frames with the discarded projector (qr([Theta1_left|U0]) / qr([Theta1_right;V0]), no augmented overlap matrices), and the Galerkin core is the projection of the evolved block, SVD-truncated. Frames are read off the evolved block so a product-state interface grows its genuine rank-2 entanglement. Because every bond is a tree node, the bond dimension grows along the whole chain (the full ballistic light cone) as a domain wall melts. The step is first order in dt with no backward substep; the validated property is the rank growth / light-cone spread. Everything stays in the U(1) block-sparse Nicole representation so the kept rank respects the charge sectors. - Replace the per-bond Trotter/KLS scheme.py with the recursive-bisection sweep.py and the symmetry-aware block_local_update in candidate.py; vendor the Krylov expv. - run(mps, mpo, opts) on a Hamiltonian MPO; Options/Summary mirror the DMRG interface. - Tests cover the full light-cone growth (peaked profile to 2**(L/2)), first-order single-step convergence, that the forward-only error does not shrink with dt, and norm / U(1) / imaginary-time conservation, all vs exact diagonalization. - Update the API docs and changelog to the MPO recursive-bisection scheme.
Replace the recursive-bisection / block-evolution scheme with a single global Basis-Update & Galerkin sweep (mirror of the Julia port): form phi = H*psi via mpo_times_mps, build augmented left/right isometries that keep psi exact and admit only the discarded part (I - U0 U0^dagger) phi per basis matrix (symmetry-blocked), then integrate one Galerkin centre tensor under the two-site effective Hamiltonian. No M/N overlap matrices, no backward substep. Exact at full bond dimension, 2nd order and convergent under truncation. - sweep.py: mpo_times_mps + k_sweep + l_sweep + global_step; the two-site effective apply and centre truncation/assembly helpers folded in. - candidate.py: deleted (block_local_update / bond_snapshot / augmenters removed). - discarded_bug.py, __init__.py: docstrings updated to the global sweep. - tests: validate vs exact diagonalisation.
Brings the rewritten discarded_bug (Lubich TTN-BUG / global discarded-projector sweep, with sweep.py + _krylov.py) onto the tdvp2 branch, replacing the old local two-site project-before candidate (candidate.py / scheme.py removed).
- two_site_bug: 'discarded' variant (project-before discarded-projector KLS update) alongside faithful KLS; pluggable local solvers - discarded_bug: global discarded-projector MPO sweep refinements - tdvp2: imaginary/real-time two-site TDVP - kls: symmetric augmentation completion + discarded candidate kernel - carry SU(2) intertwiners (intw) through to_complex / _krylov casts - opt-in Krylov-depth instrumentation (KRYLOV_LOG) in tdvp2 and two_site_bug - tests: imaginary-time groundstate, discarded variant, local solvers
k_sweep/l_sweep now admit r extra discarded-phi directions (budget = rpsi, Sulz Alg. 5) so the augmented bases can span newly reachable charge sectors; the previous maxdim-rpsi cap starved the augmentation and stalled cooling. global_step re-gauges losslessly both ways (canonical(L-1) then canonical(0)) so both halves reach the true minimal Schmidt rank at every cut instead of the right half blowing up to full Hilbert rank at larger L. maxdim truncation stays confined to the central S-step SVD (off-central re-gauge is trunc=None).
complete_column_basis/complete_row_basis padded each charge sector to its full local dimension d*r, inflating the augmented rank and defeating rank adaptivity. Replace it with the Sulz range basis orth([U0|K1]) at rank <= 2r, plus a minimal random orthonormal seed for reachable charge sectors that neither U0 nor K1 populates -- under U(1) with the opposite frame frozen, K1 stays in U0's sector so the range basis can never open a new one. The padding was doing two jobs and only the first is replaced by the seed: it also completed PARTIALLY-populated sectors to full dimension (measured: max augmented rank 16 = the full local dim, on 160/160 calls, vs 10 with 58/160 below full). That is what made the discarded variant reproduce the exact two-site evolution, and dropping it is deliberate -- padding every sector does not scale. The random fill itself is inert on this model: it fires 130 times and every firing is on a dim-1 sector, where it only picks a phase. Two tests in test_discarded_variant.py encoded the old padded behaviour and are updated to the rank-<=2r contract. Neither tolerance was loosened to pass: - test_matches_faithful_at_full_rank -> test_agrees_with_faithful_at_full_rank. Discarded and faithful span different Galerkin spaces by design, so exact agreement was never the right bar; assert close agreement (1e-8) plus a per-site <Sz_j> profile, since a vec()-based fidelity has misled before. The profile tolerance is DERIVED from the fidelity one rather than picked: a linear observable is first order in the state error while infidelity is second order, so ||dpsi|| ~ sqrt(2*infid) and the profile bound is ~1.4e-4. - test_strang_converges_second_order keeps its >8.0 threshold and simply measures in the asymptotic regime (dt 0.05/0.025 instead of 0.10/0.05). At dt=0.1 the higher-order Trotter terms contaminate the ratio: it reads ~6 at L=4, 6 and 8 alike, and rises monotonically toward 16 as dt shrinks (L=6: 6.33 -> 10.48 -> 14.63 at T=1.0/0.5/0.2) while the binding fraction stays constant. A genuine rank-projection floor would push the ratio DOWN as the Trotter error vanished, not up, so the method remains second order and the fix costs no order. Adds dense_sz_profile to conftest, built from the same _spin_matrices/_embed machinery as dense_total_sz so it is convention-exact against mps_to_vector. Full suite: 33 passed, 1 xfailed (the pre-existing kl_cutoff weight-trim xfail).
…f option Alice is the reference the other implementations are measured against, so this touches BUG code only. Nothing in network/, physics/, dmrg/ or automps was changed. Parked src/alice/algorithm/discarded_bug -> exploratory/global_sweep, with its tests. It is no longer importable from alice.algorithm, and the parked tests cannot run as-is because they import the module that moved -- they are a record, not a working suite, and the README says so. Deleted the four now-unreachable padding helpers: complete_column_basis and complete_row_basis in _kernel/linalg.py, _filter_left_aug_columns and _filter_right_aug_rows in _kernel/kls/augment.py. Verified no call sites first; only explanatory comments in symmetric_completion.py still name them, and those are worth keeping because they record why the padding was replaced. Removed kl_cutoff and kl_cutoff_min_bond from Options, run() and scheme.py. They had no callers anywhere, defaulted to off, and the only tests covering them were the TestKLCutoff class -- one of whose two live tests was a strict xfail documenting that the feature did not work. Default variant is now 'discarded': it is the canonical kernel, the one bond_update_bug! mirrors in BUG-Julia to 4.27e-11 on the L=6 Heisenberg Sz profile. 'faithful' stays available because it is the variant the XX/Heisenberg writeup validated. TWO DEVIATIONS FROM THE PLAN, both toward changing less. Kept solver / solver_substeps. The plan said to remove the non-krylov paths, but they have six passing tests covering midpoint, rk4 and trapezoid. Deleting working, tested functionality is not part of consolidating BUG. Did not split _kernel/krylov.py (579) or nicole_helpers.py (508). Ten files in this repo exceed 500 lines -- network.py is 684, dmrg/environ.py 602, automps.py 599 -- so that rule was never applied here, and splitting two of ten would be arbitrary churn in the reference implementation. Verified against a baseline run of the same suite at 1f2de1b, in a worktree: before 886 passed, 1 failed, 1 xfail -> pytest exit 1 after 867 passed, 0 failed, 0 errors -> pytest exit 0 The 19 fewer tests are accounted for: 13 in the parked global sweep, 3 in TestKLCutoff, 3 in TestGlobalSolvers, and the discarded_bug entry dropped from the imaginary-time method list. The one pre-existing failure was inside the parked suite. One mistake caught by having taken that baseline: parking discarded_bug broke COLLECTION for 13 unrelated test files, including physics and network, because src/alice/__init__.py imports it at package root so any `import alice` failed. Fixed with a one-line edit there. Without the baseline I would have read those as pre-existing.
There is now one BUG scheme, bond_update_bug (the discarded-projector K/L/S sweep), mirroring bond_update_bug! in BUG-Julia. The faithful-KLS variant, the 'variant' option, and the exploratory global-sweep BUG are deleted; TDVP2 is kept as the non-BUG comparison reference. Renames: module two_site_bug/ -> bond_update_bug/, two_site_bug.py -> bond_update_bug.py, _kernel/kls/discarded_candidate.py -> candidate.py, function _discarded_kls_local_bond_candidate -> _kls_local_bond_candidate; test dir and docs (api/bond-update-bug/) renamed to match. All docstrings reference the scheme only as bond_update_bug. Full Alice test suite green.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
bond_update_bug: a single symmetry-native Basis-Update & Galerkin time integratorAdds
alice.algorithm.bond_update_bug, the rank-adaptive Basis-Update & Galerkin(BUG) time integrator for an
MPSevolving under a nearest-neighbour Hamiltonian,plus a two-site TDVP2 integrator for comparison. There is exactly one BUG
scheme, mirrored one-to-one by
bond_update_bug!in the Julia package BUG-Julia.The method
bond_update_bugis the discarded-projector K/L/S bond update, applied as anodd/even Trotter sweep (
:strangsecond order,:liefirst). Each local update:P⊥ = I − U0 U0†to the K/Lgenerator before the exponential (non-Hermitian → symmetry-preserving
tensor Arnoldi);
Û = [U0 | Qk],V̂ = [V0 ; Ql]and acts the isometries directly inthe S-step (
Ŝ0 = Û† Θ0 V̂†) — no overlap matrices;It is symmetry-native (U(1) charge sectors of an Alice
MPS) and inverse-free.Correctness fix carried in this branch
The augmented basis is the Sulz range basis of
[U0 | K1]at rank ≤ 2r, and itis not padded to the ambient local dimension. The earlier
complete_column_basispadding was replaced by a minimal random missing-quantum-number fill (deterministic
seed): the fill still opens an empty-but-reachable charge sector, but no longer
pads a partially populated sector up to
d·r(which does not scale). Measured onL=6 Heisenberg this drops the max augmented rank from 16 (full local dim) to 10 with
no loss of accuracy. The K/L augmentation to 2r uses a symmetric two-way re-gauge.
Validation
total-Sz conservation; second-order Strang convergence; imaginary-time cooling to
the exact ground state.
to 6.33e-8, the residual attributable to the strict Sulz 2r enforcement.
Achievable order
Under the hard Sulz ≤2r bound the symmetric sweep is a clean second-order
integrator; the rank-2r Galerkin update is exact only where both frames are complete,
so composition/extrapolation cannot raise the global order past 2 (completing the
frames recovers 4th order but abandons rank-adaptivity). TDVP2 is kept as the non-BUG
comparison reference.
Companion status — BUG-Julia (real-time, symmetry-toggle)
The Julia companion
bond_update_bug!now carries the developments this PR does notneed, kept in sync as the cross-check reference:
:U1/:none. Alice is symmetry-native and always runs withU(1) sectors; the Julia side adds a dense (no-symmetry) mode for reference runs.
missing-quantum-number fill never fires, so a product initial state would freeze at
χ=1. This is closed by criterion 2 — the Ceruti–Kusch–Lubich residual of the full
2-site update
HΘ, orthogonal to the current frame and capped at Sulz 2r — which isauto-defaulted for a no-symmetry product start (padding to 2r is the fallback). Under
U(1) nothing changes: the missing-fill already grows a product state.
⟨Sz(t)⟩matches the exact free-fermionBessel profile to 8.0e-4 (U(1)) and 5.2e-5 (no-symmetry), norm = 1.000000 and
energy conserved throughout; the two modes agree to 1.8e-3 at T=30 (L=20, dt=0.05,
χ=64). This confirms the symmetric and dense paths are the same real-time propagator.
Notes for review
alice.algorithm.bond_update_bugfollowing thedmrg/tdvp2layout (
Options/Summary/run, GPL headers, docstring conventions).stable— it does not rewrite existing modules.