20Measuring energies and other quantities

An expectation value is not automatically an eigenvalue

We can now prepare a trial state and construct circuits associated with a Hamiltonian. Before building a full phase-estimation experiment, there is a more direct task: estimate the expectation value of an observable on that state.

An observable is represented by a Hermitian operator. Its eigenvalues are possible ideal measurement outcomes, while its expectation is their probability-weighted mean. For state $|\psi\rangle$, the energy expectation is

If $|\psi\rangle$ is an energy eigenstate, this equals its eigenvalue. Otherwise it averages different eigenenergies. A successful measurement of a sparse trial state therefore validates the prepared state’s energy expectation, not a claim that the full CASCI ground energy has been obtained.

Our continuing target remains neutral singlet N₂ at 1.85 Å. Its selected six-electron, six-spatial-orbital model is distinct from the original 28-orbital cc-pVDZ model and from other notebooks’ geometries. The ordinary mapped operator contains active terms; the core scalar stays outside it. Those modeling choices do not disappear when measurement begins.

The current registered expectation-estimator implementation is expectation_estimator/qdk. It uses explicit measurement circuits, a configurable circuit executor, and classical processing of counts. We inspect those implementations and tests here without running any simulation, measurement, or scientific workflow.

Why measuring occupation bits once is insufficient

A computational-basis measurement reports occupation bits in the Jordan–Wigner interpretation. For a qubit, outcome bit zero means Z eigenvalue +1 and bit one means −1. Consequently the number operator $n=(I-Z)/2$ has the expected empty/occupied outcomes. This interpretation must be transformed for parity/BK encodings; bit count is not universally electron count.

A general mapped Hamiltonian is

Its Pauli words can contain X and Y as well as Z. A Z-basis occupation measurement does not reveal the X/Y coherence needed for all of those expectations. Nor does one shot reveal an expectation value: it samples one random outcome according to the Born rule, whose probabilities are squared amplitudes in the chosen measurement basis.

A shot means preparing the same intended state again and performing the selected measurement experiment. Repeated shots approximate the corresponding probability distribution. They do not read out the entire quantum state, and the measurement generally disturbs the state rather than leaving an unknown trial state available for free reuse.

For a Pauli word, the eigenvalue is ±1. Measure each nonidentity position in its indicated basis and multiply the single-qubit ±1 outcomes. Equivalently, for the recorded bits belonging to the word, return $(-1)^{\text{number of ones}}$. Identity positions contribute +1 and need not be measured.

For X, a Hadamard before Z measurement implements the X basis. For Y, apply $S^\dagger$ then H before Z measurement. If that basis-change matrix is $B=HS^\dagger$, then $B^\dagger ZB=Y$. Reversing that order or using S instead of its adjoint changes the observable or its sign. The Qiskit fallback in the estimator implementation uses exactly sdg then h for Y.

The native Q# MeasurementBasis.qs instead calls external MResetX, MResetY, or MResetZ. Identity positions are reset without appending a result. These primitives and the runtime are external QDK contributions; the repository wrapper is what we have reviewed.

Bitstrings, masks, and omitted identity positions

QubitOperator Pauli labels place qubit 0 on the right. Measurement results contain only measured positions, compressed into a smaller register. The estimator reconstructs a mask for each term on those compressed positions; it does not directly interpret every returned bit position as the original qubit index.

For example, a group basis XIZ measures qubit 2 in X and qubit 0 in Z, leaving qubit 1 unmeasured. The result has two bits. The mask for IIZ uses the compressed bit corresponding to qubit 0; the mask for XII uses the one corresponding to qubit 2. A three-character Pauli label and a two-character result string are therefore consistent, not a missing-data error.

The code sorts measured qubit indices, assigns contiguous mask positions, parses binary strings or strings beginning with 0x as hexadecimal integers, and obtains each term’s sign from mask-and-outcome parity. It normalizes by the sum of the returned count frequencies. Unsupported formats and an empty total count raise. Its tests explicitly check compressed masks, binary/hex formats, X/Y basis construction, and zero-count failures.

A returned histogram is still only data about one measurement basis. It does not contain the original complex amplitudes. Relative phase is learned indirectly through appropriately chosen basis measurements, not extracted as an extra column next to every occupation string.

Grouping measurements: two meanings of “commuting”

If two Pauli words are qubit-wise commuting, then at every position their nonidentity letters agree. Identity can coexist with any letter. Such words share a basis made only from independent one-qubit basis changes. For example, IZ and XZ share the basis XZ and can be estimated from the same shots.

Global commutation is less restrictive. XX and YY commute because they anticommute at two positions and the two signs cancel. But one local basis cannot simultaneously be X and Y on either qubit. Their common eigenbasis requires an entangling transformation, such as a suitable Bell-basis circuit, rather than just independent H/S gates.

This distinction is decisive in the current software. The expectation estimator supports one local Pauli basis per measured group. It does not provide a general joint-diagonalization circuit for arbitrary globally commuting groups. Attaching a global-commuting partition therefore does not automatically produce a measurement saving.

The term_grouper family returns a new QubitOperator with partition metadata; it does not mutate the input or execute measurements. It rejects non-Pauli representations. Its default is commuting, which is useful for evolution but is not the universally safe default for this measurement estimator. Choose qubit_wise_commuting explicitly for a local-basis measurement workflow.

The actual bundled strategies are commuting, qubit_wise_commuting, identity, and vacuum_annihilating. The optional NetworkX plugin adds nx_commuting and nx_qubit_wise_commuting. These names denote different grouping methods, not six expectation-estimation algorithms.

How the graph-coloring methods group terms

The native commuting groupers build a noncommutation graph. Each vertex is a Pauli term, and an edge joins a pair that cannot share the requested kind of group. A graph coloring assigns different colors to connected vertices; terms with the same color form a valid group.

The native algorithm visits terms in input order and puts each into the first existing group with no conflicts, or opens a new group. This greedy first-fit rule is deterministic for a fixed ordered input, but it need not minimize the number of groups. It compares pairs, costing roughly quadratic term-count work and potentially quadratic graph storage, with Pauli-label comparisons adding a width factor.

The NetworkX variants call external greedy DSATUR coloring. DSATUR selects according to saturation degree, the number of distinct colors already present among neighbors. It can improve a grouping but remains a heuristic, not a proof of minimum colors, minimum variance, or minimum hardware runtime. The conditional plugin loader registers those variants only when its dependency is available; no installation or availability test is being performed here.

commuting and nx_commuting use global Pauli commutation. Their qubit-wise counterparts use the stricter local criterion. None of those graph-coloring decisions uses the prepared state’s covariance matrix or an optimal shot-allocation calculation. A group count is a useful structural metric, not a complete statistical cost estimate.

The identity grouper means the identity transformation of the grouping problem: it puts every term into its own singleton group. It does not select only Pauli-I terms, replace the Hamiltonian by identity, or remove small coefficients.

Vacuum-annihilating groups are for a different purpose

The vacuum_annihilating grouper has tolerance=1e-9. It groups terms by the positions they flip and by the parity of their Y count. Acting on the all-zero computational state, a term contributes

where $F$ is its X/Y support. Only terms reaching the same bitstring can cancel. Separating even/odd Y parity makes the grouped terms mutually commuting and separates real from imaginary vacuum amplitudes for real coefficients.

The implementation requires finite real coefficients and finite nonnegative tolerance. It first checks each bucket’s total cancellation, then cuts it at exactly zero running sums rather than approximate prefixes that could strand an uncancelled remainder. Diagonal I/Z terms form an unconstrained phase group. It retains every term and returns a partition; it does not subtract their expectation values from a general chemical trial state.

Its main use is preserving the all-zero workspace eigenstate in the CSWAP control construction, not measuring the nitrogen trial state, which has six active electrons and is not computational vacuum. XX and YY can be in one such group but are not qubit-wise commuting; this expectation estimator will not jointly measure them with local basis changes. The partition tests verify these distinctions, cancellation rules, invalid coefficients/tolerances, and input preservation.

If an expectation is exactly known for a specifically prepared state, it can in principle be added classically. An identity coefficient is the simplest state-independent example. A vacuum-annihilating sum has zero expectation on the vacuum, not on every state. Deleting it from an operator because of one reference state is not an operator identity, and changing the state, noise model, or observable can invalidate that shortcut.

Partition metadata and actual estimator support

FlatPartition stores groups of integer indices into a QubitOperator’s Pauli list. LayeredPartition stores group → layer → term indices for hierarchical scheduling. strategy is a descriptive string; it is not a proof that the terms have the named property. The QubitOperator constructor checks coverage of its terms exactly once, but it does not derive commutation from a label. Use nonempty valid groups; an empty group can otherwise become an invalid empty operator when materialized.

The current estimator consumes FlatPartition groups as measurement groups and verifies local-basis compatibility when constructing each basis. It does not whitelist strategy names: a globally labelled group that happens to be QWC can work, while a group containing XX and YY raises. It does not automatically recolor or split that invalid flat group. An unsupported partition subtype, including LayeredPartition, is ignored and the estimator falls back to measuring each term individually.

With no partition, it also measures each term individually. Attaching a partition is therefore a real control over circuit generation, not cosmetic metadata. A partition optimized for evolution need not be suitable for measurement.

Groupers preserve the operator coefficients, label ordering, encoding string, and fermion-mode-order field. However, the current commuting, identity, and NetworkX constructors omit an existing tapering specification when constructing their output, whereas the vacuum-annihilating grouper explicitly carries it. Keep the original sector/tapering provenance rather than assume every new object copies all metadata. The represented Pauli arrays may be unchanged while important workflow context is no longer attached.

What one estimator call does

The public contract is

energy_result, measurement_data = estimator.run(
    circuit,
    qubit_hamiltonian,
    total_shots,
    noise_model=None,
)

The circuit prepares one fixed state; the operator supplies the observable. Parameterized/unbound circuit sweeps are not implemented by this call. Despite the class name referring to energy, it can estimate other valid Hermitian Pauli sums, such as an appropriately encoded particle-number operator or a correlation observable. Their units follow their coefficients; a particle-number expectation is not measured in Hartree.

The only estimator setting is a nested circuit_executor, defaulting to AlgorithmRef("circuit_executor","qdk_sparse_state_simulator"). Other compatible executors can be selected, such as the native full-state simulator or optional Qiskit Aer backend. The estimator does not become an analytic exact expectation evaluator when the executor uses a full-state simulator: this path still obtains finite-shot counts. Backend seeds control reproducibility, not a guarantee of physical accuracy or statistical independence beyond the experiment’s assumptions.

First the estimator requires a Pauli-decomposition container and resolves measurement groups. Its own input checks do not provide a complete general Hermiticity validation: ensure the observable has finite real combined coefficients. In particular the identity-only path takes the real part of the coefficient sum, so passing an imaginary “observable” is not a supported way to estimate a non-Hermitian operator.

Next it removes identity-only groups from execution and adds their coefficient sum as an exact offset. If an identity term shares a nonidentity group, it can remain there, contributing expectation one and variance zero from the same data. The routine does not otherwise screen small coefficients or automatically detect all state-specific deterministic terms. If every group is identity-only, it returns the exact sum with zero reported variance and no circuit shots.

For each remaining group, it creates a new measurement circuit. If the supplied circuit has a Q# callable, the wrapper allocates a register of the observable’s width, applies that callable, and uses the native measurement-basis routine. It does not automatically supply extra caller-owned phase-gradient ancillas from metadata. A trial callable requiring such a suffix must be adapted with the proper ownership contract before this path can measure it.

Without a Q# callable, the estimator uses optional Qiskit conversion, checks that the circuit width equals the Pauli-label width, appends local basis rotations and a compact classical measurement register, and returns QASM. This fallback does require Qiskit, even when the eventual chosen executor is a QDK simulator. Existing measurement/classical-control structure and extra workspace need to satisfy that conversion/composition path; a rendered circuit or arbitrary QIR object is not an automatic compatibility guarantee.

total_shots is a budget request with integer rounding

Let $G$ be the number of nonidentity measurement groups. The code requires total_shots >= G and assigns

to each group. It does not distribute the remainder, weight allocation by coefficients, estimate variances adaptively, or optimize state-preparation cost. The requested executions therefore total $G\lfloor\text{total\_shots}/G\rfloor$, which can be smaller than the supplied number. This rule is independent of the executor’s random seed.

One circuit-executor call is made per group, with shots=n_g and noise=noise_model. Every shot includes fresh state preparation and the group’s basis measurement. Fewer groups can reduce separate circuit programs, but does not automatically reduce the chosen total shot budget or total uncertainty. Read measurement_data.shots_list for allocations and sum each returned count dictionary for usable observed shots.

CircuitExecutorData distinguishes clean bitstring counts, total requested/executed shots, and separate loss bitstrings. The expectation estimator currently consumes only clean bitstring_counts and does not propagate all executor metadata or loss records into its MeasurementData. Statistics normalize by the observed count sum, which can differ from the allocation after loss. That is effectively a conditional-survival estimate, with possible postselection bias; an empty clean count set raises. More shots do not automatically remove that bias.

Means, quantum variance, and variance of an estimator

For a Pauli outcome random variable $X\in\{-1,+1\}$, let its true expectation be $m$. Since $X^2=1$, its single-shot variance is $1-m^2$. For $n$ independent identically prepared shots, the sample mean $\widehat m$ has variance

The estimator computes $\widehat m$ from count-weighted parity signs and reports the plug-in estimate

This is a variance estimate of the mean, not the single-shot quantum variance. It is slightly downward-biased at finite sample size. For $n>1$, the corresponding unbiased sample-variance estimate of the mean uses $(1-\widehat m^2)/(n-1)$. The current implementation uses $n$, not $n-1$.

The intrinsic quantum energy variance is a different quantity,

An exact eigenstate has zero intrinsic energy variance. Yet separately measuring its noncommuting Pauli components can produce noisy sample means. For the ground eigenstate of the illustrative $X+Z$, each individual Pauli expectation is $-1/\sqrt2$, so each single-term measurement has variance one-half even though the state’s exact energy is sharp. Estimating the energy by separate Pauli measurements therefore still costs shots.

The current returned energy_variance is calculated as

That matches a diagonal, independent-term approximation to the variance of the energy estimator. It is not a computation of $\langle H^2\rangle-\langle H\rangle^2$, not a certified error bound, and not automatically a confidence interval.

Grouped shots introduce covariance the current result omits

When two terms use the same shots, their estimated means can be correlated. For a group, the correct variance of the weighted sample-mean estimator includes

Covariance measures whether two quantities fluctuate together or oppositely. It can increase or decrease the total depending on the state and coefficient signs. The current estimator does not calculate those cross terms, even though grouped count dictionaries contain enough joint outcome information to estimate them.

As a hand-worked counterexample, measure Z on both qubits of the Bell state $(|00\rangle+|11\rangle)/\sqrt2$, and estimate $H=Z_0+Z_1$. Suppose an illustrative 100-shot histogram has 50 00 and 50 11. Each Pauli mean is zero and its reported mean variance is 0.01, so the code’s diagonal sum is 0.02. But the per-shot energy score is always +2 or −2. Its plug-in mean variance is $4/100=0.04$. Ignoring the positive covariance understates the variance by a factor two in this example.

One way to process the raw data correctly is to assign each group outcome $b$ a score $Y_g(b)=\sum_{j\in g}h_jx_j(b)$, then compute its weighted first and second moments and divide the score variance by the number of shots. This is a mathematical post-processing prescription, not a claim that the current qdk estimator already performs it. An unbiased finite-sample version applies the appropriate sample correction. Independent groups’ estimator variances can then be summed under the independence assumptions.

The usual large-sample standard error is the square root of the estimator variance. A two-sided normal approximation such as mean ±1.96 standard errors additionally assumes adequate sampling and approximately normal mean errors; it is not guaranteed for small counts, rare events, drift, correlated shots, or systematic noise. Zero observed variance is especially deceptive: one shot or a short run with all identical outcomes does not prove the true failure probability is zero. The result’s summary printing a plus/minus square root does not remove these qualifications.

A complete toy estimate, including rounding and offsets

Consider the ungrouped illustrative operator

Request 101 shots. Identity is handled classically; the two remaining groups receive 50 shots each, so only 100 shots are allocated. Suppose X-basis outcomes are 40 plus and 10 minus, while Z-basis outcomes are 15 plus and 35 minus. Then

The code’s plug-in mean variances are $0.64/50=0.0128$ and $0.84/50=0.0168$. Because these two terms were measured in separate groups, the diagonal aggregation gives

with square root about 0.0854 in the toy energy units. If an independently retained core scalar were 0.2, total energy would be 1.02 with the same sampling variance, assuming the scalar is exact. All counts and coefficients here are illustrative, not generated by a simulation or taken from nitrogen.

This example separates four ideas: an exact identity contribution, a requested versus allocated shot count, a noisy active-operator expectation, and an independently added core scalar. It does not claim 0.0854 is a certified total chemical uncertainty.

In principle, with known per-group score variance $V_g$, equal-cost shots are allocated optimally roughly in proportion to $\sqrt{V_g}$. For independently measured individual terms this becomes $|h_j|\sqrt{1-m_j^2}$. Those variances may require pilot samples, and unequal circuit costs change the optimization. The current estimator does none of this automatically; there is no exposed adaptive-allocation setting to turn on. A graph-coloring heuristic alone cannot promise statistically optimal savings.

Returned data, units, and what survives serialization

The actual result classes are EnergyExpectationResult and MeasurementData. The first stores energy_expectation_value, energy_variance, expvals_each_term, and variances_each_term. The last two are lists of arrays organized by measured group, with entries in that group’s term order. They are not automatically flat arrays indexed like the original complete Hamiltonian.

Pauli expectations and their mean variances are dimensionless. The weighted energy expectation has the coefficients’ units, and its aggregated variance has squared units. The identity-only shortcut returns empty per-term arrays because no measurements were made. The summary’s count label should not be used to infer the original total Pauli count from the length of a group-list field.

MeasurementData pairs measured group Hamiltonians with bitstring-count dictionaries and allocated shots_list. In mixed identity/nonidentity runs, identity-only groups are omitted from that measurement list; the offset is already included in the numerical energy result, but is not a separately named MeasurementData field. Retain the original operator and core scalar if you want to reconstruct the full calculation later. The all-identity shortcut instead retains its resolved Hamiltonian groups but has empty counts/shot lists.

Both classes have JSON/HDF5 methods and tests. The numerical result stores its scalar/array statistics. MeasurementData serialization is narrower than full QubitOperator serialization: it stores Pauli labels and coefficients for each group, counts, and shot allocation, but not their full encoding/tapering/partition provenance. Its JSON path also uses raw coefficients.tolist() rather than the package’s complex-array encoding; a complex NumPy dtype can yield Python complex values that ordinary JSON dumping cannot serialize, even when imaginary parts are zero. HDF5 stores coefficient arrays directly. These are current format limits, not a reason to silently discard physical complex data. None of those file methods was called here.

The result-class tests and measurement-data tests check particular payloads and round trips. Such structural tests do not certify covariance-aware uncertainty, all dtype cases, or preserved chemical provenance.

Following the supplied workflow without overinterpreting its output

The quickstart selects the top two determinants from a CASCI wavefunction, then runs projected_multi_configuration_calculator on that support. This reoptimizes coefficients; it is not merely deleting small amplitudes and keeping the old energy. It prepares the returned sparse state, maps the matching active Hamiltonian, groups it by qubit_wise_commuting, and compares measured total energy with the sparse projected reference E_sparse.

The state-preparation notebook follows the same benchmarking logic and prints a square root of energy_variance. Its prose and plus/minus display must be read with the current covariance and finite-sample qualifications above. It is a separate source example, not newly executed nitrogen data. A classical reference already known for a small state is useful for validating state loading; the benchmark itself is not proof of quantum advantage over that classical reference calculation.

For the continuing nitrogen baseline, an unexecuted pattern is:

from qdk_chemistry.algorithms import create
from qdk_chemistry.data import AlgorithmRef

# trial_circuit prepares the intended state in the same blocked JW active basis.
measurable = create("term_grouper", "qubit_wise_commuting").run(qubit_hamiltonian)
estimator = create(
    "expectation_estimator", "qdk",
    circuit_executor=AlgorithmRef("circuit_executor", "qdk_sparse_state_simulator", seed=42),
)
result, data = estimator.run(trial_circuit, measurable, total_shots=10000)
total_energy = result.energy_expectation_value + active_hamiltonian.get_core_energy()
allocated_shots = sum(data.shots_list)
usable_shots = sum(sum(counts.values()) for counts in data.bitstring_counts if counts is not None)

The shot count is illustrative, not a recommended converged budget. result.energy_variance remains the implementation’s diagonal plug-in quantity. For grouped results, use the raw outcomes and the covariance formula before interpreting a standard error. The code measures only the supplied Pauli observable. SOS structured operators are rejected by this estimator; their separate shifted block-encoding spectral contract must not be treated as another Pauli measurement offset automatically.

What more shots can and cannot fix

Sampling uncertainty can decrease with repeated valid preparations, approximately as the inverse square root of shot count under the usual assumptions. It cannot repair an incorrect basis, a too-small active space, a DFT model used as though it were the intended interacting Hamiltonian, a wrong electron sector, lost coefficient signs, or a mismatched tapered state/operator. Those are different errors.

Noise changes the prepared/measured distribution. Survival postselection, drift, and systematic readout errors can bias the estimate even when a plug-in variance is tiny. The current estimator forwards an optional noise profile; it does not implement a general error-mitigation, confidence-calibration, or loss-correction procedure. Backend/provider behavior and richer noise treatment belong to the execution lesson.

Complete resource accounting includes state preparation for every group and every shot, basis changes, measurement/reset, classical processing, and any circuit conversion or synthesis needed by the executor. Reducing group count or reporting a small statistical variance is not the same as reducing physical-qubit overhead.

Finally, expectation estimation and phase estimation answer different questions. Pauli sampling estimates an average energy of the prepared state. Ideal phase estimation samples eigenvalue branches, with trial-state overlap and finite resolution determining their outcomes. We now have the operator, state-loading, circuit, and statistical vocabulary needed to derive the repository’s full phase-estimation methods without mistaking one kind of result for the other.