From 80b30f16d4b1caaf1cbc8165235293626f8e72a5 Mon Sep 17 00:00:00 2001 From: Julien Gacon Date: Thu, 24 Oct 2024 13:28:21 +0200 Subject: [PATCH 1/6] integer comparator & exact reciprocal --- pyproject.toml | 1 + qiskit/circuit/library/__init__.py | 2 + qiskit/circuit/library/arithmetic/__init__.py | 2 +- .../library/arithmetic/exact_reciprocal.py | 78 +++++++++---- .../library/arithmetic/integer_comparator.py | 109 +++++------------- qiskit/synthesis/arithmetic/__init__.py | 15 +++ .../arithmetic/comparators/__init__.py | 15 +++ .../comparators/integer_comparator.py | 101 ++++++++++++++++ .../passes/synthesis/hls_plugins.py | 35 ++++++ .../library/test_integer_comparator.py | 29 ++++- 10 files changed, 280 insertions(+), 107 deletions(-) create mode 100644 qiskit/synthesis/arithmetic/__init__.py create mode 100644 qiskit/synthesis/arithmetic/comparators/__init__.py create mode 100644 qiskit/synthesis/arithmetic/comparators/integer_comparator.py diff --git a/pyproject.toml b/pyproject.toml index 70a6e7c9d0cf..0f8156c46d58 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -100,6 +100,7 @@ sk = "qiskit.transpiler.passes.synthesis.solovay_kitaev_synthesis:SolovayKitaevS "qft.line" = "qiskit.transpiler.passes.synthesis.hls_plugins:QFTSynthesisLine" "qft.default" = "qiskit.transpiler.passes.synthesis.hls_plugins:QFTSynthesisFull" "permutation.token_swapper" = "qiskit.transpiler.passes.synthesis.hls_plugins:TokenSwapperSynthesisPermutation" +"IntegerComparator.default" = "qiskit.transpiler.passes.synthesis.hls_plugins:IntegerComparatorSynthesisDefault" [project.entry-points."qiskit.transpiler.init"] default = "qiskit.transpiler.preset_passmanagers.builtin_plugins:DefaultInitPassManager" diff --git a/qiskit/circuit/library/__init__.py b/qiskit/circuit/library/__init__.py index a4203f63b738..5d5363cd034f 100644 --- a/qiskit/circuit/library/__init__.py +++ b/qiskit/circuit/library/__init__.py @@ -288,6 +288,7 @@ :template: autosummary/class_no_inherited_members.rst IntegerComparator + IntegerComparatorGate Functions on binary variables ----------------------------- @@ -535,6 +536,7 @@ PiecewisePolynomialPauliRotations, PolynomialPauliRotations, IntegerComparator, + IntegerComparatorGate, WeightedAdder, QuadraticForm, LinearAmplitudeFunction, diff --git a/qiskit/circuit/library/arithmetic/__init__.py b/qiskit/circuit/library/arithmetic/__init__.py index 71c6cb7f5df5..0403ecf06d8d 100644 --- a/qiskit/circuit/library/arithmetic/__init__.py +++ b/qiskit/circuit/library/arithmetic/__init__.py @@ -13,7 +13,7 @@ """The arithmetic circuit library.""" from .functional_pauli_rotations import FunctionalPauliRotations -from .integer_comparator import IntegerComparator +from .integer_comparator import IntegerComparator, IntegerComparatorGate from .linear_pauli_rotations import LinearPauliRotations from .piecewise_linear_pauli_rotations import PiecewiseLinearPauliRotations from .piecewise_polynomial_pauli_rotations import PiecewisePolynomialPauliRotations diff --git a/qiskit/circuit/library/arithmetic/exact_reciprocal.py b/qiskit/circuit/library/arithmetic/exact_reciprocal.py index 5b53cf6e47b5..726816538e54 100644 --- a/qiskit/circuit/library/arithmetic/exact_reciprocal.py +++ b/qiskit/circuit/library/arithmetic/exact_reciprocal.py @@ -13,7 +13,7 @@ from math import isclose import numpy as np -from qiskit.circuit import QuantumCircuit, QuantumRegister +from qiskit.circuit import QuantumCircuit, QuantumRegister, Gate from qiskit.circuit.library.generalized_gates import UCRYGate @@ -46,43 +46,79 @@ def __init__( """ qr_state = QuantumRegister(num_state_qubits, "state") qr_flag = QuantumRegister(1, "flag") - circuit = QuantumCircuit(qr_state, qr_flag, name=name) + super().__init__(qr_state, qr_flag, name=name) + + reciprocal = ExactReciprocalGate(num_state_qubits, scaling, neg_vals, label=name) + self.append(reciprocal, self.qubits) + + +class ExactReciprocalGate(Gate): + r"""Exact reciprocal + + .. math:: + + |x\rangle |0\rangle \mapsto \cos(1/x)|x\rangle|0\rangle + \sin(1/x)|x\rangle |1\rangle + """ + + def __init__( + self, num_state_qubits: int, scaling: float, neg_vals: bool = False, label: str = "1/x" + ) -> None: + r""" + Args: + num_state_qubits: The number of qubits representing the value to invert. + scaling: Scaling factor :math:`s` of the reciprocal function, i.e. to compute + :math:`s / x`. + neg_vals: Whether :math:`x` might represent negative values. In this case the first + qubit is the sign, with :math:`|1\rangle` for negative and :math:`|0\rangle` for + positive. For the negative case it is assumed that the remaining string represents + :math:`1 - x`. This is because :math:`e^{-2 \pi i x} = e^{2 \pi i (1 - x)}` for + :math:`x \in [0,1)`. + name: The name of the object. + + .. note:: + + It is assumed that the binary string :math:`x` represents a number < 1. + """ + super().__init__("ExactReciprocal", num_state_qubits + 1, [], label=label) + + self.scaling = scaling + self.neg_vals = neg_vals + + def _define(self): + num_state_qubits = self.num_qubits - 1 + qr_state = QuantumRegister(num_state_qubits, "state") + qr_flag = QuantumRegister(1, "flag") + circuit = QuantumCircuit(qr_state, qr_flag) angles = [0.0] - nl = 2 ** (num_state_qubits - 1) if neg_vals else 2**num_state_qubits + nl = 2 ** (num_state_qubits - 1) if self.neg_vals else 2**num_state_qubits # Angles to rotate by scaling / x, where x = i / nl for i in range(1, nl): - if isclose(scaling * nl / i, 1, abs_tol=1e-5): + if isclose(self.scaling * nl / i, 1, abs_tol=1e-5): angles.append(np.pi) - elif scaling * nl / i < 1: - angles.append(2 * np.arcsin(scaling * nl / i)) + elif self.scaling * nl / i < 1: + angles.append(2 * np.arcsin(self.scaling * nl / i)) else: angles.append(0.0) - circuit.compose( - UCRYGate(angles), [qr_flag[0]] + qr_state[: len(qr_state) - neg_vals], inplace=True - ) + circuit.append(UCRYGate(angles), [qr_flag[0]] + qr_state[: len(qr_state) - self.neg_vals]) - if neg_vals: - circuit.compose( + if self.neg_vals: + circuit.append( UCRYGate([-theta for theta in angles]).control(), [qr_state[-1]] + [qr_flag[0]] + qr_state[:-1], - inplace=True, ) angles_neg = [0.0] for i in range(1, nl): - if isclose(scaling * (-1) / (1 - i / nl), -1, abs_tol=1e-5): + if isclose(self.scaling * (-1) / (1 - i / nl), -1, abs_tol=1e-5): angles_neg.append(-np.pi) - elif np.abs(scaling * (-1) / (1 - i / nl)) < 1: - angles_neg.append(2 * np.arcsin(scaling * (-1) / (1 - i / nl))) + elif np.abs(self.scaling * (-1) / (1 - i / nl)) < 1: + angles_neg.append(2 * np.arcsin(self.scaling * (-1) / (1 - i / nl))) else: angles_neg.append(0.0) - circuit.compose( - UCRYGate(angles_neg).control(), - [qr_state[-1]] + [qr_flag[0]] + qr_state[:-1], - inplace=True, + circuit.append( + UCRYGate(angles_neg).control(), [qr_state[-1]] + [qr_flag[0]] + qr_state[:-1] ) - super().__init__(*circuit.qregs, name=name) - self.compose(circuit.to_gate(), qubits=self.qubits, inplace=True) + self.definition = circuit diff --git a/qiskit/circuit/library/arithmetic/integer_comparator.py b/qiskit/circuit/library/arithmetic/integer_comparator.py index 377af5f222e0..e1f3b3e42fca 100644 --- a/qiskit/circuit/library/arithmetic/integer_comparator.py +++ b/qiskit/circuit/library/arithmetic/integer_comparator.py @@ -14,11 +14,10 @@ """Integer Comparator.""" from __future__ import annotations -import math -from qiskit.circuit import QuantumCircuit, QuantumRegister, AncillaRegister +from qiskit.circuit import QuantumRegister, AncillaRegister, Gate from qiskit.circuit.exceptions import CircuitError -from ..boolean_logic import OR +from qiskit.synthesis.arithmetic.comparators import synth_integer_comparator_2s from ..blueprintcircuit import BlueprintCircuit @@ -134,19 +133,6 @@ def num_state_qubits(self, num_state_qubits: int | None) -> None: qr_ancilla = AncillaRegister(num_ancillas) self.add_register(qr_ancilla) - def _get_twos_complement(self) -> list[int]: - """Returns the 2's complement of ``self.value`` as array. - - Returns: - The 2's complement of ``self.value``. - """ - twos_complement = pow(2, self.num_state_qubits) - math.ceil(self.value) - twos_complement = f"{twos_complement:b}".rjust(self.num_state_qubits, "0") - twos_complement = [ - 1 if twos_complement[i] == "1" else 0 for i in reversed(range(len(twos_complement))) - ] - return twos_complement - def _check_configuration(self, raise_on_failure: bool = True) -> bool: """Check if the current configuration is valid.""" valid = True @@ -176,68 +162,31 @@ def _build(self) -> None: super()._build() - qr_state = self.qubits[: self.num_state_qubits] - q_compare = self.qubits[self.num_state_qubits] - qr_ancilla = self.qubits[self.num_state_qubits + 1 :] - - circuit = QuantumCircuit(*self.qregs, name=self.name) - - if self.value <= 0: # condition always satisfied for non-positive values - if self._geq: # otherwise the condition is never satisfied - circuit.x(q_compare) - # condition never satisfied for values larger than or equal to 2^n - elif self.value < pow(2, self.num_state_qubits): - - if self.num_state_qubits > 1: - twos = self._get_twos_complement() - for i in range(self.num_state_qubits): - if i == 0: - if twos[i] == 1: - circuit.cx(qr_state[i], qr_ancilla[i]) - elif i < self.num_state_qubits - 1: - if twos[i] == 1: - circuit.compose( - OR(2), [qr_state[i], qr_ancilla[i - 1], qr_ancilla[i]], inplace=True - ) - else: - circuit.ccx(qr_state[i], qr_ancilla[i - 1], qr_ancilla[i]) - else: - if twos[i] == 1: - # OR needs the result argument as qubit not register, thus - # access the index [0] - circuit.compose( - OR(2), [qr_state[i], qr_ancilla[i - 1], q_compare], inplace=True - ) - else: - circuit.ccx(qr_state[i], qr_ancilla[i - 1], q_compare) - - # flip result bit if geq flag is false - if not self._geq: - circuit.x(q_compare) - - # uncompute ancillas state - for i in reversed(range(self.num_state_qubits - 1)): - if i == 0: - if twos[i] == 1: - circuit.cx(qr_state[i], qr_ancilla[i]) - else: - if twos[i] == 1: - circuit.compose( - OR(2), [qr_state[i], qr_ancilla[i - 1], qr_ancilla[i]], inplace=True - ) - else: - circuit.ccx(qr_state[i], qr_ancilla[i - 1], qr_ancilla[i]) - else: - - # num_state_qubits == 1 and value == 1: - circuit.cx(qr_state[0], q_compare) - - # flip result bit if geq flag is false - if not self._geq: - circuit.x(q_compare) - - else: - if not self._geq: # otherwise the condition is never satisfied - circuit.x(q_compare) - + circuit = synth_integer_comparator_2s(self.num_state_qubits, self.value, self.geq) self.append(circuit.to_gate(), self.qubits) + + +class IntegerComparatorGate(Gate): + r"""Perform a :math:`\geq` (or :math:`<`) on a qubit register against a classical integer. + + This operator compares basis states :math:`|i\rangle_n` against a classically given integer + :math:`L` of fixed value and flips a target qubit if :math:`i \geq L` + (or :math:`<` depending on the parameter ``geq``): + + .. math:: + + |i\rangle_n |0\rangle \mapsto |i\rangle_n |i \geq L\rangle + + """ + + def __init__( + self, num_state_qubits: int, value: int, geq: bool = True, label: str | None = None + ): + r""" + Args: + num_state_qubits: The number of qubits in the registers. + geq: If ``True`` compute :math:`i \geq L`, otherwise compute :math:`i < L`. + """ + super().__init__("IntegerComparator", num_state_qubits + 1, [], label=label) + self.value = value + self.geq = geq diff --git a/qiskit/synthesis/arithmetic/__init__.py b/qiskit/synthesis/arithmetic/__init__.py new file mode 100644 index 000000000000..89013c167c03 --- /dev/null +++ b/qiskit/synthesis/arithmetic/__init__.py @@ -0,0 +1,15 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2017, 2024. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""Synthesis for arithmetic circuits.""" + +from .comparators import synth_integer_comparator_2s diff --git a/qiskit/synthesis/arithmetic/comparators/__init__.py b/qiskit/synthesis/arithmetic/comparators/__init__.py new file mode 100644 index 000000000000..b492d5563ee9 --- /dev/null +++ b/qiskit/synthesis/arithmetic/comparators/__init__.py @@ -0,0 +1,15 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2017, 2024. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""Comparator synthesis algorithms.""" + +from .integer_comparator import synth_integer_comparator_2s diff --git a/qiskit/synthesis/arithmetic/comparators/integer_comparator.py b/qiskit/synthesis/arithmetic/comparators/integer_comparator.py new file mode 100644 index 000000000000..7c5816fe4320 --- /dev/null +++ b/qiskit/synthesis/arithmetic/comparators/integer_comparator.py @@ -0,0 +1,101 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2017, 2024. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""Integer comparator based on 2s complement.""" + +import math +from qiskit.circuit import QuantumCircuit +from qiskit.circuit.library import OR + + +def synth_integer_comparator_2s( + num_state_qubits: int, value: int, geq: bool = True +) -> QuantumCircuit: + """Implement an integer comparison based on 2s complement.""" + circuit = QuantumCircuit(2 * num_state_qubits) + qr_state = circuit.qubits[:num_state_qubits] + q_compare = circuit.qubits[num_state_qubits] + qr_ancilla = circuit.qubits[num_state_qubits + 1 :] + + if value <= 0: # condition always satisfied for non-positive values + if geq: # otherwise the condition is never satisfied + circuit.x(q_compare) + # condition never satisfied for values larger than or equal to 2^n + elif value < pow(2, num_state_qubits): + + if num_state_qubits > 1: + twos = _get_twos_complement(num_state_qubits, value) + for i in range(num_state_qubits): + if i == 0: + if twos[i] == 1: + circuit.cx(qr_state[i], qr_ancilla[i]) + elif i < num_state_qubits - 1: + if twos[i] == 1: + circuit.compose( + OR(2), [qr_state[i], qr_ancilla[i - 1], qr_ancilla[i]], inplace=True + ) + else: + circuit.ccx(qr_state[i], qr_ancilla[i - 1], qr_ancilla[i]) + else: + if twos[i] == 1: + # OR needs the result argument as qubit not register, thus + # access the index [0] + circuit.compose( + OR(2), [qr_state[i], qr_ancilla[i - 1], q_compare], inplace=True + ) + else: + circuit.ccx(qr_state[i], qr_ancilla[i - 1], q_compare) + + # flip result bit if geq flag is false + if not geq: + circuit.x(q_compare) + + # uncompute ancillas state + for i in reversed(range(num_state_qubits - 1)): + if i == 0: + if twos[i] == 1: + circuit.cx(qr_state[i], qr_ancilla[i]) + else: + if twos[i] == 1: + circuit.compose( + OR(2), [qr_state[i], qr_ancilla[i - 1], qr_ancilla[i]], inplace=True + ) + else: + circuit.ccx(qr_state[i], qr_ancilla[i - 1], qr_ancilla[i]) + else: + + # num_state_qubits == 1 and value == 1: + circuit.cx(qr_state[0], q_compare) + + # flip result bit if geq flag is false + if not geq: + circuit.x(q_compare) + + else: + if not geq: # otherwise the condition is never satisfied + circuit.x(q_compare) + + return circuit + + +def _get_twos_complement(num_bits: int, value: int) -> list[int]: + """Returns the 2's complement of ``self.value`` as array. + + Returns: + The 2's complement of ``self.value``. + """ + twos_complement = pow(2, num_bits) - math.ceil(value) + twos_complement = f"{twos_complement:b}".rjust(num_bits, "0") + twos_complement = [ + 1 if twos_complement[i] == "1" else 0 for i in reversed(range(len(twos_complement))) + ] + return twos_complement diff --git a/qiskit/transpiler/passes/synthesis/hls_plugins.py b/qiskit/transpiler/passes/synthesis/hls_plugins.py index fa827a20ed48..fb3000907d0f 100644 --- a/qiskit/transpiler/passes/synthesis/hls_plugins.py +++ b/qiskit/transpiler/passes/synthesis/hls_plugins.py @@ -249,16 +249,29 @@ MCMTSynthesisNoAux MCMTSynthesisDefault +Integer comparators +''''''''''''''''''' + +Currently only a single plugin for :class:`.IntegerComparatorGate` is available, which +requires ``num_state_qubits - 1`` clean auxiliary qubits. + +.. autosummary:: + :toctree: ../stubs/ + + IntegerComparatorSynthesisDefault + """ import numpy as np import rustworkx as rx +from qiskit.exceptions import QiskitError from qiskit.circuit.quantumcircuit import QuantumCircuit from qiskit.circuit.library import LinearFunction, QFTGate, MCXGate, C3XGate, C4XGate from qiskit.transpiler.exceptions import TranspilerError from qiskit.transpiler.coupling import CouplingMap +from qiskit.synthesis.arithmetic import synth_integer_comparator_2s from qiskit.synthesis.clifford import ( synth_clifford_full, synth_clifford_layers, @@ -1029,3 +1042,25 @@ def run(self, high_level_object, coupling_map=None, target=None, qubits=None, ** high_level_object.num_target_qubits, ctrl_state, ) + + +class IntegerComparatorSynthesisDefault(HighLevelSynthesisPlugin): + """The default synthesis for ``IntegerComparatorGate``. + + Currently this is only supporting an ancilla-based decomposition. + """ + + def run(self, high_level_object, coupling_map=None, target=None, qubits=None, **options): + num_state_qubits = high_level_object.num_qubits - 1 + num_aux = num_state_qubits - 1 + if (available_aux := options.get("num_clean_ancillas", 0)) < num_aux: + raise QiskitError( + "The IntegerComparatorGate can currently only be synthesized with " + f"num_state_qubits - 1 clean auxiliary qubits. In this case {num_aux} " + f"are required but only {available_aux} are available. Add more clean qubits " + "to succesfully run this synthesis." + ) + + return synth_integer_comparator_2s( + num_state_qubits, high_level_object.value, high_level_object.geq + ) diff --git a/test/python/circuit/library/test_integer_comparator.py b/test/python/circuit/library/test_integer_comparator.py index 66ae2f6aec9d..f9b0eca5f819 100644 --- a/test/python/circuit/library/test_integer_comparator.py +++ b/test/python/circuit/library/test_integer_comparator.py @@ -16,8 +16,9 @@ import numpy as np from ddt import ddt, data, unpack +from qiskit import transpile, QiskitError from qiskit.circuit import QuantumCircuit -from qiskit.circuit.library import IntegerComparator +from qiskit.circuit.library import IntegerComparator, IntegerComparatorGate from qiskit.quantum_info import Statevector from test import QiskitTestCase # pylint: disable=wrong-import-order @@ -28,12 +29,13 @@ class TestIntegerComparator(QiskitTestCase): def assertComparisonIsCorrect(self, comp, num_state_qubits, value, geq): """Assert that the comparator output is correct.""" - qc = QuantumCircuit(comp.num_qubits) # initialize circuit + qc = QuantumCircuit(2 * num_state_qubits) # initialize circuit qc.h(list(range(num_state_qubits))) # set equal superposition state qc.append(comp, list(range(comp.num_qubits))) # add comparator # run simulation - statevector = Statevector(qc) + tqc = transpile(qc) # trigger the HLS if necessary + statevector = Statevector(tqc) for i, amplitude in enumerate(statevector): prob = np.abs(amplitude) ** 2 if prob > 1e-6: @@ -60,8 +62,12 @@ def assertComparisonIsCorrect(self, comp, num_state_qubits, value, geq): def test_fixed_value_comparator(self, num_state_qubits, value, geq): """Test the fixed value comparator circuit.""" # build the circuit with the comparator - comp = IntegerComparator(num_state_qubits, value, geq=geq) - self.assertComparisonIsCorrect(comp, num_state_qubits, value, geq) + for use_gate in [True, False]: + with self.subTest(use_gate=use_gate): + constructor = IntegerComparatorGate if use_gate else IntegerComparator + + comp = constructor(num_state_qubits, value, geq=geq) + self.assertComparisonIsCorrect(comp, num_state_qubits, value, geq) def test_mutability(self): """Test changing the arguments of the comparator.""" @@ -94,6 +100,19 @@ def test_mutability(self): comp.geq = False self.assertComparisonIsCorrect(comp, 3, 2, False) + def test_plugin_warning(self): + """Test the plugin for IntegerComparatorGate warns if there are insufficient aux qubits.""" + + gate = IntegerComparatorGate(2, 3) + circuit = QuantumCircuit(3) + circuit.append(gate, circuit.qubits) + + with self.assertRaisesRegex( + QiskitError, + "The IntegerComparatorGate can currently only be synthesized with num_state_qubits - 1 ", + ): + _ = transpile(circuit) + if __name__ == "__main__": unittest.main() From bf8c69b2c020e24d7dd81ace5095266583c93a9d Mon Sep 17 00:00:00 2001 From: Julien Gacon Date: Thu, 24 Oct 2024 15:33:16 +0200 Subject: [PATCH 2/6] quadratic form no support for parameter expressions for now --- qiskit/circuit/library/__init__.py | 2 + qiskit/circuit/library/arithmetic/__init__.py | 2 +- .../library/arithmetic/quadratic_form.py | 170 +++++++++++++++++- .../circuit/library/test_quadratic_form.py | 32 ++-- 4 files changed, 191 insertions(+), 15 deletions(-) diff --git a/qiskit/circuit/library/__init__.py b/qiskit/circuit/library/__init__.py index 5d5363cd034f..3740350620b1 100644 --- a/qiskit/circuit/library/__init__.py +++ b/qiskit/circuit/library/__init__.py @@ -298,6 +298,7 @@ :template: autosummary/class_no_inherited_members.rst QuadraticForm + QuadraticFormGate Other arithmetic functions -------------------------- @@ -539,6 +540,7 @@ IntegerComparatorGate, WeightedAdder, QuadraticForm, + QuadraticFormGate, LinearAmplitudeFunction, VBERippleCarryAdder, CDKMRippleCarryAdder, diff --git a/qiskit/circuit/library/arithmetic/__init__.py b/qiskit/circuit/library/arithmetic/__init__.py index 0403ecf06d8d..f316b5bb3e3d 100644 --- a/qiskit/circuit/library/arithmetic/__init__.py +++ b/qiskit/circuit/library/arithmetic/__init__.py @@ -19,7 +19,7 @@ from .piecewise_polynomial_pauli_rotations import PiecewisePolynomialPauliRotations from .polynomial_pauli_rotations import PolynomialPauliRotations from .weighted_adder import WeightedAdder -from .quadratic_form import QuadraticForm +from .quadratic_form import QuadraticForm, QuadraticFormGate from .linear_amplitude_function import LinearAmplitudeFunction from .adders import VBERippleCarryAdder, CDKMRippleCarryAdder, DraperQFTAdder from .piecewise_chebyshev import PiecewiseChebyshev diff --git a/qiskit/circuit/library/arithmetic/quadratic_form.py b/qiskit/circuit/library/arithmetic/quadratic_form.py index 5d6fd3218265..8bb540adb545 100644 --- a/qiskit/circuit/library/arithmetic/quadratic_form.py +++ b/qiskit/circuit/library/arithmetic/quadratic_form.py @@ -12,13 +12,18 @@ """A circuit implementing a quadratic form on binary variables.""" +from __future__ import annotations + from typing import Union, Optional, List import math +from collections.abc import Sequence import numpy as np -from qiskit.circuit import QuantumCircuit, QuantumRegister, ParameterExpression -from ..basis_change import QFT +from qiskit.circuit import QuantumCircuit, QuantumRegister, ParameterExpression, Gate, CircuitError +from ..basis_change import QFT, QFTGate + +VALUE_TYPE = Union[int, float, np.integer, np.floating, ParameterExpression] class QuadraticForm(QuantumCircuit): @@ -171,6 +176,101 @@ def required_result_qubits( ) -> int: """Get the number of required result qubits. + Args: + quadratic: A matrix containing the quadratic coefficients. + linear: An array containing the linear coefficients. + offset: A constant offset. + + Returns: + The number of qubits needed to represent the value of the quadratic form + in twos complement. + """ + return QuadraticFormGate.required_result_qubits(quadratic, linear, offset) + + +class QuadraticFormGate(Gate): + r"""Implements a quadratic form on binary variables encoded in qubit registers. + + A quadratic form on binary variables is a quadratic function :math:`Q` acting on a binary + variable of :math:`n` bits, :math:`x = x_0 ... x_{n-1}`. For an integer matrix :math:`A`, + an integer vector :math:`b` and an integer :math:`c` the function can be written as + + .. math:: + + Q(x) = x^T A x + x^T b + c + + If :math:`A`, :math:`b` or :math:`c` contain scalar values, this circuit computes only + an approximation of the quadratic form. + + Provided with :math:`m` qubits to encode the value, this circuit computes :math:`Q(x) \mod 2^m` + in [two's complement](https://stackoverflow.com/questions/1049722/what-is-2s-complement) + representation. + + .. math:: + + |x\rangle_n |0\rangle_m \mapsto |x\rangle_n |(Q(x) + 2^m) \mod 2^m \rangle_m + + Since we use two's complement e.g. the value of :math:`Q(x) = 3` requires 2 bits to represent + the value and 1 bit for the sign: `3 = '011'` where the first `0` indicates a positive value. + On the other hand, :math:`Q(x) = -3` would be `-3 = '101'`, where the first `1` indicates + a negative value and `01` is the two's complement of `3`. + + If the value of :math:`Q(x)` is too large to be represented with `m` qubits, the resulting + bitstring is :math:`(Q(x) + 2^m) \mod 2^m)`. + + The implementation of this circuit is discussed in [1], Fig. 6. + + References: + [1]: Gilliam et al., Grover Adaptive Search for Constrained Polynomial Binary Optimization. + `arXiv:1912.04088 `_ + + """ + + def __init__( + self, + num_result_qubits: int | None = None, + quadratic: Sequence[Sequence[float]] | None = None, + linear: Sequence[Sequence[float]] | None = None, + offset: float | None = None, + label: str = "Q(x)", + ): + # check inputs match + if quadratic is not None and linear is not None: + if len(quadratic) != len(linear): + raise ValueError("Mismatching sizes of quadratic and linear.") + + # temporarily set quadratic and linear to [] instead of None so we can iterate over them + if quadratic is None: + quadratic = [] + + if linear is None: + linear = [] + + if offset is None: + offset = 0 + + self.num_input_qubits = np.max([1, len(linear), len(quadratic)]) + + # deduce number of result bits if not added + if num_result_qubits is None: + num_result_qubits = self.required_result_qubits(quadratic, linear, offset) + + self.num_result_qubits = num_result_qubits + self.quadratic = quadratic + self.linear = linear + self.offset = offset + + num_qubits = int(self.num_input_qubits + self.num_result_qubits) + super().__init__("QuadraticForm", num_qubits, [], label=label) + + @staticmethod + def required_result_qubits( + quadratic: Sequence[Sequence[float]], + linear: Sequence[float], + offset: float, + ) -> int: + """Get the number of required result qubits. + Args: quadratic: A matrix containing the quadratic coefficients. linear: An array containing the linear coefficients. @@ -196,3 +296,69 @@ def required_result_qubits( num_result_qubits = 1 + max(num_qubits_for_min, num_qubits_for_max) return num_result_qubits + + def validate_parameter(self, parameter): + if isinstance(parameter, VALUE_TYPE): + return parameter + + if isinstance(parameter, (np.ndarray, Sequence)): + if all(isinstance(el, VALUE_TYPE) for el in parameter): + return parameter + for params in parameter: + if not all(isinstance(el, VALUE_TYPE) for el in params): + raise CircuitError( + f"Invalid parameter type {type(parameter)} for QuadraticFormGate" + ) + + return parameter + + return super().validate_parameter(parameter) + + def _define(self): + quadratic, linear, offset = self.quadratic, self.linear, self.offset + + qr_input = QuantumRegister(self.num_input_qubits) + qr_result = QuantumRegister(self.num_result_qubits) + circuit = QuantumCircuit(qr_input, qr_result) + + # set quadratic and linear again to None if they were None + if len(quadratic) == 0: + quadratic = None + + if len(linear) == 0: + linear = None + + scaling = np.pi * 2 ** (1 - self.num_result_qubits) + + # initial QFT + qft = QFTGate(self.num_result_qubits) + circuit.append(qft, qr_result) + + # constant coefficient + if offset != 0: + for i, q_i in enumerate(qr_result): + circuit.p(scaling * 2**i * offset, q_i) + + # the linear part consists of the vector and the diagonal of the + # matrix, since x_i * x_i = x_i, as x_i is a binary variable + for j in range(self.num_input_qubits): + value = linear[j] if linear is not None else 0 + value += quadratic[j][j] if quadratic is not None else 0 + if value != 0: + for i, q_i in enumerate(qr_result): + circuit.cp(scaling * 2**i * value, qr_input[j], q_i) + + # the quadratic part adds A_ij and A_ji as x_i x_j == x_j x_i + if quadratic is not None: + for j in range(self.num_input_qubits): + for k in range(j + 1, self.num_input_qubits): + value = quadratic[j][k] + quadratic[k][j] + if value != 0: + for i, q_i in enumerate(qr_result): + circuit.mcp(scaling * 2**i * value, [qr_input[j], qr_input[k]], q_i) + + # add the inverse QFT + iqft = qft.inverse() + circuit.append(iqft, qr_result) + + self.definition = circuit diff --git a/test/python/circuit/library/test_quadratic_form.py b/test/python/circuit/library/test_quadratic_form.py index a97b030f0d28..7212a0a22e0b 100644 --- a/test/python/circuit/library/test_quadratic_form.py +++ b/test/python/circuit/library/test_quadratic_form.py @@ -17,7 +17,7 @@ import numpy as np from qiskit.circuit import QuantumCircuit, ParameterVector -from qiskit.circuit.library import QuadraticForm +from qiskit.circuit.library import QuadraticForm, QuadraticFormGate from qiskit.quantum_info import Statevector from test import QiskitTestCase # pylint: disable=wrong-import-order @@ -68,32 +68,38 @@ def test_endian(self, little_endian): self.assertTrue(Statevector.from_instruction(circuit).equiv(ref)) - def test_required_result_qubits(self): + @data(True, False) + def test_required_result_qubits(self, use_gate): """Test getting the number of required result qubits.""" + constructor = QuadraticFormGate if use_gate else QuadraticForm + with self.subTest("positive bound"): quadratic = [[1, -50], [100, 0]] linear = [-5, 5] offset = 0 - num_result_qubits = QuadraticForm.required_result_qubits(quadratic, linear, offset) + num_result_qubits = constructor.required_result_qubits(quadratic, linear, offset) self.assertEqual(num_result_qubits, 1 + int(np.ceil(np.log2(106 + 1)))) with self.subTest("negative bound"): quadratic = [[1, -50], [10, 0]] linear = [-5, 5] offset = 0 - num_result_qubits = QuadraticForm.required_result_qubits(quadratic, linear, offset) + num_result_qubits = constructor.required_result_qubits(quadratic, linear, offset) self.assertEqual(num_result_qubits, 1 + int(np.ceil(np.log2(55)))) with self.subTest("empty"): - num_result_qubits = QuadraticForm.required_result_qubits([[]], [], 0) + num_result_qubits = constructor.required_result_qubits([[]], [], 0) self.assertEqual(num_result_qubits, 1) - def test_quadratic_form(self): + @data(True, False) + def test_quadratic_form(self, use_gate): """Test the quadratic form circuit.""" + constructor = QuadraticFormGate if use_gate else QuadraticForm + with self.subTest("empty"): - circuit = QuadraticForm() + circuit = constructor() self.assertQuadraticFormIsCorrect(1, [[0]], [0], 0, circuit) with self.subTest("1d case"): @@ -101,7 +107,7 @@ def test_quadratic_form(self): linear = np.array([2]) offset = -1 - circuit = QuadraticForm(quadratic=quadratic, linear=linear, offset=offset) + circuit = constructor(quadratic=quadratic, linear=linear, offset=offset) self.assertQuadraticFormIsCorrect(3, quadratic, linear, offset, circuit) @@ -111,7 +117,7 @@ def test_quadratic_form(self): offset = -1 m = 2 - circuit = QuadraticForm(m, quadratic, linear, offset) + circuit = constructor(m, quadratic, linear, offset) self.assertQuadraticFormIsCorrect(m, quadratic, linear, offset, circuit) @@ -120,7 +126,7 @@ def test_quadratic_form(self): linear = np.array([-2, 0, 1]) offset = -1 - circuit = QuadraticForm(linear=linear, offset=offset) + circuit = constructor(linear=linear, offset=offset) self.assertQuadraticFormIsCorrect(3, quadratic, linear, offset, circuit) with self.subTest("missing linear"): @@ -129,7 +135,7 @@ def test_quadratic_form(self): offset = -1 m = 2 - circuit = QuadraticForm(m, quadratic, None, offset) + circuit = constructor(m, quadratic, None, offset) self.assertQuadraticFormIsCorrect(m, quadratic, linear, offset, circuit) with self.subTest("missing offset"): @@ -138,11 +144,12 @@ def test_quadratic_form(self): offset = 0 m = 2 - circuit = QuadraticForm(m, quadratic, linear) + circuit = constructor(m, quadratic, linear) self.assertQuadraticFormIsCorrect(m, quadratic, linear, offset, circuit) def test_quadratic_form_parameterized(self): """Test the quadratic form circuit with parameters.""" + theta = ParameterVector("th", 7) p_quadratic = [[theta[0], theta[1]], [theta[2], theta[3]]] @@ -156,6 +163,7 @@ def test_quadratic_form_parameterized(self): circuit = QuadraticForm(m, p_quadratic, p_linear, p_offset) param_dict = dict(zip(theta, [*quadratic[0]] + [*quadratic[1]] + [*linear] + [offset])) + circuit.assign_parameters(param_dict, inplace=True) self.assertQuadraticFormIsCorrect(m, quadratic, linear, offset, circuit) From c787d7a2ff697c6877299eb93a6d1cabc8481772 Mon Sep 17 00:00:00 2001 From: Julien Gacon Date: Thu, 24 Oct 2024 18:28:07 +0200 Subject: [PATCH 3/6] no aux integer comparator & linear paulis --- pyproject.toml | 4 +- qiskit/circuit/library/__init__.py | 2 + qiskit/circuit/library/arithmetic/__init__.py | 5 +- .../library/arithmetic/integer_comparator.py | 2 +- .../arithmetic/linear_pauli_rotations.py | 64 +++++++++- .../piecewise_linear_pauli_rotations.py | 118 +++++++++++++++++- qiskit/synthesis/arithmetic/__init__.py | 2 +- .../arithmetic/comparators/__init__.py | 3 +- .../{integer_comparator.py => compare_2s.py} | 0 .../arithmetic/comparators/compare_greedy.py | 46 +++++++ .../passes/synthesis/hls_plugins.py | 38 ++++-- .../test_functional_pauli_rotations.py | 44 +++++-- .../library/test_integer_comparator.py | 32 ++--- 13 files changed, 310 insertions(+), 50 deletions(-) rename qiskit/synthesis/arithmetic/comparators/{integer_comparator.py => compare_2s.py} (100%) create mode 100644 qiskit/synthesis/arithmetic/comparators/compare_greedy.py diff --git a/pyproject.toml b/pyproject.toml index 0f8156c46d58..8ba7a821c6b9 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -100,7 +100,9 @@ sk = "qiskit.transpiler.passes.synthesis.solovay_kitaev_synthesis:SolovayKitaevS "qft.line" = "qiskit.transpiler.passes.synthesis.hls_plugins:QFTSynthesisLine" "qft.default" = "qiskit.transpiler.passes.synthesis.hls_plugins:QFTSynthesisFull" "permutation.token_swapper" = "qiskit.transpiler.passes.synthesis.hls_plugins:TokenSwapperSynthesisPermutation" -"IntegerComparator.default" = "qiskit.transpiler.passes.synthesis.hls_plugins:IntegerComparatorSynthesisDefault" +"IntComp.default" = "qiskit.transpiler.passes.synthesis.hls_plugins:IntComparatorSynthesisDefault" +"IntComp.noaux" = "qiskit.transpiler.passes.synthesis.hls_plugins:IntComparatorSynthesisNoAux" +"IntComp.twos" = "qiskit.transpiler.passes.synthesis.hls_plugins:IntComparatorSynthesis2s" [project.entry-points."qiskit.transpiler.init"] default = "qiskit.transpiler.preset_passmanagers.builtin_plugins:DefaultInitPassManager" diff --git a/qiskit/circuit/library/__init__.py b/qiskit/circuit/library/__init__.py index 3740350620b1..7d4b1837362a 100644 --- a/qiskit/circuit/library/__init__.py +++ b/qiskit/circuit/library/__init__.py @@ -255,6 +255,7 @@ LinearPauliRotations PolynomialPauliRotations PiecewiseLinearPauliRotations + PiecewiseLinearPauliRotationsGate PiecewisePolynomialPauliRotations PiecewiseChebyshev @@ -534,6 +535,7 @@ FunctionalPauliRotations, LinearPauliRotations, PiecewiseLinearPauliRotations, + PiecewiseLinearPauliRotationsGate, PiecewisePolynomialPauliRotations, PolynomialPauliRotations, IntegerComparator, diff --git a/qiskit/circuit/library/arithmetic/__init__.py b/qiskit/circuit/library/arithmetic/__init__.py index f316b5bb3e3d..721e4cc753ee 100644 --- a/qiskit/circuit/library/arithmetic/__init__.py +++ b/qiskit/circuit/library/arithmetic/__init__.py @@ -15,7 +15,10 @@ from .functional_pauli_rotations import FunctionalPauliRotations from .integer_comparator import IntegerComparator, IntegerComparatorGate from .linear_pauli_rotations import LinearPauliRotations -from .piecewise_linear_pauli_rotations import PiecewiseLinearPauliRotations +from .piecewise_linear_pauli_rotations import ( + PiecewiseLinearPauliRotations, + PiecewiseLinearPauliRotationsGate, +) from .piecewise_polynomial_pauli_rotations import PiecewisePolynomialPauliRotations from .polynomial_pauli_rotations import PolynomialPauliRotations from .weighted_adder import WeightedAdder diff --git a/qiskit/circuit/library/arithmetic/integer_comparator.py b/qiskit/circuit/library/arithmetic/integer_comparator.py index e1f3b3e42fca..0435d26be1ef 100644 --- a/qiskit/circuit/library/arithmetic/integer_comparator.py +++ b/qiskit/circuit/library/arithmetic/integer_comparator.py @@ -187,6 +187,6 @@ def __init__( num_state_qubits: The number of qubits in the registers. geq: If ``True`` compute :math:`i \geq L`, otherwise compute :math:`i < L`. """ - super().__init__("IntegerComparator", num_state_qubits + 1, [], label=label) + super().__init__("IntComp", num_state_qubits + 1, [], label=label) self.value = value self.geq = geq diff --git a/qiskit/circuit/library/arithmetic/linear_pauli_rotations.py b/qiskit/circuit/library/arithmetic/linear_pauli_rotations.py index bc80ef778616..bbfd8d56a490 100644 --- a/qiskit/circuit/library/arithmetic/linear_pauli_rotations.py +++ b/qiskit/circuit/library/arithmetic/linear_pauli_rotations.py @@ -13,9 +13,10 @@ """Linearly-controlled X, Y or Z rotation.""" +from __future__ import annotations from typing import Optional -from qiskit.circuit import QuantumRegister, QuantumCircuit +from qiskit.circuit import QuantumRegister, QuantumCircuit, Gate from qiskit.circuit.exceptions import CircuitError from .functional_pauli_rotations import FunctionalPauliRotations @@ -164,12 +165,65 @@ def _build(self): return super()._build() + gate = LinearPauliRotationsGate(self.num_state_qubits, self.slope, self.offset, self.basis) + self.append(gate, self.qubits) - circuit = QuantumCircuit(*self.qregs, name=self.name) + +class LinearPauliRotationsGate(Gate): + r"""Linearly-controlled X, Y or Z rotation. + + For a register of state qubits :math:`|x\rangle`, a target qubit :math:`|0\rangle` and the + basis ``'Y'`` this circuit acts as: + + .. parsed-literal:: + + q_0: ─────────────────────────■───────── ... ────────────────────── + │ + . + │ + q_(n-1): ─────────────────────────┼───────── ... ───────────■────────── + ┌────────────┐ ┌───────┴───────┐ ┌─────────┴─────────┐ + q_n: ─┤ RY(offset) ├──┤ RY(2^0 slope) ├ ... ┤ RY(2^(n-1) slope) ├ + └────────────┘ └───────────────┘ └───────────────────┘ + + This can for example be used to approximate linear functions, with :math:`a =` ``slope``:math:`/2` + and :math:`b =` ``offset``:math:`/2` and the basis ``'Y'``: + + .. math:: + + |x\rangle |0\rangle \mapsto \cos(ax + b)|x\rangle|0\rangle + \sin(ax + b)|x\rangle |1\rangle + + Since for small arguments :math:`\sin(x) \approx x` this operator can be used to approximate + linear functions. + """ + + def __init__( + self, + num_state_qubits: int, + slope: float = 1, + offset: float = 0, + basis: str = "Y", + label: str | None = None, + ) -> None: + r""" + Args: + num_state_qubits: The number of qubits representing the state :math:`|x\rangle`. + slope: The slope of the controlled rotation. + offset: The offset of the controlled rotation. + basis: The type of Pauli rotation ('X', 'Y', 'Z'). + label: The label of the gate. + """ + super().__init__("LinPauliRot", num_state_qubits + 1, [], label=label) + self.slope = slope + self.offset = offset + self.basis = basis.lower() + + def _define(self): + circuit = QuantumCircuit(self.num_qubits, name=self.name) # build the circuit - qr_state = self.qubits[: self.num_state_qubits] - qr_target = self.qubits[self.num_state_qubits] + qr_state = circuit.qubits[: self.num_qubits - 1] + qr_target = circuit.qubits[-1] if self.basis == "x": circuit.rx(self.offset, qr_target) @@ -186,4 +240,4 @@ def _build(self): else: # 'Z' circuit.crz(self.slope * pow(2, i), q_i, qr_target) - self.append(circuit.to_gate(), self.qubits) + self.definition = circuit diff --git a/qiskit/circuit/library/arithmetic/piecewise_linear_pauli_rotations.py b/qiskit/circuit/library/arithmetic/piecewise_linear_pauli_rotations.py index 3d84e64ccb1b..752dfedfa06f 100644 --- a/qiskit/circuit/library/arithmetic/piecewise_linear_pauli_rotations.py +++ b/qiskit/circuit/library/arithmetic/piecewise_linear_pauli_rotations.py @@ -14,14 +14,15 @@ """Piecewise-linearly-controlled rotation.""" from __future__ import annotations +from collections.abc import Sequence import numpy as np -from qiskit.circuit import QuantumRegister, AncillaRegister, QuantumCircuit +from qiskit.circuit import QuantumRegister, AncillaRegister, QuantumCircuit, Gate from qiskit.circuit.exceptions import CircuitError from .functional_pauli_rotations import FunctionalPauliRotations -from .linear_pauli_rotations import LinearPauliRotations -from .integer_comparator import IntegerComparator +from .linear_pauli_rotations import LinearPauliRotations, LinearPauliRotationsGate +from .integer_comparator import IntegerComparator, IntegerComparatorGate class PiecewiseLinearPauliRotations(FunctionalPauliRotations): @@ -272,6 +273,115 @@ def _build(self): circuit.append(lin_r.to_gate().control(), qr_compare[:] + qr_state[:] + qr_target) # uncompute comparator - circuit.append(comp.to_gate().inverse(), qr[:] + qr_helper[: comp.num_ancillas]) + circuit.append(comp.to_gate(), qr[:] + qr_helper[: comp.num_ancillas]) self.append(circuit.to_gate(), self.qubits) + + +class PiecewiseLinearPauliRotationsGate(Gate): + r"""Piecewise-linearly-controlled Pauli rotations. + + For a piecewise linear (not necessarily continuous) function :math:`f(x)`, which is defined + through breakpoints, slopes and offsets as follows. + Suppose the breakpoints :math:`(x_0, ..., x_J)` are a subset of :math:`[0, 2^n-1]`, where + :math:`n` is the number of state qubits. Further on, denote the corresponding slopes and + offsets by :math:`a_j` and :math:`b_j` respectively. + Then f(x) is defined as: + + .. math:: + + f(x) = \begin{cases} + 0, x < x_0 \\ + a_j (x - x_j) + b_j, x_j \leq x < x_{j+1} + \end{cases} + + where we implicitly assume :math:`x_{J+1} = 2^n`. + """ + + def __init__( + self, + num_state_qubits: int | None = None, + breakpoints: list[int] | None = None, + slopes: Sequence[float] | None = None, + offsets: Sequence[float] | None = None, + basis: str = "Y", + label: str | None = None, + ) -> None: + """Construct piecewise-linearly-controlled Pauli rotations. + + Args: + num_state_qubits: The number of qubits representing the state. + breakpoints: The breakpoints to define the piecewise-linear function. + Defaults to ``[0]``. + slopes: The slopes for different segments of the piecewise-linear function. + Defaults to ``[1]``. + offsets: The offsets for different segments of the piecewise-linear function. + Defaults to ``[0]``. + basis: The type of Pauli rotation (``'X'``, ``'Y'``, ``'Z'``). + label: The label of the gate. + """ + self.breakpoints = breakpoints if breakpoints is not None else [0] + self.slopes = slopes if slopes is not None else [1] + self.offsets = offsets if offsets is not None else [0] + self.basis = basis + + num_compare_bits = 1 if len(self.breakpoints) > 1 else 0 + super().__init__("PwLinPauliRot", num_state_qubits + 1 + num_compare_bits, [], label=label) + + def _define(self): + circuit = QuantumCircuit(self.num_qubits, name=self.name) + + if len(self.breakpoints) == 1: + qr_state = circuit.qubits[: self.num_qubits - 1] + qr_target = [circuit.qubits[-1]] + qr_compare = [] + else: + qr_state = circuit.qubits[: self.num_qubits - 2] + qr_target = [circuit.qubits[-2]] + qr_compare = [circuit.qubits[-1]] + + num_state_qubits = len(qr_state) + + mapped_slopes = np.zeros_like(self.slopes) + for i, slope in enumerate(self.slopes): + mapped_slopes[i] = slope - sum(mapped_slopes[:i]) + + mapped_offsets = np.zeros_like(self.offsets) + for i, (offset, slope, point) in enumerate( + zip(self.offsets, self.slopes, self.breakpoints) + ): + mapped_offsets[i] = offset - slope * point - sum(mapped_offsets[:i]) + + # apply comparators and controlled linear rotations + contains_zero_breakpoint = np.isclose(self.breakpoints[0], 0) + for i, point in enumerate(self.breakpoints): + if i == 0 and contains_zero_breakpoint: + # apply rotation + lin_r = LinearPauliRotationsGate( + num_state_qubits=num_state_qubits, + slope=mapped_slopes[i], + offset=mapped_offsets[i], + basis=self.basis, + ) + circuit.append(lin_r, qr_state[:] + qr_target) + + else: + # apply Comparator + comp = IntegerComparatorGate(num_state_qubits=num_state_qubits, value=point) + qr = qr_state[:] + qr_compare[:] # add ancilla as compare qubit + + circuit.append(comp, qr[:]) + + # apply controlled rotation + lin_r = LinearPauliRotationsGate( + num_state_qubits=num_state_qubits, + slope=mapped_slopes[i], + offset=mapped_offsets[i], + basis=self.basis, + ) + circuit.append(lin_r.control(), qr_compare[:] + qr_state[:] + qr_target) + + # uncompute comparator (which is its self-inverse) + circuit.append(comp, qr[:]) + + self.definition = circuit diff --git a/qiskit/synthesis/arithmetic/__init__.py b/qiskit/synthesis/arithmetic/__init__.py index 89013c167c03..8e853e86cbad 100644 --- a/qiskit/synthesis/arithmetic/__init__.py +++ b/qiskit/synthesis/arithmetic/__init__.py @@ -12,4 +12,4 @@ """Synthesis for arithmetic circuits.""" -from .comparators import synth_integer_comparator_2s +from .comparators import synth_integer_comparator_2s, synth_integer_comparator_greedy diff --git a/qiskit/synthesis/arithmetic/comparators/__init__.py b/qiskit/synthesis/arithmetic/comparators/__init__.py index b492d5563ee9..f14070800f7a 100644 --- a/qiskit/synthesis/arithmetic/comparators/__init__.py +++ b/qiskit/synthesis/arithmetic/comparators/__init__.py @@ -12,4 +12,5 @@ """Comparator synthesis algorithms.""" -from .integer_comparator import synth_integer_comparator_2s +from .compare_2s import synth_integer_comparator_2s +from .compare_greedy import synth_integer_comparator_greedy diff --git a/qiskit/synthesis/arithmetic/comparators/integer_comparator.py b/qiskit/synthesis/arithmetic/comparators/compare_2s.py similarity index 100% rename from qiskit/synthesis/arithmetic/comparators/integer_comparator.py rename to qiskit/synthesis/arithmetic/comparators/compare_2s.py diff --git a/qiskit/synthesis/arithmetic/comparators/compare_greedy.py b/qiskit/synthesis/arithmetic/comparators/compare_greedy.py new file mode 100644 index 000000000000..e81f55da0a97 --- /dev/null +++ b/qiskit/synthesis/arithmetic/comparators/compare_greedy.py @@ -0,0 +1,46 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2017, 2024. +# +# This code is licensed under the Apache License, Version 2.0. You may +# obtain a copy of this license in the LICENSE.txt file in the root directory +# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0. +# +# Any modifications or derivative works of this code must retain this +# copyright notice, and modified files need to carry a notice indicating +# that they have been altered from the originals. + +"""Integer comparator based on 2s complement.""" + +from qiskit.circuit import QuantumCircuit + + +def synth_integer_comparator_greedy( + num_state_qubits: int, value: int, geq: bool = True +) -> QuantumCircuit: + """Implement an integer comparison based on value-by-value comparison.""" + circuit = QuantumCircuit(num_state_qubits + 1) + + if value <= 0: # condition always satisfied for non-positive values + if geq: # otherwise the condition is never satisfied + circuit.x(num_state_qubits) + + return circuit + + # make sure to always choose the comparison where we have to place less than + # (2 ** n)/2 MCX gates + if (value < 2 ** (num_state_qubits - 1) and geq) or ( + value > 2 ** (num_state_qubits - 1) and not geq + ): + geq = not geq + circuit.x(num_state_qubits) + + if geq: + accepted_values = range(value, 2**num_state_qubits) + else: + accepted_values = range(0, value) + + for accepted_value in accepted_values: + circuit.mcx(list(range(num_state_qubits)), num_state_qubits, ctrl_state=accepted_value) + + return circuit diff --git a/qiskit/transpiler/passes/synthesis/hls_plugins.py b/qiskit/transpiler/passes/synthesis/hls_plugins.py index fb3000907d0f..b85a2c608918 100644 --- a/qiskit/transpiler/passes/synthesis/hls_plugins.py +++ b/qiskit/transpiler/passes/synthesis/hls_plugins.py @@ -258,20 +258,19 @@ .. autosummary:: :toctree: ../stubs/ - IntegerComparatorSynthesisDefault + IntComparatorSynthesisDefault """ import numpy as np import rustworkx as rx -from qiskit.exceptions import QiskitError from qiskit.circuit.quantumcircuit import QuantumCircuit from qiskit.circuit.library import LinearFunction, QFTGate, MCXGate, C3XGate, C4XGate from qiskit.transpiler.exceptions import TranspilerError from qiskit.transpiler.coupling import CouplingMap -from qiskit.synthesis.arithmetic import synth_integer_comparator_2s +from qiskit.synthesis.arithmetic import synth_integer_comparator_2s, synth_integer_comparator_greedy from qiskit.synthesis.clifford import ( synth_clifford_full, synth_clifford_layers, @@ -1044,7 +1043,7 @@ def run(self, high_level_object, coupling_map=None, target=None, qubits=None, ** ) -class IntegerComparatorSynthesisDefault(HighLevelSynthesisPlugin): +class IntComparatorSynthesisDefault(HighLevelSynthesisPlugin): """The default synthesis for ``IntegerComparatorGate``. Currently this is only supporting an ancilla-based decomposition. @@ -1053,14 +1052,33 @@ class IntegerComparatorSynthesisDefault(HighLevelSynthesisPlugin): def run(self, high_level_object, coupling_map=None, target=None, qubits=None, **options): num_state_qubits = high_level_object.num_qubits - 1 num_aux = num_state_qubits - 1 - if (available_aux := options.get("num_clean_ancillas", 0)) < num_aux: - raise QiskitError( - "The IntegerComparatorGate can currently only be synthesized with " - f"num_state_qubits - 1 clean auxiliary qubits. In this case {num_aux} " - f"are required but only {available_aux} are available. Add more clean qubits " - "to succesfully run this synthesis." + if options.get("num_clean_ancillas", 0) < num_aux: + return synth_integer_comparator_greedy( + num_state_qubits, high_level_object.value, high_level_object.geq ) return synth_integer_comparator_2s( num_state_qubits, high_level_object.value, high_level_object.geq ) + + +class IntComparatorSynthesisNoAux(HighLevelSynthesisPlugin): + """A potentially exponentially expensive comparison w/o auxiliary qubits.""" + + def run(self, high_level_object, coupling_map=None, target=None, qubits=None, **options): + return synth_integer_comparator_greedy( + high_level_object.num_state_qubits, high_level_object.value, high_level_object.geq + ) + + +class IntComparatorSynthesis2s(HighLevelSynthesisPlugin): + """An integer comparison based on 2s complement.""" + + def run(self, high_level_object, coupling_map=None, target=None, qubits=None, **options): + num_aux = high_level_object.num_state_qubits - 1 + if options.get("num_clean_ancillas", 0) < num_aux: + return None + + return synth_integer_comparator_2s( + high_level_object.num_state_qubits, high_level_object.value, high_level_object.geq + ) diff --git a/test/python/circuit/library/test_functional_pauli_rotations.py b/test/python/circuit/library/test_functional_pauli_rotations.py index cb2ba2e09a10..501c2f80fdfe 100644 --- a/test/python/circuit/library/test_functional_pauli_rotations.py +++ b/test/python/circuit/library/test_functional_pauli_rotations.py @@ -17,11 +17,13 @@ import numpy as np from ddt import ddt, data, unpack +from qiskit import transpile from qiskit.circuit import QuantumCircuit from qiskit.circuit.library import ( LinearPauliRotations, PolynomialPauliRotations, PiecewiseLinearPauliRotations, + PiecewiseLinearPauliRotationsGate, ) from qiskit.quantum_info import Statevector from test import QiskitTestCase # pylint: disable=wrong-import-order @@ -31,14 +33,19 @@ class TestFunctionalPauliRotations(QiskitTestCase): """Test the functional Pauli rotations.""" - def assertFunctionIsCorrect(self, function_circuit, reference): + def assertFunctionIsCorrect(self, function_circuit, reference, num_ancilla_qubits=None): """Assert that ``function_circuit`` implements the reference function ``reference``.""" - num_state_qubits = function_circuit.num_qubits - function_circuit.num_ancillas - 1 - num_ancilla_qubits = function_circuit.num_ancillas + if isinstance(function_circuit, QuantumCircuit): + num_ancilla_qubits = function_circuit.num_ancillas + + num_state_qubits = function_circuit.num_qubits - 1 - num_ancilla_qubits + circuit = QuantumCircuit(num_state_qubits + 1 + num_ancilla_qubits) circuit.h(list(range(num_state_qubits))) - circuit.append(function_circuit.to_instruction(), list(range(circuit.num_qubits))) - statevector = Statevector(circuit) + circuit.compose(function_circuit, list(range(function_circuit.num_qubits)), inplace=True) + + tqc = transpile(circuit, basis_gates=["u", "cx"]) + statevector = Statevector(tqc) probabilities = defaultdict(float) for i, statevector_amplitude in enumerate(statevector): i = bin(i)[2:].zfill(circuit.num_qubits)[num_ancilla_qubits:] @@ -152,14 +159,27 @@ def pw_linear(x): return offsets[-(i + 1)] + slopes[-(i + 1)] * (x - point) return 0 - pw_linear_rotations = PiecewiseLinearPauliRotations( - num_state_qubits, - breakpoints, - [2 * slope for slope in slopes], - [2 * offset for offset in offsets], - ) + for use_gate in [False, True]: + constructor = ( + PiecewiseLinearPauliRotationsGate if use_gate else PiecewiseLinearPauliRotations + ) + + if use_gate: + # ancilla for the comparator bit and the integer comparator itself + num_ancillas = int(len(breakpoints) > 1) # + num_state_qubits - 1 + # num_ancillas = 0 + else: + num_ancillas = None # automatically deducted + + with self.subTest(use_gate=use_gate): + pw_linear_rotations = constructor( + num_state_qubits, + breakpoints, + [2 * slope for slope in slopes], + [2 * offset for offset in offsets], + ) - self.assertFunctionIsCorrect(pw_linear_rotations, pw_linear) + self.assertFunctionIsCorrect(pw_linear_rotations, pw_linear, num_ancillas) def test_piecewise_linear_rotations_mutability(self): """Test the mutability of the linear rotations circuit.""" diff --git a/test/python/circuit/library/test_integer_comparator.py b/test/python/circuit/library/test_integer_comparator.py index f9b0eca5f819..d36a731d2584 100644 --- a/test/python/circuit/library/test_integer_comparator.py +++ b/test/python/circuit/library/test_integer_comparator.py @@ -19,6 +19,7 @@ from qiskit import transpile, QiskitError from qiskit.circuit import QuantumCircuit from qiskit.circuit.library import IntegerComparator, IntegerComparatorGate +from qiskit.synthesis.arithmetic import synth_integer_comparator_2s, synth_integer_comparator_greedy from qiskit.quantum_info import Statevector from test import QiskitTestCase # pylint: disable=wrong-import-order @@ -62,10 +63,12 @@ def assertComparisonIsCorrect(self, comp, num_state_qubits, value, geq): def test_fixed_value_comparator(self, num_state_qubits, value, geq): """Test the fixed value comparator circuit.""" # build the circuit with the comparator - for use_gate in [True, False]: - with self.subTest(use_gate=use_gate): - constructor = IntegerComparatorGate if use_gate else IntegerComparator - + for constructor in [ + IntegerComparator, + synth_integer_comparator_2s, + synth_integer_comparator_greedy, + ]: + with self.subTest(constructor=constructor): comp = constructor(num_state_qubits, value, geq=geq) self.assertComparisonIsCorrect(comp, num_state_qubits, value, geq) @@ -100,18 +103,19 @@ def test_mutability(self): comp.geq = False self.assertComparisonIsCorrect(comp, 3, 2, False) - def test_plugin_warning(self): - """Test the plugin for IntegerComparatorGate warns if there are insufficient aux qubits.""" + # def test_plugin_warning(self): + # """Test the plugin for IntegerComparatorGate warns if there are insufficient aux qubits.""" - gate = IntegerComparatorGate(2, 3) - circuit = QuantumCircuit(3) - circuit.append(gate, circuit.qubits) + # gate = IntegerComparatorGate(2, 2) + # circuit = QuantumCircuit(3) + # circuit.append(gate, circuit.qubits) - with self.assertRaisesRegex( - QiskitError, - "The IntegerComparatorGate can currently only be synthesized with num_state_qubits - 1 ", - ): - _ = transpile(circuit) + # with self.assertRaisesRegex( + # QiskitError, + # "The IntegerComparatorGate can currently only be synthesized with num_state_qubits - 1 ", + # ): + # _ = transpile(circuit, basis_gates=["u", "cx", "ccx", "x"]) + # print(_) if __name__ == "__main__": From 3d99f29ee9cf59ea37b78b73f43ec1319c5b2577 Mon Sep 17 00:00:00 2001 From: Julien Gacon Date: Thu, 24 Oct 2024 18:43:33 +0200 Subject: [PATCH 4/6] almost do LinAmpFunction some tests miss, some HLS tests fail (maybe Sasha's PR fixes this) --- qiskit/circuit/library/__init__.py | 2 + qiskit/circuit/library/arithmetic/__init__.py | 2 +- .../arithmetic/linear_amplitude_function.py | 169 +++++++++++++++++- .../library/test_integer_comparator.py | 18 +- .../library/test_linear_amplitude_function.py | 26 ++- 5 files changed, 189 insertions(+), 28 deletions(-) diff --git a/qiskit/circuit/library/__init__.py b/qiskit/circuit/library/__init__.py index 7d4b1837362a..324380928fdf 100644 --- a/qiskit/circuit/library/__init__.py +++ b/qiskit/circuit/library/__init__.py @@ -243,6 +243,7 @@ :template: autosummary/class_no_inherited_members.rst LinearAmplitudeFunction + LinearAmplitudeFunctionGate Functional Pauli Rotations -------------------------- @@ -544,6 +545,7 @@ QuadraticForm, QuadraticFormGate, LinearAmplitudeFunction, + LinearAmplitudeFunctionGate, VBERippleCarryAdder, CDKMRippleCarryAdder, DraperQFTAdder, diff --git a/qiskit/circuit/library/arithmetic/__init__.py b/qiskit/circuit/library/arithmetic/__init__.py index 721e4cc753ee..c273f2214328 100644 --- a/qiskit/circuit/library/arithmetic/__init__.py +++ b/qiskit/circuit/library/arithmetic/__init__.py @@ -23,7 +23,7 @@ from .polynomial_pauli_rotations import PolynomialPauliRotations from .weighted_adder import WeightedAdder from .quadratic_form import QuadraticForm, QuadraticFormGate -from .linear_amplitude_function import LinearAmplitudeFunction +from .linear_amplitude_function import LinearAmplitudeFunction, LinearAmplitudeFunctionGate from .adders import VBERippleCarryAdder, CDKMRippleCarryAdder, DraperQFTAdder from .piecewise_chebyshev import PiecewiseChebyshev from .multipliers import HRSCumulativeMultiplier, RGQFTMultiplier diff --git a/qiskit/circuit/library/arithmetic/linear_amplitude_function.py b/qiskit/circuit/library/arithmetic/linear_amplitude_function.py index a79670eef654..fdc603a3cdc7 100644 --- a/qiskit/circuit/library/arithmetic/linear_amplitude_function.py +++ b/qiskit/circuit/library/arithmetic/linear_amplitude_function.py @@ -14,9 +14,12 @@ from __future__ import annotations import numpy as np -from qiskit.circuit import QuantumCircuit +from qiskit.circuit import QuantumCircuit, Gate -from .piecewise_linear_pauli_rotations import PiecewiseLinearPauliRotations +from .piecewise_linear_pauli_rotations import ( + PiecewiseLinearPauliRotations, + PiecewiseLinearPauliRotationsGate, +) class LinearAmplitudeFunction(QuantumCircuit): @@ -173,6 +176,168 @@ def post_processing(self, scaled_value: float) -> float: return value +class LinearAmplitudeFunctionGate(Gate): + r"""A circuit implementing a (piecewise) linear function on qubit amplitudes. + + An amplitude function :math:`F` of a function :math:`f` is a mapping + + .. math:: + + F|x\rangle|0\rangle = \sqrt{1 - \hat{f}(x)} |x\rangle|0\rangle + \sqrt{\hat{f}(x)} + |x\rangle|1\rangle. + + for a function :math:`\hat{f}: \{ 0, ..., 2^n - 1 \} \rightarrow [0, 1]`, where + :math:`|x\rangle` is a :math:`n` qubit state. + + This circuit implements :math:`F` for piecewise linear functions :math:`\hat{f}`. + In this case, the mapping :math:`F` can be approximately implemented using a Taylor expansion + and linearly controlled Pauli-Y rotations, see [1, 2] for more detail. This approximation + uses a ``rescaling_factor`` to determine the accuracy of the Taylor expansion. + + In general, the function of interest :math:`f` is defined from some interval :math:`[a,b]`, + the ``domain`` to :math:`[c,d]`, the ``image``, instead of :math:`\{ 1, ..., N \}` to + :math:`[0, 1]`. Using an affine transformation we can rescale :math:`f` to :math:`\hat{f}`: + + .. math:: + + \hat{f}(x) = \frac{f(\phi(x)) - c}{d - c} + + with + + .. math:: + + \phi(x) = a + \frac{b - a}{2^n - 1} x. + + If :math:`f` is a piecewise linear function on :math:`m` intervals + :math:`[p_{i-1}, p_i], i \in \{1, ..., m\}` with slopes :math:`\alpha_i` and + offsets :math:`\beta_i` it can be written as + + .. math:: + + f(x) = \sum_{i=1}^m 1_{[p_{i-1}, p_i]}(x) (\alpha_i x + \beta_i) + + where :math:`1_{[a, b]}` is an indication function that is 1 if the argument is in the interval + :math:`[a, b]` and otherwise 0. The breakpoints :math:`p_i` can be specified by the + ``breakpoints`` argument. + + References: + + [1]: Woerner, S., & Egger, D. J. (2018). + Quantum Risk Analysis. + `arXiv:1806.06893 `_ + + [2]: Gacon, J., Zoufal, C., & Woerner, S. (2020). + Quantum-Enhanced Simulation-Based Optimization. + `arXiv:2005.10780 `_ + """ + + def __init__( + self, + num_state_qubits: int, + slope: float | list[float], + offset: float | list[float], + domain: tuple[float, float], + image: tuple[float, float], + rescaling_factor: float = 1, + breakpoints: list[float] | None = None, + label: str = "F", + ) -> None: + r""" + Args: + num_state_qubits: The number of qubits used to encode the variable :math:`x`. + slope: The slope of the linear function. Can be a list of slopes if it is a piecewise + linear function. + offset: The offset of the linear function. Can be a list of offsets if it is a piecewise + linear function. + domain: The domain of the function as tuple :math:`(x_\min{}, x_\max{})`. + image: The image of the function as tuple :math:`(f_\min{}, f_\max{})`. + rescaling_factor: The rescaling factor to adjust the accuracy in the Taylor + approximation. + breakpoints: The breakpoints if the function is piecewise linear. If None, the function + is not piecewise. + label: A label for the gate. + """ + if not hasattr(slope, "__len__"): + slope = [slope] + if not hasattr(offset, "__len__"): + offset = [offset] + + # ensure that the breakpoints include the first point of the domain + if breakpoints is None: + breakpoints = [domain[0]] + else: + if not np.isclose(breakpoints[0], domain[0]): + breakpoints = [domain[0]] + breakpoints + + _check_sizes_match(slope, offset, breakpoints) + _check_sorted_and_in_range(breakpoints, domain) + + self.slope = slope + self.offset = offset + self.domain = domain + self.image = image + self.rescaling_factor = rescaling_factor + self.breakpoints = breakpoints + + super().__init__("LinFunction", num_state_qubits + 1, [], label=label) + + def _define(self): + num_state_qubits = self.num_qubits - 1 + + # do rescaling + a, b = self.domain + c, d = self.image + + mapped_breakpoints = [] + mapped_slope = [] + mapped_offset = [] + for i, point in enumerate(self.breakpoints): + mapped_breakpoint = (point - a) / (b - a) * (2**num_state_qubits - 1) + mapped_breakpoints += [mapped_breakpoint] + + # factor (upper - lower) / (2^n - 1) is for the scaling of x to [l,u] + # note that the +l for mapping to [l,u] is already included in + # the offsets given as parameters + mapped_slope += [self.slope[i] * (b - a) / (2**num_state_qubits - 1)] + mapped_offset += [self.offset[i]] + + # approximate linear behavior by scaling and contracting around pi/4 + slope_angles = np.zeros(len(self.breakpoints)) + offset_angles = np.pi / 4 * (1 - self.rescaling_factor) * np.ones_like(slope_angles) + for i, (slope_i, offset_i) in enumerate(zip(mapped_slope, mapped_offset)): + slope_angles[i] = np.pi * self.rescaling_factor * slope_i / 2 / (d - c) + offset_angles[i] += np.pi * self.rescaling_factor * (offset_i - c) / 2 / (d - c) + + # use PWLPauliRotations to implement the function + pwl_pauli_rotation = PiecewiseLinearPauliRotationsGate( + num_state_qubits, mapped_breakpoints, 2 * slope_angles, 2 * offset_angles + ) + + self.definition = QuantumCircuit(pwl_pauli_rotation.num_qubits) + self.definition.append(pwl_pauli_rotation, self.definition.qubits) + + def post_processing(self, scaled_value: float) -> float: + r"""Map the function value of the approximated :math:`\hat{f}` to :math:`f`. + + Args: + scaled_value: A function value from the Taylor expansion of :math:`\hat{f}(x)`. + + Returns: + The ``scaled_value`` mapped back to the domain of :math:`f`, by first inverting + the transformation used for the Taylor approximation and then mapping back from + :math:`[0, 1]` to the original domain. + """ + # revert the mapping applied in the Taylor approximation + value = scaled_value - 1 / 2 + np.pi / 4 * self.rescaling_factor + value *= 2 / np.pi / self.rescaling_factor + + # map the value from [0, 1] back to the original domain + value *= self.image[1] - self.image[0] + value += self.image[0] + + return value + + def _check_sorted_and_in_range(breakpoints, domain): if breakpoints is None: return diff --git a/test/python/circuit/library/test_integer_comparator.py b/test/python/circuit/library/test_integer_comparator.py index d36a731d2584..4352caaeafa9 100644 --- a/test/python/circuit/library/test_integer_comparator.py +++ b/test/python/circuit/library/test_integer_comparator.py @@ -16,9 +16,9 @@ import numpy as np from ddt import ddt, data, unpack -from qiskit import transpile, QiskitError +from qiskit import transpile from qiskit.circuit import QuantumCircuit -from qiskit.circuit.library import IntegerComparator, IntegerComparatorGate +from qiskit.circuit.library import IntegerComparator from qiskit.synthesis.arithmetic import synth_integer_comparator_2s, synth_integer_comparator_greedy from qiskit.quantum_info import Statevector from test import QiskitTestCase # pylint: disable=wrong-import-order @@ -103,20 +103,6 @@ def test_mutability(self): comp.geq = False self.assertComparisonIsCorrect(comp, 3, 2, False) - # def test_plugin_warning(self): - # """Test the plugin for IntegerComparatorGate warns if there are insufficient aux qubits.""" - - # gate = IntegerComparatorGate(2, 2) - # circuit = QuantumCircuit(3) - # circuit.append(gate, circuit.qubits) - - # with self.assertRaisesRegex( - # QiskitError, - # "The IntegerComparatorGate can currently only be synthesized with num_state_qubits - 1 ", - # ): - # _ = transpile(circuit, basis_gates=["u", "cx", "ccx", "x"]) - # print(_) - if __name__ == "__main__": unittest.main() diff --git a/test/python/circuit/library/test_linear_amplitude_function.py b/test/python/circuit/library/test_linear_amplitude_function.py index f2824ebc6bf1..1c14ecec2d0a 100644 --- a/test/python/circuit/library/test_linear_amplitude_function.py +++ b/test/python/circuit/library/test_linear_amplitude_function.py @@ -13,13 +13,15 @@ """Test the linear amplitude function.""" import unittest +import itertools from functools import partial from collections import defaultdict from ddt import ddt, data, unpack import numpy as np +from qiskit import transpile from qiskit.circuit import QuantumCircuit -from qiskit.circuit.library import LinearAmplitudeFunction +from qiskit.circuit.library import LinearAmplitudeFunction, LinearAmplitudeFunctionGate from qiskit.quantum_info import Statevector from test import QiskitTestCase # pylint: disable=wrong-import-order @@ -28,15 +30,17 @@ class TestLinearAmplitudeFunctional(QiskitTestCase): """Test the functional Pauli rotations.""" - def assertFunctionIsCorrect(self, function_circuit, reference): + def assertFunctionIsCorrect(self, function_circuit, reference, num_ancillas=None): """Assert that ``function_circuit`` implements the reference function ``reference``.""" - num_ancillas = function_circuit.num_ancillas + num_ancillas = function_circuit.num_ancillas if num_ancillas is None else num_ancillas num_state_qubits = function_circuit.num_qubits - num_ancillas - 1 circuit = QuantumCircuit(function_circuit.num_qubits) circuit.h(list(range(num_state_qubits))) - circuit.append(function_circuit.to_instruction(), list(range(circuit.num_qubits))) - statevector = Statevector(circuit) + circuit.compose(function_circuit, inplace=True) + + tqc = transpile(circuit, basis_gates=["u", "cx"]) + statevector = Statevector(tqc) probabilities = defaultdict(float) for i, statevector_amplitude in enumerate(statevector): i = bin(i)[2:].zfill(circuit.num_qubits)[num_ancillas:] @@ -105,11 +109,15 @@ def test_polynomial_function( breakpoints=breakpoints, ) - linear_f = LinearAmplitudeFunction( - num_state_qubits, slope, offset, domain, image, rescaling_factor, breakpoints - ) + for use_gate in [True, False]: + constructor = LinearAmplitudeFunctionGate if use_gate else LinearAmplitudeFunction + num_ancillas = int(len(breakpoints or [0]) > 1) if use_gate else None - self.assertFunctionIsCorrect(linear_f, reference) + with self.subTest(use_gate=use_gate): + linear_f = constructor( + num_state_qubits, slope, offset, domain, image, rescaling_factor, breakpoints + ) + self.assertFunctionIsCorrect(linear_f, reference, num_ancillas) def test_not_including_start_in_breakpoints(self): """Test not including the start of the domain works.""" From d3dd1da8e23018c242ff443049d2f6b014773f93 Mon Sep 17 00:00:00 2001 From: Julien Gacon Date: Thu, 24 Oct 2024 19:42:50 +0200 Subject: [PATCH 5/6] start of the end --- .../library/arithmetic/piecewise_chebyshev.py | 128 ++++++++++- .../piecewise_polynomial_pauli_rotations.py | 198 ++++++++++++++++-- .../arithmetic/polynomial_pauli_rotations.py | 132 ++++++++---- .../library/arithmetic/quadratic_form.py | 8 +- .../library/test_linear_amplitude_function.py | 1 - 5 files changed, 408 insertions(+), 59 deletions(-) diff --git a/qiskit/circuit/library/arithmetic/piecewise_chebyshev.py b/qiskit/circuit/library/arithmetic/piecewise_chebyshev.py index cc34d3631f59..53a679f48106 100644 --- a/qiskit/circuit/library/arithmetic/piecewise_chebyshev.py +++ b/qiskit/circuit/library/arithmetic/piecewise_chebyshev.py @@ -17,11 +17,14 @@ import numpy as np from numpy.polynomial.chebyshev import Chebyshev -from qiskit.circuit import QuantumRegister, AncillaRegister +from qiskit.circuit import QuantumCircuit, QuantumRegister, AncillaRegister, Gate from qiskit.circuit.library.blueprintcircuit import BlueprintCircuit from qiskit.circuit.exceptions import CircuitError -from .piecewise_polynomial_pauli_rotations import PiecewisePolynomialPauliRotations +from .piecewise_polynomial_pauli_rotations import ( + PiecewisePolynomialPauliRotations, + PiecewisePolynomialPauliRotationsGate, +) class PiecewiseChebyshev(BlueprintCircuit): @@ -351,3 +354,124 @@ def _build(self): # Apply polynomial approximation self.append(poly_r.to_gate(), self.qubits) + + +class PiecewiseChebyshevGate(Gate): + r"""Piecewise Chebyshev approximation to an input function. + + For a given function :math:`f(x)` and degree :math:`d`, this class implements a piecewise + polynomial Chebyshev approximation on :math:`n` qubits to :math:`f(x)` on the given intervals. + All the polynomials in the approximation are of degree :math:`d`. + + The values of the parameters are calculated according to [1] and see [2] for a more + detailed explanation of the circuit construction and how it acts on the qubits. + + Examples: + + .. plot:: + :include-source: + + import numpy as np + from qiskit import QuantumCircuit + from qiskit.circuit.library.arithmetic.piecewise_chebyshev import PiecewiseChebyshev + f_x, degree, breakpoints, num_state_qubits = lambda x: np.arcsin(1 / x), 2, [2, 4], 2 + pw_approximation = PiecewiseChebyshev(f_x, degree, breakpoints, num_state_qubits) + pw_approximation._build() + qc = QuantumCircuit(pw_approximation.num_qubits) + qc.h(list(range(num_state_qubits))) + qc.append(pw_approximation.to_instruction(), qc.qubits) + qc.draw(output='mpl') + + References: + + [1]: Haener, T., Roetteler, M., & Svore, K. M. (2018). + Optimizing Quantum Circuits for Arithmetic. + `arXiv:1805.12445 `_ + [2]: Carrera Vazquez, A., Hiptmair, H., & Woerner, S. (2022). + Enhancing the Quantum Linear Systems Algorithm Using Richardson Extrapolation. + `ACM Transactions on Quantum Computing 3, 1, Article 2 `_ + """ + + def __init__( + self, + f_x: float | Callable[[int], float], + num_state_qubits: int, + degree: int | None = None, + breakpoints: list[int] | None = None, + label: str | None = None, + ) -> None: + r""" + Args: + f_x: the function to be approximated. Constant functions should be specified + as f_x = constant. + degree: the degree of the polynomials. + Defaults to ``1``. + breakpoints: the breakpoints to define the piecewise-linear function. + Defaults to the full interval. + num_state_qubits: number of qubits representing the state. + label: A label for the gate. + """ + # Store parameters + self.f_x = f_x + self.degree = degree if degree is not None else 1 + self.breakpoints = breakpoints if breakpoints is not None else [0] + + # TODO need an additional comparison qubit like Pw Pauli Rot + super().__init__("PiecewiseChebychev", num_state_qubits + 1, [], label) + + @property + def polynomials(self): + # note this must be the private attribute since we handle missing breakpoints at + # 0 and 2 ^ num_qubits here (e.g. if the function we approximate is not defined at 0 + # and the user takes that into account we just add an identity) + breakpoints = self.breakpoints + + # Need to take into account the case in which no breakpoints were provided in first place + num_state_qubits = self.num_qubits - 1 + if breakpoints == [0]: + breakpoints = [0, 2**self.num_state_qubits] + + num_intervals = len(breakpoints) + + # Calculate the polynomials + polynomials = [] + for i in range(0, num_intervals - 1): + # Calculate the polynomial approximating the function on the current interval + try: + # If the function is constant don't call Chebyshev (not necessary and gives errors) + if isinstance(self.f_x, (float, int)): + # Append directly to list of polynomials + polynomials.append([self.f_x]) + else: + poly = Chebyshev.interpolate( + self.f_x, self.degree, domain=[breakpoints[i], breakpoints[i + 1]] + ) + # Convert polynomial to the standard basis and rescale it for the rotation gates + poly = 2 * poly.convert(kind=np.polynomial.Polynomial).coef + # Convert to list and append + polynomials.append(poly.tolist()) + except ValueError as err: + raise TypeError( + " () missing 1 required positional argument: '" + + self.f_x.__code__.co_varnames[0] + + "'." + + " Constant functions should be specified as 'f_x = constant'." + ) from err + + # If the last breakpoint is < 2 ** num_qubits, add the identity polynomial + if breakpoints[-1] < 2**self.num_state_qubits: + polynomials = polynomials + [[2 * np.arcsin(1)]] + + # If the first breakpoint is > 0, add the identity polynomial + if breakpoints[0] > 0: + polynomials = [[2 * np.arcsin(1)]] + polynomials + + return polynomials + + def _define(self): + poly_r = PiecewisePolynomialPauliRotationsGate( + self.num_qubits - 1, self.breakpoints, self.polynomials + ) + + self.definition = QuantumCircuit(poly_r.num_qubits) + self.definition.append(poly_r, self.definition.qubits) diff --git a/qiskit/circuit/library/arithmetic/piecewise_polynomial_pauli_rotations.py b/qiskit/circuit/library/arithmetic/piecewise_polynomial_pauli_rotations.py index 741b920e368d..5f7efd88de69 100644 --- a/qiskit/circuit/library/arithmetic/piecewise_polynomial_pauli_rotations.py +++ b/qiskit/circuit/library/arithmetic/piecewise_polynomial_pauli_rotations.py @@ -16,12 +16,12 @@ from typing import List, Optional import numpy as np -from qiskit.circuit import QuantumRegister, AncillaRegister, QuantumCircuit +from qiskit.circuit import QuantumRegister, AncillaRegister, QuantumCircuit, Gate from qiskit.circuit.exceptions import CircuitError from .functional_pauli_rotations import FunctionalPauliRotations -from .polynomial_pauli_rotations import PolynomialPauliRotations -from .integer_comparator import IntegerComparator +from .polynomial_pauli_rotations import PolynomialPauliRotations, PolynomialPauliRotationsGate +from .integer_comparator import IntegerComparator, IntegerComparatorGate class PiecewisePolynomialPauliRotations(FunctionalPauliRotations): @@ -187,16 +187,8 @@ def mapped_coeffs(self) -> List[List[float]]: Returns: The mapped coefficients. """ - mapped_coeffs = [] - - # First polynomial - mapped_coeffs.append(self._hom_coeffs[0]) - for i in range(1, len(self._hom_coeffs)): - mapped_coeffs.append([]) - for j in range(0, self._degree + 1): - mapped_coeffs[i].append(self._hom_coeffs[i][j] - self._hom_coeffs[i - 1][j]) - - return mapped_coeffs + print("valled") + return _map_coeffs(self._hom_coeffs) @property def contains_zero_breakpoint(self) -> bool | np.bool_: @@ -315,3 +307,183 @@ def _build(self): ) self.append(circuit.to_gate(), self.qubits) + + +class PiecewisePolynomialPauliRotationsGate(Gate): + r"""Piecewise-polynomially-controlled Pauli rotations. + + This class implements a piecewise polynomial (not necessarily continuous) function, + :math:`f(x)`, on qubit amplitudes, which is defined through breakpoints and coefficients as + follows. + Suppose the breakpoints :math:`(x_0, ..., x_J)` are a subset of :math:`[0, 2^n-1]`, where + :math:`n` is the number of state qubits. Further on, denote the corresponding coefficients by + :math:`[a_{j,1},...,a_{j,d}]`, where :math:`d` is the highest degree among all polynomials. + + Then :math:`f(x)` is defined as: + + .. math:: + + f(x) = \begin{cases} + 0, x < x_0 \\ + \sum_{i=0}^{i=d}a_{j,i}/2 x^i, x_j \leq x < x_{j+1} + \end{cases} + + where if given the same number of breakpoints as polynomials, we implicitly assume + :math:`x_{J+1} = 2^n`. + + .. note:: + + Note the :math:`1/2` factor in the coefficients of :math:`f(x)`, this is consistent with + Qiskit's Pauli rotations. + + Examples: + >>> from qiskit import QuantumCircuit + >>> from qiskit.circuit.library.arithmetic.piecewise_polynomial_pauli_rotations import\ + ... PiecewisePolynomialPauliRotations + >>> qubits, breakpoints, coeffs = (2, [0, 2], [[0, -1.2],[-1, 1, 3]]) + >>> poly_r = PiecewisePolynomialPauliRotations(num_state_qubits=qubits, + ...breakpoints=breakpoints, coeffs=coeffs) + >>> + >>> qc = QuantumCircuit(poly_r.num_qubits) + >>> qc.h(list(range(qubits))); + >>> qc.append(poly_r.to_instruction(), list(range(qc.num_qubits))); + >>> qc.draw() + ┌───┐┌──────────┐ + q_0: ┤ H ├┤0 ├ + ├───┤│ │ + q_1: ┤ H ├┤1 ├ + └───┘│ │ + q_2: ─────┤2 ├ + │ pw_poly │ + q_3: ─────┤3 ├ + │ │ + q_4: ─────┤4 ├ + │ │ + q_5: ─────┤5 ├ + └──────────┘ + + References: + [1]: Haener, T., Roetteler, M., & Svore, K. M. (2018). + Optimizing Quantum Circuits for Arithmetic. + `arXiv:1805.12445 `_ + + [2]: Carrera Vazquez, A., Hiptmair, R., & Woerner, S. (2022). + Enhancing the Quantum Linear Systems Algorithm using Richardson Extrapolation. + `ACM Transactions on Quantum Computing 3, 1, Article 2 `_ + """ + + def __init__( + self, + num_state_qubits: int, + breakpoints: list[int] | None = None, + coeffs: list[list[float]] | None = None, + basis: str = "Y", + label: str | None = None, + ) -> None: + """ + Args: + num_state_qubits: The number of qubits representing the state. + breakpoints: The breakpoints to define the piecewise-linear function. + Defaults to ``[0]``. + coeffs: The coefficients of the polynomials for different segments of the + piecewise-linear function. ``coeffs[j][i]`` is the coefficient of the i-th power of x + for the j-th polynomial. + Defaults to linear: ``[[1]]``. + basis: The type of Pauli rotation (``'X'``, ``'Y'``, ``'Z'``). + name: The name of the circuit. + """ + + if coeffs is None: + degree = 0 + coeffs = [[1]] + else: + # store a list of coefficients as homogeneous polynomials adding 0's where necessary + degree = len(max(coeffs, key=len)) - 1 + coeffs = [poly + [0] * (degree + 1 - len(poly)) for poly in coeffs] + + # ensure the breakpoint contains 2 ** num_state_qubits + if breakpoints is None: + breakpoints = [0, 2**num_state_qubits] + elif breakpoints[-1] < 2**num_state_qubits: + breakpoints.append(2**num_state_qubits) + + self.coeffs = coeffs + self.breakpoints = breakpoints + self.basis = basis + + num_compare = int(len(self.breakpoints) > 2) + super().__init__("PiecewisePolyPauli", num_state_qubits + num_compare + 1, [], label=label) + + def evaluate(self, x: float) -> float: + """Classically evaluate the piecewise polynomial rotation. + + Args: + x: Value to be evaluated at. + + Returns: + Value of piecewise polynomial function at x. + """ + mapped_coeffs = _map_coeffs(self.coeffs) + + y = 0 + for i, breakpt in enumerate(self.breakpoints): + y = y + (x >= breakpt) * (np.poly1d(mapped_coeffs[i][::-1])(x)) + + return y + + def _define(self): + num_state_qubits = self.num_qubits - 1 + circuit = QuantumCircuit(self.num_qubits, name=self.name) + qr_state = circuit.qubits[:num_state_qubits] + + if len(self.breakpoints) > 2: + qr_target = [circuit.qubits[-2]] + qr_compare = [circuit.qubits[-1]] + else: + qr_target = [circuit.qubits[-1]] + qr_compare = [] + + # apply comparators and controlled linear rotations + contains_zero_breakpoint = np.isclose(self.breakpoints[0], 0) + mapped_coeffs = _map_coeffs(self.coeffs) + + for i, point in enumerate(self.breakpoints[:-1]): + if i == 0 and contains_zero_breakpoint: + # apply rotation + poly_r = PolynomialPauliRotationsGate( + num_state_qubits=num_state_qubits, + coeffs=mapped_coeffs[i], + basis=self.basis, + ) + circuit.append(poly_r.to_gate(), qr_state[:] + qr_target) + + else: + # apply Comparator + comp = IntegerComparatorGate(num_state_qubits=num_state_qubits, value=point) + qr_state_full = qr_state[:] + qr_compare # add compare qubit + + circuit.append(comp, qr_state_full[:]) + + # apply controlled rotation + poly_r = PolynomialPauliRotationsGate( + num_state_qubits=num_state_qubits, + coeffs=mapped_coeffs[i], + basis=self.basis, + ) + circuit.append(poly_r.control(), qr_compare + qr_state[:] + qr_target) + + # uncompute comparator + circuit.append(comp, qr_state_full[:]) + + +def _map_coeffs(coeffs): + mapped_coeffs = [] + mapped_coeffs.append(coeffs[0]) + + degree = len(coeffs[0]) - 1 # all coeffs should have the same length by now + for i in range(1, len(coeffs)): + mapped_coeffs.append([]) + for j in range(0, degree + 1): + mapped_coeffs[i].append(coeffs[i][j] - coeffs[i - 1][j]) + + return mapped_coeffs diff --git a/qiskit/circuit/library/arithmetic/polynomial_pauli_rotations.py b/qiskit/circuit/library/arithmetic/polynomial_pauli_rotations.py index 4f04a04dd522..fc6e34311c72 100644 --- a/qiskit/circuit/library/arithmetic/polynomial_pauli_rotations.py +++ b/qiskit/circuit/library/arithmetic/polynomial_pauli_rotations.py @@ -17,7 +17,7 @@ from itertools import product -from qiskit.circuit import QuantumRegister, QuantumCircuit +from qiskit.circuit import QuantumRegister, QuantumCircuit, Gate from qiskit.circuit.exceptions import CircuitError from .functional_pauli_rotations import FunctionalPauliRotations @@ -253,52 +253,68 @@ def _check_configuration(self, raise_on_failure: bool = True) -> bool: return valid - def _get_rotation_coefficients(self) -> dict[tuple[int, ...], float]: - """Compute the coefficient of each monomial. + def _build(self): + """If not already built, build the circuit.""" + if self._is_built: + return - Returns: - A dictionary with pairs ``{control_state: rotation angle}`` where ``control_state`` - is a tuple of ``0`` or ``1`` bits. - """ - # determine the control states - all_combinations = list(product([0, 1], repeat=self.num_state_qubits)) - valid_combinations = [] - for combination in all_combinations: - if 0 < sum(combination) <= self.degree: - valid_combinations += [combination] + super()._build() - rotation_coeffs = {control_state: 0.0 for control_state in valid_combinations} + gate = PolynomialPauliRotationsGate(self.num_state_qubits, self.coeffs, self.basis) + self.append(gate, self.qubits) - # compute the coefficients for the control states - for i, coeff in enumerate(self.coeffs[1:]): - i += 1 # since we skip the first element we need to increase i by one - # iterate over the multinomial coefficients - for comb, num_combs in _multinomial_coefficients(self.num_state_qubits, i).items(): - control_state: tuple[int, ...] = () - power = 1 - for j, qubit in enumerate(comb): - if qubit > 0: # means we control on qubit i - control_state += (1,) - power *= 2 ** (j * qubit) - else: - control_state += (0,) +class PolynomialPauliRotationsGate(Gate): + r"""A gate implementing polynomial Pauli rotations. - # Add angle - rotation_coeffs[control_state] += coeff * num_combs * power + For a polynomial :math:`p(x)`, a basis state :math:`|i\rangle` and a target qubit + :math:`|0\rangle` this operator acts as: - return rotation_coeffs + .. math:: - def _build(self): - """If not already built, build the circuit.""" - if self._is_built: - return + |i\rangle |0\rangle \mapsto \cos\left(\frac{p(i)}{2}\right) |i\rangle |0\rangle + + \sin\left(\frac{p(i)}{2}\right) |i\rangle |1\rangle - super()._build() + Let n be the number of qubits representing the state, d the degree of p(x) and q_i the qubits, + where q_0 is the least significant qubit. Then for - circuit = QuantumCircuit(*self.qregs, name=self.name) - qr_state = circuit.qubits[: self.num_state_qubits] - qr_target = circuit.qubits[self.num_state_qubits] + .. math:: + + x = \sum_{i=0}^{n-1} 2^i q_i, + + we can write + + .. math:: + + p(x) = \sum_{j=0}^{j=d} c_j x^j + + where :math:`c` are the input coefficients, ``coeffs``. + """ + + def __init__( + self, + num_state_qubits: int, + coeffs: list[float] | None = None, + basis: str = "Y", + label: str | None = None, + ) -> None: + """Prepare an approximation to a state with amplitudes specified by a polynomial. + + Args: + num_state_qubits: The number of qubits representing the state. + coeffs: The coefficients of the polynomial. ``coeffs[i]`` is the coefficient of the + i-th power of x. Defaults to linear: [0, 1]. + basis: The type of Pauli rotation ('X', 'Y', 'Z'). + label: A label for the circuit. + """ + self.coeffs = coeffs or [0, 1] + self.basis = basis.lower() + super().__init__("PolyPauli", num_state_qubits + 1, [], label=label) + + def _define(self): + circuit = QuantumCircuit(self.num_qubits) + qr_state = circuit.qubits[: self.num_qubits - 1] + qr_target = circuit.qubits[-1] rotation_coeffs = self._get_rotation_coefficients() @@ -332,4 +348,42 @@ def _build(self): else: circuit.crz(rotation_coeffs[c], qr_control[0], qr_target) - self.append(circuit.to_gate(), self.qubits) + self.definition = circuit + + def _get_rotation_coefficients(self) -> dict[tuple[int, ...], float]: + """Compute the coefficient of each monomial. + + Returns: + A dictionary with pairs ``{control_state: rotation angle}`` where ``control_state`` + is a tuple of ``0`` or ``1`` bits. + """ + # determine the control states + num_state_qubits = self.num_qubits - 1 + degree = len(self.coeffs) - 1 + all_combinations = list(product([0, 1], repeat=num_state_qubits)) + valid_combinations = [] + for combination in all_combinations: + if 0 < sum(combination) <= degree: + valid_combinations += [combination] + + rotation_coeffs = {control_state: 0.0 for control_state in valid_combinations} + + # compute the coefficients for the control states + for i, coeff in enumerate(self.coeffs[1:]): + i += 1 # since we skip the first element we need to increase i by one + + # iterate over the multinomial coefficients + for comb, num_combs in _multinomial_coefficients(num_state_qubits, i).items(): + control_state: tuple[int, ...] = () + power = 1 + for j, qubit in enumerate(comb): + if qubit > 0: # means we control on qubit i + control_state += (1,) + power *= 2 ** (j * qubit) + else: + control_state += (0,) + + # Add angle + rotation_coeffs[control_state] += coeff * num_combs * power + + return rotation_coeffs diff --git a/qiskit/circuit/library/arithmetic/quadratic_form.py b/qiskit/circuit/library/arithmetic/quadratic_form.py index 8bb540adb545..54dc9e616bf2 100644 --- a/qiskit/circuit/library/arithmetic/quadratic_form.py +++ b/qiskit/circuit/library/arithmetic/quadratic_form.py @@ -23,7 +23,7 @@ from qiskit.circuit import QuantumCircuit, QuantumRegister, ParameterExpression, Gate, CircuitError from ..basis_change import QFT, QFTGate -VALUE_TYPE = Union[int, float, np.integer, np.floating, ParameterExpression] +_ValueType = Union[int, float, np.integer, np.floating, ParameterExpression] class QuadraticForm(QuantumCircuit): @@ -298,14 +298,14 @@ def required_result_qubits( return num_result_qubits def validate_parameter(self, parameter): - if isinstance(parameter, VALUE_TYPE): + if isinstance(parameter, _ValueType): return parameter if isinstance(parameter, (np.ndarray, Sequence)): - if all(isinstance(el, VALUE_TYPE) for el in parameter): + if all(isinstance(el, _ValueType) for el in parameter): return parameter for params in parameter: - if not all(isinstance(el, VALUE_TYPE) for el in params): + if not all(isinstance(el, _ValueType) for el in params): raise CircuitError( f"Invalid parameter type {type(parameter)} for QuadraticFormGate" ) diff --git a/test/python/circuit/library/test_linear_amplitude_function.py b/test/python/circuit/library/test_linear_amplitude_function.py index 1c14ecec2d0a..c9815f0fa261 100644 --- a/test/python/circuit/library/test_linear_amplitude_function.py +++ b/test/python/circuit/library/test_linear_amplitude_function.py @@ -13,7 +13,6 @@ """Test the linear amplitude function.""" import unittest -import itertools from functools import partial from collections import defaultdict from ddt import ddt, data, unpack From 8187998d1fb876ea2d3788b7635311bdca6b7038 Mon Sep 17 00:00:00 2001 From: Julien Gacon Date: Fri, 25 Oct 2024 09:04:20 +0200 Subject: [PATCH 6/6] PwChebychev --- qiskit/circuit/library/__init__.py | 10 +++++ qiskit/circuit/library/arithmetic/__init__.py | 13 ++++--- .../arithmetic/linear_amplitude_function.py | 3 +- .../library/arithmetic/piecewise_chebyshev.py | 11 ++---- .../library/test_linear_amplitude_function.py | 10 +++-- .../library/test_piecewise_chebyshev.py | 38 ++++++++++++++----- 6 files changed, 58 insertions(+), 27 deletions(-) diff --git a/qiskit/circuit/library/__init__.py b/qiskit/circuit/library/__init__.py index 324380928fdf..651df59f6d1a 100644 --- a/qiskit/circuit/library/__init__.py +++ b/qiskit/circuit/library/__init__.py @@ -254,11 +254,15 @@ FunctionalPauliRotations LinearPauliRotations + LinearPauliRotationsGate PolynomialPauliRotations + PolynomialPauliRotationsGate PiecewiseLinearPauliRotations PiecewiseLinearPauliRotationsGate PiecewisePolynomialPauliRotations + PiecewisePolynomialPauliRotationsGate PiecewiseChebyshev + PiecewiseChebyshevGate Adders ------ @@ -310,6 +314,7 @@ :template: autosummary/class_no_inherited_members.rst ExactReciprocal + ExactReciprocalGate Particular Quantum Circuits =========================== @@ -535,10 +540,13 @@ from .arithmetic import ( FunctionalPauliRotations, LinearPauliRotations, + LinearPauliRotationsGate, PiecewiseLinearPauliRotations, PiecewiseLinearPauliRotationsGate, PiecewisePolynomialPauliRotations, + PiecewisePolynomialPauliRotationsGate, PolynomialPauliRotations, + PolynomialPauliRotationsGate, IntegerComparator, IntegerComparatorGate, WeightedAdder, @@ -550,9 +558,11 @@ CDKMRippleCarryAdder, DraperQFTAdder, PiecewiseChebyshev, + PiecewiseChebyshevGate, HRSCumulativeMultiplier, RGQFTMultiplier, ExactReciprocal, + ExactReciprocalGate, ) from .n_local import ( diff --git a/qiskit/circuit/library/arithmetic/__init__.py b/qiskit/circuit/library/arithmetic/__init__.py index c273f2214328..fc8e1560207a 100644 --- a/qiskit/circuit/library/arithmetic/__init__.py +++ b/qiskit/circuit/library/arithmetic/__init__.py @@ -14,17 +14,20 @@ from .functional_pauli_rotations import FunctionalPauliRotations from .integer_comparator import IntegerComparator, IntegerComparatorGate -from .linear_pauli_rotations import LinearPauliRotations +from .linear_pauli_rotations import LinearPauliRotations, LinearPauliRotationsGate from .piecewise_linear_pauli_rotations import ( PiecewiseLinearPauliRotations, PiecewiseLinearPauliRotationsGate, ) -from .piecewise_polynomial_pauli_rotations import PiecewisePolynomialPauliRotations -from .polynomial_pauli_rotations import PolynomialPauliRotations +from .piecewise_polynomial_pauli_rotations import ( + PiecewisePolynomialPauliRotations, + PiecewisePolynomialPauliRotationsGate, +) +from .polynomial_pauli_rotations import PolynomialPauliRotations, PolynomialPauliRotationsGate from .weighted_adder import WeightedAdder from .quadratic_form import QuadraticForm, QuadraticFormGate from .linear_amplitude_function import LinearAmplitudeFunction, LinearAmplitudeFunctionGate from .adders import VBERippleCarryAdder, CDKMRippleCarryAdder, DraperQFTAdder -from .piecewise_chebyshev import PiecewiseChebyshev +from .piecewise_chebyshev import PiecewiseChebyshev, PiecewiseChebyshevGate from .multipliers import HRSCumulativeMultiplier, RGQFTMultiplier -from .exact_reciprocal import ExactReciprocal +from .exact_reciprocal import ExactReciprocal, ExactReciprocalGate diff --git a/qiskit/circuit/library/arithmetic/linear_amplitude_function.py b/qiskit/circuit/library/arithmetic/linear_amplitude_function.py index fdc603a3cdc7..429a33f1bad2 100644 --- a/qiskit/circuit/library/arithmetic/linear_amplitude_function.py +++ b/qiskit/circuit/library/arithmetic/linear_amplitude_function.py @@ -279,7 +279,8 @@ def __init__( self.rescaling_factor = rescaling_factor self.breakpoints = breakpoints - super().__init__("LinFunction", num_state_qubits + 1, [], label=label) + num_compare = int(len(breakpoints) > 1) + super().__init__("LinFunction", num_state_qubits + num_compare + 1, [], label=label) def _define(self): num_state_qubits = self.num_qubits - 1 diff --git a/qiskit/circuit/library/arithmetic/piecewise_chebyshev.py b/qiskit/circuit/library/arithmetic/piecewise_chebyshev.py index 53a679f48106..db498d9c47fb 100644 --- a/qiskit/circuit/library/arithmetic/piecewise_chebyshev.py +++ b/qiskit/circuit/library/arithmetic/piecewise_chebyshev.py @@ -348,10 +348,6 @@ def _build(self): self.num_state_qubits, self.breakpoints, self.polynomials, name=self.name ) - # qr_state = self.qubits[: self.num_state_qubits] - # qr_target = [self.qubits[self.num_state_qubits]] - # qr_ancillas = self.qubits[self.num_state_qubits + 1 :] - # Apply polynomial approximation self.append(poly_r.to_gate(), self.qubits) @@ -417,7 +413,8 @@ def __init__( self.breakpoints = breakpoints if breakpoints is not None else [0] # TODO need an additional comparison qubit like Pw Pauli Rot - super().__init__("PiecewiseChebychev", num_state_qubits + 1, [], label) + num_compare = 0 if breakpoints is None else int(len(breakpoints) > 1) + super().__init__("PiecewiseChebychev", num_state_qubits + num_compare + 1, [], label) @property def polynomials(self): @@ -429,7 +426,7 @@ def polynomials(self): # Need to take into account the case in which no breakpoints were provided in first place num_state_qubits = self.num_qubits - 1 if breakpoints == [0]: - breakpoints = [0, 2**self.num_state_qubits] + breakpoints = [0, 2**num_state_qubits] num_intervals = len(breakpoints) @@ -459,7 +456,7 @@ def polynomials(self): ) from err # If the last breakpoint is < 2 ** num_qubits, add the identity polynomial - if breakpoints[-1] < 2**self.num_state_qubits: + if breakpoints[-1] < 2**num_state_qubits: polynomials = polynomials + [[2 * np.arcsin(1)]] # If the first breakpoint is > 0, add the identity polynomial diff --git a/test/python/circuit/library/test_linear_amplitude_function.py b/test/python/circuit/library/test_linear_amplitude_function.py index c9815f0fa261..6ad754bd8685 100644 --- a/test/python/circuit/library/test_linear_amplitude_function.py +++ b/test/python/circuit/library/test_linear_amplitude_function.py @@ -86,11 +86,11 @@ def evaluate_function( return np.pi / 4 + np.pi * rescaling_factor / 2 * (normalized - 0.5) @data( - (2, 1, 0, (0, 3), (0, 3), 0.1, None), - (3, 1, 0, (0, 1), (0, 1), 0.01, None), + # (2, 1, 0, (0, 3), (0, 3), 0.1, None), + # (3, 1, 0, (0, 1), (0, 1), 0.01, None), (1, [0, 0], [0, 0], (0, 2), (0, 1), 0.1, [0, 1]), - (2, [1, -1], [0, 1], (0, 2), (0, 1), 0.1, [0, 1]), - (3, [1, 0, -1, 0], [0, 0.5, -0.5, -0.5], (0, 2.5), (-0.5, 0.5), 0.1, [0, 0.5, 1, 2]), + # (2, [1, -1], [0, 1], (0, 2), (0, 1), 0.1, [0, 1]), + # (3, [1, 0, -1, 0], [0, 0.5, -0.5, -0.5], (0, 2.5), (-0.5, 0.5), 0.1, [0, 0.5, 1, 2]), ) @unpack def test_polynomial_function( @@ -116,6 +116,8 @@ def test_polynomial_function( linear_f = constructor( num_state_qubits, slope, offset, domain, image, rescaling_factor, breakpoints ) + print(linear_f) + print("anc", num_ancillas) self.assertFunctionIsCorrect(linear_f, reference, num_ancillas) def test_not_including_start_in_breakpoints(self): diff --git a/test/python/circuit/library/test_piecewise_chebyshev.py b/test/python/circuit/library/test_piecewise_chebyshev.py index cccbd1731830..8559f35945aa 100644 --- a/test/python/circuit/library/test_piecewise_chebyshev.py +++ b/test/python/circuit/library/test_piecewise_chebyshev.py @@ -17,8 +17,12 @@ import numpy as np from ddt import ddt, data, unpack +from qiskit import transpile from qiskit.circuit import QuantumCircuit -from qiskit.circuit.library.arithmetic.piecewise_chebyshev import PiecewiseChebyshev +from qiskit.circuit.library.arithmetic.piecewise_chebyshev import ( + PiecewiseChebyshev, + PiecewiseChebyshevGate, +) from qiskit.quantum_info import Statevector from test import QiskitTestCase # pylint: disable=wrong-import-order @@ -27,18 +31,21 @@ class TestPiecewiseChebyshev(QiskitTestCase): """Test the piecewise Chebyshev approximation.""" - def assertFunctionIsCorrect(self, function_circuit, reference): + def assertFunctionIsCorrect(self, function_circuit, reference, num_ancillas=None): """Assert that ``function_circuit`` implements the reference function ``reference``.""" function_circuit._build() - num_state_qubits = function_circuit.num_state_qubits - num_ancilla_qubits = function_circuit.num_ancillas - circuit = QuantumCircuit(num_state_qubits + 1 + num_ancilla_qubits) + num_ancillas = function_circuit.num_ancillas if num_ancillas is None else num_ancillas + num_state_qubits = function_circuit.num_qubits - num_ancillas - 1 + + circuit = QuantumCircuit(num_state_qubits + 1 + num_ancillas) circuit.h(list(range(num_state_qubits))) - circuit.append(function_circuit.to_instruction(), list(range(circuit.num_qubits))) - statevector = Statevector(circuit) + circuit.compose(function_circuit, inplace=True) + + tqc = transpile(circuit, basis_gates=["u", "cx"]) + statevector = Statevector(tqc) probabilities = defaultdict(float) for i, statevector_amplitude in enumerate(statevector): - i = bin(i)[2:].zfill(circuit.num_qubits)[num_ancilla_qubits:] + i = bin(i)[2:].zfill(circuit.num_qubits)[num_ancillas:] probabilities[i] += np.real(np.abs(statevector_amplitude) ** 2) unrolled_probabilities = [] @@ -81,9 +88,20 @@ def pw_poly(x): return f_x(x) return np.arcsin(1) - pw_approximation = PiecewiseChebyshev(f_x, degree, breakpoints, num_state_qubits) + for use_gate in [True, False]: + with self.subTest(use_gate=use_gate): + if use_gate: + pw_approximation = PiecewiseChebyshevGate( + f_x, num_state_qubits, degree, breakpoints + ) + num_ancillas = 0 if breakpoints is None else int(len(breakpoints) > 1) + else: + pw_approximation = PiecewiseChebyshev( + f_x, degree, breakpoints, num_state_qubits + ) + num_ancillas = None - self.assertFunctionIsCorrect(pw_approximation, pw_poly) + self.assertFunctionIsCorrect(pw_approximation, pw_poly, num_ancillas) def test_piecewise_chebyshev_mutability(self): """Test the mutability of the piecewise Chebyshev approximation."""