Double Factorization

English · 简体中文

Category C2 · Module oracq.algorithms.input_model.lowrank · Stage V1

Overview

The double factorization (DF) of the quantum chemistry electronic Hamiltonian compresses the two-body integral tensor into a low-rank form

\[ H = \mathrm{scalar}\cdot I + \sum_r U_r\, \mathrm{diag}(g_r)\, U_r^\dagger , \]

each rank term being a Hermitian matrix whose diagonal spectrum \(g_r\) is conjugated by the rotation \(U_r\). The implementation follows the “low-rank tensor → LCU → BE” pipeline of Berry et al. 2019 and von Burg et al. 2021 (literature citations in the module docstring): the integral tensor is given as classical input data, this module performs no real chemistry integral computation, and the resulting block encoding can enter qubitization directly.

Interface and input model

DoubleFactorization(scalar, rotations, spectra)
DoubleFactorization.from_symmetric(scalar, rotations, factors)
diagonalize_symmetric(matrix, *, tolerance=1e-12, max_sweeps=100)
double_factorized_encoding(df, *, name=None)

API entry points: DoubleFactorization, diagonalize_symmetric, double_factorized_encoding

• DoubleFactorization: the DF input model (input model CP — the tensor data is given directly as classical parameters, i.e. the “low-rank tensor (HAM) → BE” pipeline in algorithm-coverage). rotations are explicit small unitary matrices (dimensions powers of two within 2..32, unitarity checked by column orthogonality, tolerance 1e-9), and spectra are the corresponding real spectra; the dataclass attributes width gives the number of target qubits and rank the number of rank terms.

• from_symmetric: the classical preprocessing entry point of physical DF — takes explicit unitaries \(U_r\) and real symmetric \(G_r\), and folds the eigenvector matrix of \(G_r = V_r\,\mathrm{diag}(g_r)\,V_r^{\mathsf T}\) into the rotation (\(U_r V_r\)).

• diagonalize_symmetric: Jacobi eigendecomposition of a real symmetric matrix, returning (eigenvalues, eigenvector matrix); non-symmetric input raises ValidationError at generation time.

• double_factorized_encoding: assembles the LCU block encoding (the name parameter is retained but currently takes no part in naming).

Returns a BlockEncoding with module attributes:

Attribute	Meaning
be_alpha	\(\lvert\mathrm{scalar}\rvert + \sum_r \lVert g_r\rVert_1\) (sum of the spectra’s 1-norms)
be_form	"double_factorization"
df_rank / lcu_terms	number of rank terms / number of nonzero terms of the outer LCU
df_lambda	a λ report under the same accounting, equal in value to be_alpha

Implementation notes

The outer LCU’s PREPARE runs over the rank index \(r\) with weights \(\propto \lVert g_r\rVert_1\); a nonzero scalar joins the same LCU as (scalar, identity), and all-zero coefficients (including all-zero spectra) raise ValidationError. Each rank term’s term block encoding is \(U_r^\dagger\cdot\mathrm{diag}\cdot U_r\): the diagonal part applies a controlled \(R_y(\theta_t)\) on a single-bit signal for every basis state \(\lvert t\rangle\), with \(\cos(\theta_t/2) = g_t/\alpha_r\) and \(\alpha_r = \sum_p\lvert g_p\rvert\) taken as the spectral 1-norm, so the (0,0) block is exactly \(\mathrm{diag}(g_r)/\alpha_r\); the synthesized \(U_r\) is applied on the two sides to complete the conjugation.

\(U_r\) is synthesized by two-level decomposition: columns are eliminated into diagonal phases and replayed in reverse order; each two-level unitary is realized by a ZYZ decomposition (\(e^{i\varphi}R_z(\alpha)R_y(\beta)R_z(\gamma)\)) plus multi-controlled X transpositions along a Gray path. Register layout: the term block encoding is target(n) | signal(1), and the outer LCU’s signal is the selector (as many bits as the term count requires) concatenated with work(1).

Applicability boundary: the explicit-matrix path supports only validation-sized instances of dimension 2..32 (at most 5 qubits); at larger scales the spectra and rotation angle tables of \(U_r\)/\(G_r\) should switch to QRAM data binding (see the three-layer paradigm of prepare_select). In terms of normalization, the DF α of a single-bit positive definite \(g\) equals \(\mathrm{tr}(g)\) and introduces no off-diagonal redundancy from the Pauli expansion (tightness witness in the validation section).

Validation approach

Category C2 (approximately continuous semantics; the acceptance criteria are in ../development/validation-plan.md §2): the block encoding’s zero-signal block must agree with the target Hamiltonian within tolerance. Three layers of evidence:

• Structure: tests/core/test_lowrank.py:DiagonalizeSymmetricTests.test_reconstructs_factor (\(V\,\mathrm{diag}(\lambda)\,V^{\mathsf T}\) reconstruction and eigenvector orthogonality, places = 10) and test_rejects_non_symmetric; the construction attribute assertions of DoubleFactorizationTests (be_form, df.rank == 1 / df.width == 1); ThcTests.test_invalid_inputs_fail_at_generation covers the DF negative cases (empty rank terms, mismatched spectrum lengths, non-unitary rotations, non-numeric scalar, inconsistent counts in from_symmetric).

• Numerical: DoubleFactorizationTests.test_single_rank_block_equals_hamiltonian — single-rank DF cross-checked column by column against \(G\) via assert_block_equals (places = 9); test_scalar_and_rotated_ranks — a scalar and a Hadamard-rotated second rank, with the expected matrix assembled independently inside the test before cross-checking, and α cross-checked against \(\lvert\mathrm{scalar}\rvert + \sum_r\lVert g_r\rVert_1\) (places = 10); test_feeds_qubitization_walk — after the product enters transforms.qubitization_walk there are no unbound slots and the registers are ["target", "signal"].

• Binding: DoubleFactorizationTests.test_matches_pauli_encoding_block — the DF encoding and the Pauli LCU encoding of the same Hamiltonian are cross-checked column by column across BEs via the block_column helper (places = 9), covering the consistency of two independent encoding paths.

α tightness witnesses (added in V1): test_alpha_matches_closed_form_eigenvalues independently computes \(\sum\lvert\lambda\rvert\) from the 2×2 closed-form eigenvalues \(\lambda_\pm = (t \pm \sqrt{t^2-4d})/2\) (\(t\) the trace, \(d\) the determinant) and cross-checks it against be.alpha (places = 10, measured 2.0); test_df_alpha_tighter_than_pauli takes an off-diagonal-dominant positive definite \(g\) and measures DF α = 1.4 ≤ Pauli LCU α = 1.5 (theoretical condition: on single-bit PSD matrices Pauli α − DF α = \(\lvert g_{01}\rvert + \lvert\Delta/2\rvert - \mathrm{tr}/2\), positive in the off-diagonal-dominant case).

Known gaps and planned stages

No known gaps; the stage V1 witnesses are complete (block cross-checks + α tightness against the closed-form eigenvalues and the Pauli comparison + cross-encoding binding consistency + the qubitization hookup).

Related links

• Source: src/oracq/algorithms/input_model/lowrank.py

• Algorithms in the same module: THC block encoding

• API reference: Chemical low-rank decomposition block encodings

• Concepts: Oracles and operator representations

• Validation matrix: Validation coverage matrix

Numerical validation

Paper-grade numerical validation script: tests/verification/verify_blockencoding.py (real-backend execution, no mock substitutes, no skips; all 73 cases passed on 2026-09-16), artifact out/verification/blockencoding.json. This page corresponds to diagonalize-symmetric-*, df-encoding-*, and diagonal-be-*, 11 cases in total (the diagonal block encoding _diagonal_encoding is exactly the kernel of a DF rank term).

Experiment design:

1. Diagonal BE across scales (_diagonal_encoding, the diagonal BE named in the paper): three spectrum groups (mixed signs, containing zeros, generated from trigonometric functions) with dim ∈ {2, 4, 8} go through OriginIR-ext + UniQC to_matrix full-unitary extraction (the (0,0) block compared against \(\mathrm{diag}(g)/\alpha\), with leakage into the failure branch reported), and dim ∈ {16, 32, 64} go through reference + rir-pysparq column by column; there is also a public assembly-consistency case — the same diagonal spectrum yields element-wise identical blocks via double_factorized_encoding (identity rotations, two ranks with the same spectrum) and via the direct diagonal BE.

2. Jacobi eigendecomposition: random real symmetric matrices with dim ∈ {2, 4, 8, 16, 32}, eigenvalues compared against numpy.linalg.eigh, plus checks of the reconstruction \(V\mathrm{diag}(\lambda)V^{\mathsf T}\) and orthogonality.

3. DF block encoding across ranks and scales: 2×2 double rank + scalar (from_symmetric preprocessing, full-unitary extraction), 4×4 three ranks with random orthogonal rotations (full-unitary extraction), 8×8 double rank (reference + rir); the expected matrices are assembled independently with numpy, and α is compared against the closed form \(|\mathrm{scalar}|+\sum_r\lVert g_r\rVert_1\).

Case	Scale	Paths	Metric	Value
diagonal-be-unitary-d2/d4/d8	dim 2/4/8, α=2.0/3.25/4.0	originir-ext + to_matrix, three-backend cross	max_error	5.6e-17 / 5.6e-17 / 1.1e-16
diagonal-be-wide-d16/d32/d64	dim 16/32/64	reference + rir	max_error	1.1e-16 / 1.0e-16 / 1.1e-16
diagonal-be-public-df-d4	dim 4 (two paths)	reference	max_error	5.6e-17
diagonalize-symmetric-d2..d32	5 groups	independent classical cross-check	reconstruction / eigenvalue / orthogonality error	≤6.2e-13 / ≤5.7e-14 / ≤7.3e-15
df-encoding-rank2-scalar-d2	rank 2 + scalar, α=5.05	to_matrix + three paths	max_error / α − closed form	4.4e-16 / −8.9e-16
df-encoding-rank3-d4	rank 3, α=8.24	to_matrix + reference	max_error	1.6e-15
df-encoding-rank2-d8	rank 2 + scalar 0.7	reference + rir	max_error	1.2e-15

Reproduction command:

PYTHONPATH=src /path/to/backend/python tests/verification/verify_blockencoding.py

Artifacts: out/verification/blockencoding.json.