09Describing interacting electrons classically

Optimize the many-electron coefficients, not only the orbitals

At the stretched N₂ geometry, one optimized determinant can remain an inadequate description even after its SCF residual is small and its allowed stability tests pass. Orbitals describe one-electron shapes. Correlation also requires flexibility in how several electrons occupy those shapes together.

We now hold a chosen orbital representation fixed and write

$D_I$ is a Slater determinant, $\mathcal D$ is the chosen determinant collection, and $c_I$ is its amplitude. Unlike the MO coefficients $C_{\mu p}$, these amplitudes describe a many-electron superposition. Configuration interaction (CI) optimizes them by solving the Hamiltonian eigenvalue problem in that determinant space.

The same principle underlies full CI, CASCI, selected CI, and the projected solver in this repository. They differ primarily in which determinant space they provide and how they search or diagonalize it. They do not all optimize the orbitals; that is the MCSCF question in Lesson 11.

Our running molecular model remains neutral singlet N₂ at 1.85 Å with 28 spatial orbitals in cc-pVDZ. Nothing in this lesson claims a newly executed calculation. Small matrices are illustrative, and separate molecular examples are explicitly identified.

Configuration records need an orbital basis and a consistent sector

For an electron spin-half Configuration, each spatial-orbital slot has one of the values 0, u, d, or 2: empty, alpha-only, beta-only, or doubly occupied. Configuration.from_spin_half_string("2du0") records occupations, not spatial orbital shapes or amplitudes. A determinant acquires its physical meaning only when paired with a particular ordered orbital basis and fermionic ordering convention.

A ConfigurationSet associates a collection of configurations with Orbitals and a single-particle sector name such as "electrons". That sector label identifies the kind of particles/basis being represented; it is not itself a constraint fixing electron number or total spin. For active-space calculations, configuration positions refer to the active-orbital ordering, not necessarily to consecutive positions in the full original MO matrix. Inactive filled orbitals are represented through the orbital metadata and Hamiltonian, rather than repeated in every active configuration string.

The packed representation uses two bits per spin-half spatial slot, or one bit per slot for a generic bitstring configuration. capacity() reflects byte-padded storage; ConfigurationSet.num_modes() derives the meaningful non-padded count from its orbitals. Here a spin-half slot can represent two electron spin modes; it is not automatically one occupation qubit. Conversion to MACIS bitsets places alpha occupation in the low half and beta occupation in the high half. A string’s displayed length, packed capacity, active spatial-orbital count, and qubit count must not be treated as synonyms.

The concrete ConfigurationSet validation checks consistent capacities and, for a nonempty defined active window, sufficient capacity and no occupied overhang. Despite broader header wording, it deliberately omits particle-number and spin-conservation checks. It is also a vector-like collection, not a mathematical set that automatically removes duplicates. Supply unique determinants in the desired electron-number sector yourself.

For fixed alpha/beta counts, the determinant dimension in $K$ spatial orbitals is $\binom K{N_\alpha}\binom K{N_\beta}$. We already counted 400 for three alpha and three beta electrons in six spatial orbitals. The full nitrogen count is enormous; a compact bit record does not remove the combinatorial number of such records.

Fixing $N_\alpha-N_\beta$ fixes a spin projection. It does not uniquely request a total-spin eigenvalue. Even a full determinant calculation in one projection sector can contain several total-spin sectors, and an incomplete selected determinant space may not be closed under total-spin operations. These interfaces take active alpha/beta counts, not a separate guaranteed singlet-targeting instruction. Spin interpretation must be checked rather than inferred only from equal counts.

Slater–Condon connectivity keeps the matrix structured

Define the CI matrix by $H_{IJ}=\langle D_I|\hat H|D_J\rangle$. The diagonal is the energy of determinant $I$; off-diagonal entries couple different occupations. The stationary coefficients satisfy

A one-electron operator can remove an electron from one occupied spin orbital and place it in one unoccupied spin orbital: a single excitation. A two-electron operator can change at most two occupations: a double excitation. Consequently, a Hamiltonian containing only one- and two-electron terms has zero matrix elements between determinants differing by more than two excitations. These are the central Slater–Condon rules. The diagonal and connected entries can be evaluated from the molecular integrals without integrating the full many-electron coordinate wavefunction afresh.

The nonzero values include fermionic signs. Creation and annihilation operators are algebraic operations for adding or removing occupation in a specified ordered list of spin orbitals. Passing an operator through an occupied mode changes sign. That bookkeeping expresses antisymmetry; omitting it changes interference and energies, not merely the appearance of the coefficient vector. The bundled MACIS determinant operations and excitation contribution code implement these connectivity and sign operations.

An explicitly illustrative two-determinant matrix is

The first determinant alone has energy −1.0 Hartree. Optimizing both amplitudes gives a lowest eigenvalue $-0.8-\sqrt{0.05}\approx-1.0236068$ Hartree. The lowering comes from coherent mixing through the off-diagonal coupling. This is not a molecular result or a repeat of the different toy matrix in Lesson 2.

For an exact Hermitian Hamiltonian and an orthonormal determinant subspace, the lowest projected eigenvalue is variational: enlarging a nested determinant space cannot raise its exact minimum. In practice matrix-element screening and incomplete eigensolver convergence introduce separate numerical approximations. A smaller selected space may also be replaced, rather than enlarged, during selection refinement, so successive reported energies need not decrease monotonically.

FCI, CASCI, and selected CI reduce different things

Full configuration interaction (FCI) includes every determinant with the prescribed particle counts in the finite orbital space. It is exact within that finite model and the specified sector, not exact continuum chemistry. It still carries basis, relativistic/ECP, frozen-core, and other modeling assumptions.

Complete active space CI (CASCI) includes every allowed determinant in a selected active orbital space while holding inactive orbitals filled and excluded orbitals outside the correlated treatment. The word complete refers to the active determinant space. Selecting six active orbitals and then solving all 400 determinants for three alpha/three beta electrons is different from selecting 400 determinants spread over 28 active orbitals.

Selected CI keeps only chosen determinants from a potentially much larger orbital space. This can target important configurations economically, but sacrifices completeness. ASCI, adaptive sampling configuration interaction, selects and refines that determinant list using estimates of missing amplitudes; its name does not mean that this native solver executes a quantum sampling experiment.

The actual multi-configuration interface is run(hamiltonian, n_active_alpha_electrons, n_active_beta_electrons), returning (energy, wavefunction). Those counts are active electrons, not automatically the molecule’s full fourteen. The Hamiltonian has an active one-body operator, two-body integrals, and a constant core energy. Native solvers diagonalize the active electronic matrix and add hamiltonian.get_core_energy() to the reported result. Do not add that constant a second time or assume it is only physical nuclear repulsion when inactive orbitals or ECPs are present.

The native complete-space solver and its real limits

The default implementation is multi_configuration_calculator/macis_cas. The adapter requires a restricted Hamiltonian and identical alpha/beta active orbital indices. Restricted means a common spatial-orbital representation; unequal active alpha/beta counts can still describe an open-shell sector. An unrestricted HF result is not accepted directly just because a common API produced it.

MACIS generates alpha occupation combinations, beta combinations, and their Cartesian product. The adapter constructs a SortedDoubleLoopHamiltonianGenerator from restricted one- and two-electron integrals. Its sorting exploits determinant connectivity when building matrix elements; it does not change the physical determinant space.

For one determinant the result is its diagonal matrix element with coefficient one. For a determinant dimension at or below iterative_solver_dimension_cutoff=2000, the adapter builds a sparse CSR Hamiltonian, converts it to a dense matrix, and calls LAPACK’s symmetric eigensolver. CSR, compressed sparse row storage, records nonzero values and their column positions grouped by row. Dense diagonalization stores all matrix entries and computes all eigenpairs even though this API returns only the lowest state.

For larger spaces the complete determinant list still exists, and the selected-CI diagonalization machinery builds a sparse Hamiltonian and runs Davidson. This is not a magic matrix-free FCI method that avoids enumerating all determinants. Memory includes the determinant list, sparse nonzero entries, vectors, and eigensolver workspace.

The true bitset dispatch in macis_base.hpp supports active-orbital representation widths up to 2048 spatial slots through successively larger template bitsets. Some older method comments say 128. Neither number is a feasible FCI problem-size promise: a representation may fit one determinant while the required number of determinants is impossible to store.

The tests exercise one-determinant cases and force the iterative path on a 441-determinant water/STO-3G space by lowering the cutoff. Those are regression cases, not demonstrations that the full 28-orbital nitrogen problem is cheap.

Davidson solves the coefficient problem in a small expanding subspace

Given a trial coefficient vector $c$, its energy estimate is the Rayleigh quotient, $E=c^{\dagger}Hc/(c^{\dagger}c)$. The residual

measures failure to satisfy the eigenvalue equation. A small energy change alone is not the same as a small residual.

Davidson maintains several orthonormal trial vectors, projects $H$ into their span, diagonalizes that small projected matrix, and constructs an improved full-space vector. It then approximately preconditions the residual using the Hamiltonian diagonal:

Small denominators are regularized; the new direction is orthogonalized against the current subspace and added if useful. This trades dense full-space diagonalization for repeated sparse matrix–vector products plus small dense eigenproblems.

The bundled implementation uses double reorthogonalization, a default minimum denominator magnitude of 1e-12, and an explicit residual-norm stopping test. It calls sparsexx for sparse matrix/vector work and lobpcgxx’s Rayleigh–Ritz helper for the projected problem. That helper use does not mean the exposed CI solver is a separate LOBPCG algorithm. BLAS/LAPACK provide underlying numerical algebra.

ci_residual_tolerance defaults to 1e-6; its declared range is zero to one, but zero should not be read as a useful guarantee of an exactly zero floating-point residual. ci_matel_tol defaults to machine epsilon, approximately 2.2e-16, and controls sparsification of CI Hamiltonian entries, not selection of determinant amplitudes.

A naming subtlety matters: the QDK setting max_solver_iterations=200 is mapped to MACIS ci_max_subspace, the maximum Davidson subspace dimension. In the current non-restarted Davidson implementation this also limits the sequence of expansions, but it is not an independently configurable general iteration budget. The limit is clamped to matrix dimension; failure to converge throws Davidson Did Not Converge!. Increasing it uses more full coefficient-sized vectors, not just a few scalar bookkeeping entries. Dense CAS/PMC cases do not use this residual stopping loop.

MACIS is bundled source with Lawrence Berkeley National Laboratory/University of California origins and Microsoft contributions recorded in its headers. It supplies determinant algorithms and solvers, not a remote service. Its tree contains MPI paths, but the enclosing QDK build configuration disables MACIS MPI in the configured integration, and selected incremental/RDM paths contain additional MPI restrictions. Do not infer that every visible parallel library or every build option participates in every solve.

A small complete CI handoff, without attempting full nitrogen FCI

The repository’s MC example uses H₂ at 0.740848 Å in STO-3G. Here is that separate small-system pattern with density outputs requested explicitly:

from qdk_chemistry.algorithms import create
from qdk_chemistry.data import Ansatz, Structure

h2_structure = Structure.from_xyz("""2
Separate small CI example
H 0.0 0.0 0.0
H 0.0 0.0 0.740848
""")
_, h2_reference = create("scf_solver", "qdk").run(h2_structure, 0, 1, "sto-3g")
h2_hamiltonian = create("hamiltonian_constructor", "qdk").run(h2_reference.get_orbitals())
ci = create(
    "multi_configuration_calculator", "macis_cas",
    calculate_one_rdm=True, calculate_two_rdm=True,
)
h2_energy, h2_ci = ci.run(h2_hamiltonian, 1, 1)
energy_from_density = Ansatz(h2_hamiltonian, h2_ci).calculate_energy()

These calls were not executed. In this finite two-spatial-orbital example there are four determinants in the fixed one-alpha/one-beta sector. The final line requests a consistency check using the returned densities, not another optimization. Its energy should agree with the CI result within the applicable numerical conventions and tolerances. We will select a tractable nitrogen active space in Lesson 10 instead of pretending this H₂ calculation solved stretched N₂.

Projected MC optimizes coefficients on a supplied determinant list

The implementation projected_multi_configuration_calculator/macis_pmc receives a Hamiltonian and an explicit list of Configuration objects. It builds and diagonalizes the Hamiltonian restricted to those determinants, using the same one-determinant/dense/iterative size distinction as CAS. The list is a free input: it may come from a scientific heuristic, a sampled support, or a deliberately constructed excitation set.

For the separate H₂ example, a two-configuration request would be:

from qdk_chemistry.data import Configuration

pmc = create("projected_multi_configuration_calculator", "macis_pmc")
selected = [Configuration.from_spin_half_string("20"),
            Configuration.from_spin_half_string("02")]
projected_energy, projected_state = pmc.run(h2_hamiltonian, selected)

The result optimizes both coefficients anew. In contrast, wavefunction.truncate(k) keeps the largest-magnitude amplitudes and renormalizes the retained vector; it does not rediagonalize the projected Hamiltonian. For a fixed valid subspace, coefficient reoptimization can do no worse than any normalized truncated candidate in that same subspace. The projected solver is therefore not merely a sophisticated spelling of truncation.

The current projected adapter rejects unrestricted Hamiltonians, insufficient configuration capacity and empty determinant lists. Its list validation is not a complete physical-sector or uniqueness validator. Provide unique determinants with consistent active indexing and particle counts. Duplicate records are not independent orthonormal basis functions, even if an input path accepts them. The PMC tests establish the projected energy/coefficient interface; they should not be read as permission to ignore those mathematical requirements.

ASCI estimates which missing determinants are worth adding

In macis_asci, the variational determinant list evolves. Suppose a current approximation has energy $E$ and amplitudes $c_I$. Rearranging the CI equation suggests a missing determinant $J$ might have amplitude approximately

where $\mathcal C$ is a selected set of important core determinants. Here core means search seeds, not frozen core electrons or core orbitals. The denominator measures the candidate diagonal energy relative to the current eigenvalue; near-degenerate candidates can become important. Signed contributions from several parents must be accumulated before ranking because they can reinforce or cancel.

The contribution kernels generate connected single and double excitations, include fermionic signs, compute diagonal estimates efficiently, and store numerator/denominator contributions. The search accumulates repeated candidates and keeps those with large ratio magnitudes. These ratios guide selection; they are not normalized final amplitudes. The chosen list is then rediagonalized to obtain the variational coefficients.

The native adapter starts from the canonical HF determinant. Growth increases the target determinant count; refinement repeatedly replaces/searches the list and rediagonalizes at approximately fixed size. If ntdets_max is strictly greater than the computed full-CI dimension, the adapter instead invokes the complete-space path. Equality is not that shortcut’s condition. The CAS/PMC dense-size cutoff is not the adaptive ASCI solver’s general dispatch rule.

The important default controls are:

Percentage selection can require far more than 100 core determinants when weight is diffuse. Thus ncdets_max=100 does not make the default search cheap in a strongly multiconfigurational state. just_singles=False includes doubles; setting it true changes which connected configurations are searched and can omit essential correlation. Search constraint levels partition work; they are not a total-spin setting or a declaration that only that excitation rank exists.

Growth uses a minimum factor of 1.01, backs off with growth_backoff_rate=0.5 if the requested size is not achieved, and recovers toward the requested factor using growth_recovery_rate=1.1. It can stop when no further growth is achieved. These growth semantics are more informative than treating ntdets_max as a guaranteed determinant count.

Refinement checks energy changes but retains a separately converged CI solve at each step. The current refinement routine can detect an oscillation and diagonalize the union of two determinant sets, extending the iteration budget by a bounded amount. That union can exceed the nominal trial target, so the target is not a hard peak-memory bound. If refinement does not meet its criterion, it throws rather than silently labeling the last list converged.

Warm-starting Davidson is enabled by default. Old coefficients are mapped by determinant identity into the new list and normalized if their retained vector norm meets min_warm_start_overlap=0.5. This is an L2 norm threshold, not a fifty-percent retained probability statement. Incremental Hamiltonian construction reuses sufficiently overlapping lists with min_patch_overlap=0.3; otherwise it rebuilds. grow_ci_residual_tolerance=0 means use the refinement tolerance, and positive values allow a looser growth solve. taper_grow_factor=0 disables the optional reduced final expansion factor.

The matrix-construction strategies sorted_double_loop, residue_arrays, and dynamic_bit_masking organize detection of connected determinants, not different physical Hamiltonians. The default is the sorted double loop. Residue arrays use excitation-deletion patterns and incur per-determinant electron-pair memory; the adapter falls back to the sorted method above sixty active electrons. Dynamic bit masking filters candidate connections using bit patterns; dynamic_bit_masking_num_masks=0 uses its internal default. Unknown strategy strings currently fall back to the sorted method, so inspect the chosen settings rather than infer a new algorithm from a typo.

pair_size_max=500000000 concerns stored search contributions and can imply substantial memory; it is not the total determinant target or a process-memory budget. More generous selection thresholds and larger determinant spaces need a deliberate resource plan. Conversely, a sparse state with large retained weight can still miss an energetically important collection of small amplitudes. Selection and refinement thresholds are not certified chemical-error bounds.

Several pt2_* controls are exposed from the bundled ASCI settings, but this wrapper’s return path does not execute an additional perturbative PT2 energy-correction stage. Its returned energy is the selected variational result plus the core constant. Also leave grow_with_rot=False in the present course workflow: the bundled optional route rotates integral arrays during growth, while this adapter does not return the accumulated orbital transformation with its original orbital object. It should not be presented as a provenance-safe drop-in orbital optimization. Use explicit orbital transformations and consistent Hamiltonian reconstruction instead.

Reduced density matrices retain what few-electron observables need

A full many-electron state contains vastly more information than is required for a one-electron observable such as density. Reducing, or tracing out, the unobserved particles means summing over their possibilities while retaining coherence in the coordinates or orbitals we still observe. It is not choosing the single most likely configuration.

For a normalized $N$-electron wavefunction, a coordinate-space 1-RDM can be expressed as

where $x$ includes position and spin and $X$ stands for the other electron coordinates. Its trace is $N$, not one, because it represents electron counts. A two-particle RDM analog retains two electron coordinates with a factor $N(N-1)$, counting ordered electron pairs. In a finite orbital representation, these become matrices/tensors that can be built by summing amplitude products and operator matrix elements between connected determinants.

Use a consistent convention. For the real molecular path here, the one-body density convention can be written $\gamma_{pq}=\langle a_q^{\dagger}a_p\rangle$. The QDK two-body convention is

The symbols $a^{\dagger}$ and $a$ are the creation/removal operators already used for fermionic signs. A diagonal one-body element is an occupation; a two-body entry describes a pair occupation or pair coherence. For one doubly occupied spatial orbital, the opposite-spin $\Gamma^{\alpha\beta}_{iiii}$ is one and its spin-traced ordered-pair value is two. In a normalized fixed-N full active-space density, $\sum_{pr}\Gamma(p,p,r,r)=N(N-1)$.

The spin-dependent accessors return alpha and beta 1-RDM blocks and the aaaa, aabb, bbbb 2-RDM blocks. Spin tracing adds alpha/beta one-body blocks; for two-body data it includes both opposite-spin orientations, with an index-pair transpose where needed. It is not always just doubling one supplied tensor. Public two-body arrays are flattened with the last index fastest: $pK^3+qK^2+rK+s$. These conventions matter as much as array shape.

For a restricted Hamiltonian, the actual Ansatz.calculate_energy contracts

Here $g$ is the stored chemist-ordered two-electron tensor, and the explicit pair permutation is the implementation’s indexing convention. For unrestricted storage, same-spin terms carry one-half factors and the combined opposite-spin contribution carries factor one. Follow the matching producer and consumer’s explicit index contractions: a legacy beta-first comment is not a universal storage convention. Lesson 12 traces the native alpha-alpha|beta-beta production and reconstruction paths. Reuse the provided contraction or verify the exact convention before writing a guessed einsum. This particular native energy-contraction implementation expects real RDM variants; complex coefficient storage elsewhere does not imply universal complex-property support.

There is a genuine interface-normalization caveat. CAS/ASCI’s shared wavefunction builder multiplies MACIS spin-dependent two-RDMs by two to convert the bundled half-weighted convention to the QDK convention. PMC’s separate eager RDM path stores the raw spin-traced MACIS two-RDM without the same conversion. Therefore eager PMC two-RDM data must not be assumed interchangeable with CAS/ASCI for the above energy formula; verify pair traces and energy consistency before consuming it. Our RDM example deliberately uses CAS. This conclusion comes from the builder, the projected adapter, and MACIS RDM contribution factors, not a newly run calculation.

Ask for expensive density and entropy data deliberately

The common MC flags default false: calculate_one_rdm, calculate_two_rdm, calculate_single_orbital_entropies, calculate_two_orbital_entropies, and calculate_mutual_information. CAS/ASCI’s builder supports the corresponding requested products. A one-RDM block has $K^2$ entries; a two-RDM block has $K^4$, and several spin blocks can be required. Computing them also involves determinant-pair connectivity and amplitude products. “The energy converged” does not make all of those properties free.

PMC has a different eager path: requesting either one- or two-RDM allocates and computes both spin-traced densities. It does not use the same shared entropy builder. The normalization warning above and its limited spin-resolved payload are further reasons to inspect capabilities rather than rely on a shared settings class alone.

A single spatial orbital has the four local possibilities from Lesson 8, so its entropy can be computed from spin-resolved 1/2-RDM information. Two-orbital entropy retains the joint subsystem; it generally needs more information than blindly reading a few pair-density entries. The MC entropy path evaluates it from the determinant expansion. The container can also derive mutual information from already available entropy data using $I_{ij}=s_i+s_j-s_{ij}$, with no extra one-half factor in this implementation. Natural logarithms are used. Request only the properties needed for the next scientific decision.

A Wavefunction container is a capability contract

StateVectorContainer stores determinant records and one coefficient per record, with real or complex coefficients supported as data. Its constructor checks matching coefficient/determinant counts; it does not generally normalize arbitrary user amplitudes automatically. Native MACIS solvers produce real coefficient vectors for their restricted real-integral path.

Single determinants can generate compatible occupation-based RDMs lazily. A general multi-determinant state does not automatically generate every missing RDM when an accessor is called; the required densities need to have been supplied or computed by the producer. Capability methods such as has_one_rdm_spin_dependent() and has_two_rdm_spin_traced() matter. When a localizer supplies a non-determinantal 1-RDM alongside an Aufbau carrier, the container avoids pretending that missing higher densities can always be manufactured from that carrier. This preserves the caveat from Lesson 8.

The occupation accessors have another concrete nuance: for stored nontrivial RDM data, get_active_orbital_occupations() diagonalizes each real spin 1-RDM and returns descending eigenvalues. It does not rotate the attached coefficients at that moment. These spectral occupations must not automatically be assigned to the current non-natural orbital columns as their diagonal populations. For a specific orbital basis use the appropriate 1-RDM diagonal, or explicitly transform to a consistent natural-orbital representation.

AmplitudeContainer, used for MP2/CC-type results, instead stores a reference and excitation-amplitude blocks. It is pure storage: it does not expand those amplitudes into CI determinants or generate RDMs. Determinant-coefficient, occupation, norm, overlap, and RDM-based operations are not generally available on that amplitude state. contains_determinant checks the reference support, not an expansion of the correlated exponential. An amplitude-type tag identifies Møller–Plesset or coupled-cluster meaning; it does not turn the amplitudes into normalized state-vector probabilities. We will use this distinction explicitly in Lesson 11.

For normalized, aligned determinant expansions, the overlap is $\langle\Psi|\Phi\rangle=\sum_I c_I^*d_I$, and pure-state fidelity is its squared magnitude. If inputs are not normalized, divide by the product of squared norms. The current state-vector overlap method requires the same determinant count, compatible particle count, and the same orbital object, then takes a coefficient inner product. It does not reconcile different determinant orderings/supports for you. Align the determinants explicitly before treating that number as physical overlap; equal vector lengths alone are insufficient.

Mathematically, truncating a normalized unique-determinant expansion to support of weight $w=\sum_{I\in\mathcal S}|c_I|^2$ and renormalizing gives fidelity $w$ with the original. That is a retained-state overlap, not an energy-error bound. The current truncate implementation creates a new state without carrying old correlated density payloads. Do not use the size-restricted overlap API on unequal supports as if it performed that derivation automatically.

Serialization preserves the defining records and orbital association, but object identity and physical equivalence remain distinct: ConfigurationSet equality uses orbital-pointer identity rather than a complete deep comparison. Full provenance and migration treatment comes later; for now retain the model, active ordering, coefficients, and requested property capabilities together.

Classical diagonalization after qubit mapping is still classical

The separate qubit_hamiltonian_solver interface takes a QubitOperator and returns (ground_energy, numpy_eigenvector), not a molecular Wavefunction. Its implementations are qdk_dense_matrix_solver and default qdk_sparse_matrix_solver.

The Python code materializes the operator matrix. The dense implementation calls a LAPACK-backed full symmetric eigensolver. The sparse implementation constructs a sparse matrix, converts it to sparsexx through the binding, and runs Davidson with tol=1e-8 and max_m=20. It is not a matrix-free operator-action simulation. For $q$ qubits the vector has $2^q$ entries, and the dense matrix has $4^q$ entries; sparsity helps matrix storage without eliminating vector-size growth.

Two limitations are essential. Both implementations currently take the real part of the matrix before solving. They are therefore appropriate only when discarding the imaginary part does not change the Hamiltonian; general complex Hermitian operators are not faithfully handled by this path. Also, these solvers do not expose an electron-number or spin-sector projection argument. They find the lowest state of the supplied qubit operator’s entire represented space, which need not be the target molecular particle sector unless that restriction has been supplied by some other valid construction.

The tests cover small real-matrix qubit operators and compare energies/vectors. They establish a useful small-system reference, not quantum execution or advantage. A mapped-operator eigenvector also needs its encoding and mode order before it can be interpreted as electron occupations.

A correlated approximation can guide the next reduction

We now have several distinct tools. CASCI gives a complete coefficient solution inside a finite chosen active space. ASCI trades completeness for adaptive determinant selection. PMC optimizes on externally supplied support. RDMs and entropy products summarize specific information from those states. Classical qubit diagonalization provides another small-model check after mapping, with its own matrix and sector limitations.

None requires that our production strategy first solve the impossible full nitrogen FCI problem merely to decide what to keep. A carefully chosen preliminary space, a bounded selected-CI calculation, or another justified approximate correlated density can supply diagnostic evidence. Its quality and omissions must be assessed rather than assumed exact. Lesson 10 turns that evidence into an explicit active-orbital and active-electron model, separating omission error from coefficient-solver error and keeping the original molecular question visible.