Skip to content

diagnose() at d >= 8: the exponential action now propagates the trajectory — D18 and the resolvent/non-normality layers are the remaining bottlenecks #178

Open

Description

@marcohost33-maker

Context

PR #176 (issue #162) moves relaxation-trajectory propagation to the exponential action (expm_multiply) for n = d^2 >= 64. That speeds up the propagation itself 30–46× at n = 256–1024. It does not make diagnose() fast at those sizes, and the PR says so.

Measured (single-threaded BLAS, damped transverse-field Ising chain, bootstrap_B=20)

d diagnose() total relaxation layer of which propagation non-normality resolvent transient
8 5.9 s 2.11 s (fits/bootstrap) 0.05 s 2.00 s 1.17 s 0.60 s
16 > 25 min (killed) — — — — —

The profile at d = 8 was also noisy with ResolventConvergenceWarnings: the power iteration in resolvent_norm did not converge in 200 steps at many frequencies.

Candidates, ordered by value/effort

  1. D18, lep.initial_state_sensitivity. It materialises sla.expm(L t_eval) once and then applies it only to n_samples vectors. This is the same pattern Trajectory backend: use matrix-exponential action instead of materialising expm #162 removed. A block action expm_multiply(L*t_eval, B) with B = [vec(rho_k)] fits directly.
    • One caveat: for n0 columns, the deterministic exact-norm branch of SciPy requires ||A||_1 <= 63.36 / n0 (condition 3.13). The sub-stepping from numerics.propagation has to scale its threshold by n0, or the randomised onenormest path, which touches global np.random state, comes back.
  2. Resolvent / non-normality layers. They dominate at d = 8 and are the likely reason d = 16 does not finish. This needs its own profile at d = 16 (which frequencies, how many resolvent_norm calls, whether the power iteration hits its cap) before choosing a remedy.
  3. transient.py. It computes svdvals(expm(L t)) per grid point. The propagator norm genuinely needs more than one action, so the action does not apply directly here. A Krylov or randomised norm estimator would be a methodology change and belongs in its own issue.

Acceptance for item 1

  • D18 is unchanged bit-for-bit for n < 64.
  • For n >= 64, dense and action agree within the declared tolerance.
  • Global RNG state stays untouched: D18 draws its samples from its own default_rng, and onenormest must not touch the global state either.
  • The backend is recorded, as in RelaxationResult.trajectory_backend.

Refs #162, #176.

Activity

  1. added 17 commits that reference this issue on Sep 29, 2026
  2. marcohost33-maker commented on Sep 29, 2026

    @marcohost33-maker
    OwnerAuthor

    2026-09-30 D18 implementation / crossover falsification

    Draft PR #181 implements item 1 on top of merged #176: D18 now has a single-time block exponential primitive and can evolve the Haar RHS ensemble with expm_multiply(L*t, B) without materialising the propagator.

    Two details from the implementation changed the naive plan:

    1. Condition (3.13) must scale with RHS width. SciPy's exact-norm eligibility uses n0 = B.shape[1]; using feat(#162): propagate relaxation trajectories by the exponential action #176's one-vector substep threshold for a 10-column Haar block would re-enter randomised onenormest and consume global np.random. perf(#178): propagate D18 Haar ensemble by deterministic block action #181 therefore uses the same safety margin divided by n_rhs, with a regression probe on every block call plus an explicit global-RNG invariance test.
    2. D18 cannot reuse Trajectory backend: use matrix-exponential action instead of materialising expm #162's n=64 auto crossover. D18 forms exactly one dense expm and applies it to a small block; the trajectory benchmark formed many exponentials. An exploratory check on the same damped-Ising family showed that block action can win at n=64, lose at n=256, and win clearly again by n=1024 on the test environment. That falsifies a simple monotone transfer of the trajectory threshold. perf(#178): propagate D18 Haar ensemble by deterministic block action #181 therefore keeps dense auto below a conservative BLOCK_ACTION_MIN_DIM=512; forced action remains available for differential evidence.

    The PR also records LepResult.initial_state_backend, pins small-system legacy arithmetic, dense/action numerical agreement, report-JSON persistence, and keeps D18 advisory (no classifier change).

    Current exact-head gates on #181: Quality Contract, QuTiP and CodeQL are green; full Python matrix and reusable pytest are still executing. Do not treat this comment as merge evidence until those exact-head runs finish.

    After item 1 lands, item 2 should start from a fresh d=16 profile rather than assuming the d=8 proportions remain representative: record frequency count, non-converged resolvent_norm calls, iteration counts/cap hits and time by non-normality/resolvent subroutine before choosing an algorithmic replacement.

  3. added 8 commits that reference this issue on Sep 30, 2026
  4. marcohost33-maker commented on Sep 30, 2026

    @marcohost33-maker
    OwnerAuthor

    Item 1 (D18) ist als Draft-PR #182 umgesetzt. Kleine Systeme n=d^2<64 behalten den historischen Dense-Pfad exakt. Ab n>=64 werden die Haar-Zustände in derselben RNG-Reihenfolge als Block mit expm_multiply propagiert; die Condition-(3.13)-Schwelle wird durch n_rhs geteilt, damit SciPy auf dem deterministischen Exact-Norm-Pfad bleibt und globales np.random unberührt bleibt. LepResult.initial_state_backend dokumentiert den Backend-Wechsel additiv. Neue Tests vergleichen gegen die Dense-Referenz und pinnen RNG-Isolation. Ruff- und mypy-Befunde der ersten CI-Runden wurden bereits behoben; Resolvent/Non-Normality bleiben bewusst separater nächster Profiling-Slice.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Metadata

Metadata

Assignees

No one assigned

    Labels

    enhancementNew feature or request

    Projects

    No projects

      Milestone

      No milestone

      Relationships

      None yet

      Development

      No branches or pull requests

      Issue actions