diff --git a/pyqpanda-algorithm/example/QAlgBase/testeg_QPE.py b/pyqpanda-algorithm/example/QAlgBase/testeg_QPE.py new file mode 100644 index 00000000..dea3688c --- /dev/null +++ b/pyqpanda-algorithm/example/QAlgBase/testeg_QPE.py @@ -0,0 +1,33 @@ +"""Minimal Quantum Phase Estimation example. + +The Z eigenstate |1> has eigenvalue -1 = exp(2*pi*i*0.5), so three precision +qubits should identify phase 0.5 exactly. +""" + +from pyqpanda3.core import QCircuit, X, Z + +from pyqpanda_alg import QPE + + +def prepare_one(qubits): + circuit = QCircuit() + circuit << X(qubits[0]) + return circuit + + +def z_unitary(qubits): + circuit = QCircuit() + circuit << Z(qubits[0]) + return circuit + + +if __name__ == "__main__": + result = QPE.run_qpe( + z_unitary, + target_width=1, + precision_bits=3, + prepare_eigenstate=prepare_one, + ) + print("phase:", result.phase) + print("maximum bitstring:", result.bitstring) + print("probabilities:", result.probability_dict()) diff --git a/pyqpanda-algorithm/pyqpanda_alg/QPE/README.md b/pyqpanda-algorithm/pyqpanda_alg/QPE/README.md new file mode 100644 index 00000000..f027f73a --- /dev/null +++ b/pyqpanda-algorithm/pyqpanda_alg/QPE/README.md @@ -0,0 +1,92 @@ +# Quantum Phase Estimation (QPE) + +This package adds a standalone Quantum Phase Estimation implementation to +pyqpanda-algorithm. + +QPE estimates the eigenphase phi of a unitary U on an eigenstate |psi>: + + U |psi> = exp(2*pi*i*phi) |psi>, 0 <= phi < 1. + +The implementation is split into two layers so algorithm behavior is useful +even when a PyQPanda3 runtime is not present: + +1. Pure reference helpers compute the exact ideal finite-register QPE + probability distribution, decode the maximum-likelihood grid point, and + compare phases on the unit circle. +2. The PyQPanda3 layer builds the standard controlled-power circuit and runs it + on CPUQVM. Imports are lazy, so importing the pure helpers does not require + PyQPanda3. + +## Public API + + from pyqpanda_alg import QPE + + probabilities = QPE.qpe_reference_probabilities( + phase=0.375, + precision_bits=4, + ) + estimate = QPE.decode_phase(probabilities, 4) + +For a runtime circuit, provide a factory returning the unitary QCircuit and, +optionally, a factory preparing its eigenstate: + + from pyqpanda_alg import QPE + from pyqpanda3.core import QCircuit, X, Z + + def prepare_one(qubits): + circuit = QCircuit() + circuit << X(qubits[0]) + return circuit + + def z_unitary(qubits): + circuit = QCircuit() + circuit << Z(qubits[0]) + return circuit + + result = QPE.run_qpe( + z_unitary, + target_width=1, + precision_bits=3, + prepare_eigenstate=prepare_one, + ) + print(result.phase) # 0.5 + print(result.bitstring) # counting-register maximum + +Z|1> = -|1> = exp(2*pi*i*0.5)|1>, so this example has an exactly +representable phase. + +## Conventions + +- Phases are fractions of one full turn in the half-open interval [0, 1). +- Counting qubit i controls U raised to 2**i. +- The inverse QFT and returned bit-string integer convention follow the same + ordering used by the existing pyqpanda_alg.QAE implementation. +- The circuit builder adds no measurement operations, so callers can compose + it into larger programs. +- target_width can exceed one; the unitary factory receives the complete target + register and may return any compatible QCircuit. +- precision_bits is capped at 24 in the pure helpers to avoid accidental + exponential allocation of a probability vector. + +## Reference distribution + +For N = 2**m and integer output y, QPE assigns + + P(y) = |(1/N) * sum_{k=0}^{N-1} + exp(2*pi*i*k*(phi-y/N))|**2. + +The implementation evaluates the corresponding Dirichlet-kernel expression, +including the removable singularity at exact grid points. The returned vector +is normalized with math.fsum. + +## Failure behavior + +Invalid precision, non-finite phases/probabilities, negative probability mass, +wrong distribution widths, empty runtime mappings, overlapping registers, and +invalid circuit factories fail with ValueError instead of silently producing a +phase estimate. + +## Scope + +This contribution is a general algorithm primitive. It does not change QAE, +Grover, Simon, SPSA, QKMeans, QARM, QSVR, QWalk, or QEC behavior. diff --git a/pyqpanda-algorithm/pyqpanda_alg/QPE/__init__.py b/pyqpanda-algorithm/pyqpanda_alg/QPE/__init__.py new file mode 100644 index 00000000..89f28b67 --- /dev/null +++ b/pyqpanda-algorithm/pyqpanda_alg/QPE/__init__.py @@ -0,0 +1,25 @@ +"""Quantum Phase Estimation and deterministic reference helpers.""" + +from .qpe import ( + QPEResult, + canonical_phase, + circular_phase_distance, + decode_phase, + decode_phase_mapping, + nearest_phase_grid_point, + qpe_circuit, + qpe_reference_probabilities, + run_qpe, +) + +__all__ = [ + "QPEResult", + "canonical_phase", + "circular_phase_distance", + "decode_phase", + "decode_phase_mapping", + "nearest_phase_grid_point", + "qpe_circuit", + "qpe_reference_probabilities", + "run_qpe", +] diff --git a/pyqpanda-algorithm/pyqpanda_alg/QPE/qpe.py b/pyqpanda-algorithm/pyqpanda_alg/QPE/qpe.py new file mode 100644 index 00000000..a315ba9c --- /dev/null +++ b/pyqpanda-algorithm/pyqpanda_alg/QPE/qpe.py @@ -0,0 +1,325 @@ +"""Quantum Phase Estimation for PyQPanda3. + +The module exposes two layers: + +* Pure reference helpers for the exact ideal QPE distribution and deterministic + phase decoding. These helpers do not import PyQPanda3. +* A lazy PyQPanda3 circuit builder/runner for arbitrary unitary circuits and + caller-supplied eigenstate preparation. + +Phase convention +---------------- +Eigenphases are fractions in [0, 1): U|psi> = exp(2*pi*i*phase)|psi>. +Counting-register powers follow the same ordering used by pyqpanda_alg.QAE: +counting qubit i controls U**(2**i), followed by inverse QFT over the counting +register. Runtime bit strings are decoded using the integer convention returned +by PyQPanda3's probability dictionary. +""" + +from __future__ import annotations + +from dataclasses import dataclass +import math +from collections.abc import Callable, Iterable, Mapping, Sequence +from typing import Any + + +CircuitFactory = Callable[[Sequence[Any]], Any] + + +def _validate_precision_bits(precision_bits: int) -> int: + if not isinstance(precision_bits, int) or isinstance(precision_bits, bool): + raise ValueError("precision_bits must be an integer") + if precision_bits <= 0: + raise ValueError("precision_bits must be positive") + if precision_bits > 24: + raise ValueError("precision_bits must be <= 24") + return precision_bits + + +def canonical_phase(phase: float) -> float: + """Return a finite phase reduced into the half-open interval [0, 1).""" + try: + value = float(phase) + except (TypeError, ValueError) as exc: + raise ValueError("phase must be a finite real number") from exc + if not math.isfinite(value): + raise ValueError("phase must be a finite real number") + value %= 1.0 + # Avoid returning 1.0 through floating-point edge cases. + return 0.0 if value == 1.0 else value + + +def circular_phase_distance(left: float, right: float) -> float: + """Return the shortest distance between two phases on the unit circle.""" + a = canonical_phase(left) + b = canonical_phase(right) + delta = abs(a - b) + return min(delta, 1.0 - delta) + + +def qpe_reference_probabilities( + phase: float, + precision_bits: int, +) -> tuple[float, ...]: + """Return the ideal QPE counting-register distribution. + + For N = 2**precision_bits and output integer y, the probability is + + |(1/N) sum_k exp(2*pi*i*k*(phase-y/N))|**2. + + The implementation evaluates the closed-form Dirichlet-kernel expression, + with an exact-grid branch that avoids the removable zero/zero singularity. + """ + bits = _validate_precision_bits(precision_bits) + phi = canonical_phase(phase) + size = 1 << bits + probabilities: list[float] = [] + + for outcome in range(size): + # math.remainder gives the numerically closest periodic displacement. + delta = math.remainder(phi - outcome / size, 1.0) + if abs(delta) <= 1e-15: + probability = 1.0 + else: + denominator = math.sin(math.pi * delta) + if abs(denominator) <= 1e-15: + probability = 1.0 + else: + numerator = math.sin(math.pi * size * delta) + probability = (numerator / (size * denominator)) ** 2 + probabilities.append(max(0.0, float(probability))) + + total = math.fsum(probabilities) + if not math.isfinite(total) or total <= 0.0: + raise ValueError("reference distribution is not normalizable") + return tuple(value / total for value in probabilities) + + +def _validate_probability_vector( + probabilities: Sequence[float], + precision_bits: int, +) -> tuple[float, ...]: + bits = _validate_precision_bits(precision_bits) + expected = 1 << bits + if len(probabilities) != expected: + raise ValueError( + f"expected {expected} probabilities for {bits} precision bits" + ) + values: list[float] = [] + for probability in probabilities: + try: + value = float(probability) + except (TypeError, ValueError) as exc: + raise ValueError("probabilities must be finite numeric values") from exc + if not math.isfinite(value): + raise ValueError("probabilities must be finite") + if value < 0.0: + raise ValueError("probabilities must be non-negative") + values.append(value) + total = math.fsum(values) + if total <= 0.0: + raise ValueError("probabilities must contain positive mass") + return tuple(value / total for value in values) + + +def decode_phase( + probabilities: Sequence[float], + precision_bits: int, +) -> float: + """Decode the maximum-likelihood grid phase from a QPE distribution.""" + values = _validate_probability_vector(probabilities, precision_bits) + outcome = max(range(len(values)), key=lambda index: values[index]) + return outcome / (1 << precision_bits) + + +def decode_phase_mapping( + probabilities: Mapping[str, int | float], + precision_bits: int, +) -> float: + """Decode a PyQPanda-style binary probability/count mapping. + + Keys are interpreted exactly as conventional binary integers, matching the + existing QAE runner in pyqpanda-algorithm. + """ + bits = _validate_precision_bits(precision_bits) + if not probabilities: + raise ValueError("probability mapping must not be empty") + + winner: str | None = None + winner_weight = -1.0 + for raw_key, raw_weight in probabilities.items(): + key = str(raw_key) + if len(key) != bits or any(ch not in "01" for ch in key): + raise ValueError( + f"every probability key must be a {bits}-bit binary string" + ) + try: + weight = float(raw_weight) + except (TypeError, ValueError) as exc: + raise ValueError("probability weights must be finite numbers") from exc + if not math.isfinite(weight) or weight < 0.0: + raise ValueError("probability weights must be finite and non-negative") + if weight > winner_weight or ( + weight == winner_weight and (winner is None or key < winner) + ): + winner = key + winner_weight = weight + + if winner is None or winner_weight <= 0.0: + raise ValueError("probability mapping must contain positive mass") + return int(winner, 2) / (1 << bits) + + +def nearest_phase_grid_point(phase: float, precision_bits: int) -> float: + """Return the nearest representable QPE phase grid point.""" + bits = _validate_precision_bits(precision_bits) + phi = canonical_phase(phase) + size = 1 << bits + return (int(math.floor(phi * size + 0.5)) % size) / size + + +def _resolve_registers( + precision_bits: int, + target_width: int, + counting_qubits: Iterable[Any] | None, + target_qubits: Iterable[Any] | None, +) -> tuple[list[Any], list[Any]]: + bits = _validate_precision_bits(precision_bits) + if not isinstance(target_width, int) or isinstance(target_width, bool): + raise ValueError("target_width must be an integer") + if target_width <= 0: + raise ValueError("target_width must be positive") + + counting = ( + list(range(bits)) + if counting_qubits is None + else list(counting_qubits) + ) + target = ( + list(range(bits, bits + target_width)) + if target_qubits is None + else list(target_qubits) + ) + if len(counting) != bits: + raise ValueError(f"expected exactly {bits} counting qubits") + if len(target) != target_width: + raise ValueError(f"expected exactly {target_width} target qubits") + labels = [str(qubit) for qubit in counting + target] + if len(set(labels)) != len(labels): + raise ValueError("counting and target qubits must be distinct") + return counting, target + + +def qpe_circuit( + unitary: CircuitFactory, + *, + target_width: int, + precision_bits: int, + prepare_eigenstate: CircuitFactory | None = None, + counting_qubits: Iterable[Any] | None = None, + target_qubits: Iterable[Any] | None = None, +): + """Build a standard Quantum Phase Estimation circuit. + + unitary(target_qubits) must return a QCircuit implementing U. + prepare_eigenstate(target_qubits), when provided, must return a circuit + preparing an eigenstate of U from |0...0>. + + The circuit does not add measurements, allowing callers to compose it into + larger programs or inspect exact probabilities. + """ + if not callable(unitary): + raise ValueError("unitary must be callable") + if prepare_eigenstate is not None and not callable(prepare_eigenstate): + raise ValueError("prepare_eigenstate must be callable") + + counting, target = _resolve_registers( + precision_bits, + target_width, + counting_qubits, + target_qubits, + ) + + from pyqpanda3.core import H, QCircuit + from ..plugin import QFT + + circuit = QCircuit() + if prepare_eigenstate is not None: + circuit << prepare_eigenstate(target) + + for qubit in counting: + circuit << H(qubit) + + # Match the power ordering already used by pyqpanda_alg.QAE. + for power_index, control in enumerate(counting): + repetitions = 1 << power_index + for _ in range(repetitions): + operation = unitary(target) + if operation is None: + raise ValueError("unitary factory returned None") + circuit << operation.control([control]) + + circuit << QFT(counting).dagger() + return circuit + + +@dataclass(frozen=True) +class QPEResult: + """Result of a CPUQVM QPE execution.""" + + phase: float + bitstring: str + probabilities: tuple[tuple[str, float], ...] + + def probability_dict(self) -> dict[str, float]: + return dict(self.probabilities) + + +def run_qpe( + unitary: CircuitFactory, + *, + target_width: int, + precision_bits: int = 6, + prepare_eigenstate: CircuitFactory | None = None, + shots: int = 1024, +) -> QPEResult: + """Execute QPE on CPUQVM and return the maximum-likelihood eigenphase. + + The arbitrary unitary is supplied as a circuit factory so the same public + API supports one- and multi-qubit eigenstates without matrix conversion. + """ + bits = _validate_precision_bits(precision_bits) + if not isinstance(shots, int) or isinstance(shots, bool) or shots <= 0: + raise ValueError("shots must be a positive integer") + + counting = list(range(bits)) + target = list(range(bits, bits + target_width)) + + from pyqpanda3.core import CPUQVM, QProg + + program = QProg() + program << qpe_circuit( + unitary, + target_width=target_width, + precision_bits=bits, + prepare_eigenstate=prepare_eigenstate, + counting_qubits=counting, + target_qubits=target, + ) + + machine = CPUQVM() + machine.run(program, shots) + raw = machine.result().get_prob_dict(counting) + normalized = {str(key): float(value) for key, value in raw.items()} + phase = decode_phase_mapping(normalized, bits) + winner = max( + normalized, + key=lambda key: (normalized[key], -int(key, 2)), + ) + ordered = tuple(sorted(normalized.items(), key=lambda item: item[0])) + return QPEResult( + phase=phase, + bitstring=winner, + probabilities=ordered, + ) diff --git a/test/QPE/Test_qpe.py b/test/QPE/Test_qpe.py new file mode 100644 index 00000000..5dceefbd --- /dev/null +++ b/test/QPE/Test_qpe.py @@ -0,0 +1,91 @@ +import math + +import pytest + +from pyqpanda_alg import QPE + + +def test_exact_grid_reference_is_one_hot_and_decodes_exactly(): + for bits in range(1, 7): + size = 1 << bits + for outcome in range(size): + phase = outcome / size + probabilities = QPE.qpe_reference_probabilities(phase, bits) + assert sum(probabilities) == pytest.approx(1.0) + assert probabilities[outcome] == pytest.approx(1.0, abs=1e-12) + assert QPE.decode_phase(probabilities, bits) == phase + + +def test_non_grid_phase_decodes_to_nearest_high_probability_bin(): + for bits in (3, 5, 7): + resolution = 1.0 / (1 << bits) + for phase in (0.071, 0.247, 0.499, 0.731, 0.997): + probabilities = QPE.qpe_reference_probabilities(phase, bits) + estimate = QPE.decode_phase(probabilities, bits) + assert QPE.circular_phase_distance(estimate, phase) <= resolution + assert QPE.circular_phase_distance( + QPE.nearest_phase_grid_point(phase, bits), + phase, + ) <= resolution / 2 + 1e-15 + + +def test_phase_periodicity_and_circular_distance(): + assert QPE.canonical_phase(1.25) == pytest.approx(0.25) + assert QPE.canonical_phase(-0.25) == pytest.approx(0.75) + assert QPE.qpe_reference_probabilities(0.125, 4) == pytest.approx( + QPE.qpe_reference_probabilities(1.125, 4) + ) + assert QPE.circular_phase_distance(0.99, 0.01) == pytest.approx(0.02) + + +def test_mapping_decoder_uses_binary_integer_convention(): + mapping = { + "000": 0.0, + "001": 0.1, + "100": 0.8, + "111": 0.1, + } + assert QPE.decode_phase_mapping(mapping, 3) == pytest.approx(0.5) + + +def test_validation_rejects_invalid_inputs(): + for bad in (0, -1, 25, 1.5, True): + with pytest.raises(ValueError): + QPE.qpe_reference_probabilities(0.5, bad) + + for bad_phase in (math.nan, math.inf, -math.inf, "not-a-phase"): + with pytest.raises(ValueError): + QPE.canonical_phase(bad_phase) + + with pytest.raises(ValueError, match="expected 8"): + QPE.decode_phase([1.0, 0.0], 3) + with pytest.raises(ValueError, match="non-negative"): + QPE.decode_phase([1.0, -1.0, 0.0, 0.0], 2) + with pytest.raises(ValueError, match="positive mass"): + QPE.decode_phase([0.0] * 4, 2) + with pytest.raises(ValueError, match="3-bit"): + QPE.decode_phase_mapping({"10": 1.0}, 3) + + +def test_cpuqvm_z_eigenphase_when_pyqpanda3_is_available(): + pytest.importorskip("pyqpanda3") + from pyqpanda3.core import QCircuit, X, Z + + def prepare_one(qubits): + circuit = QCircuit() + circuit << X(qubits[0]) + return circuit + + def z_unitary(qubits): + circuit = QCircuit() + circuit << Z(qubits[0]) + return circuit + + result = QPE.run_qpe( + z_unitary, + target_width=1, + precision_bits=3, + prepare_eigenstate=prepare_one, + ) + assert QPE.circular_phase_distance(result.phase, 0.5) < 1e-9 + assert result.bitstring == "100"