diff --git a/qualtran/bloqs/arithmetic/addition.py b/qualtran/bloqs/arithmetic/addition.py index 882722ec7..47e2e1fd1 100644 --- a/qualtran/bloqs/arithmetic/addition.py +++ b/qualtran/bloqs/arithmetic/addition.py @@ -129,8 +129,19 @@ def on_classical_vals( ) -> Dict[str, 'ClassicalValT']: unsigned = isinstance(self.a_dtype, (QUInt, QMontgomeryUInt)) b_bitsize = self.b_dtype.bitsize - N = 2**b_bitsize if unsigned else 2 ** (b_bitsize - 1) - return {'a': a, 'b': int(math.fmod(a + b, N))} + N = 2**b_bitsize + if unsigned: + return {'a': a, 'b': int((a + b) % N)} + + # Addition of signed integers can result in overflow. In most classical programming languages (e.g. C++) + # what happens when an overflow happens is left as an implementation detail for compiler designers. + # However for quantum subtraction the operation should be unitary and that means that the unitary of + # the bloq should be a permutation matrix. + # If we hold `a` constant then the valid range of values of `b` [-N/2, N/2) gets shifted forward or backwards + # by `a`. to keep the operation unitary overflowing values wrap around. this is the same as moving the range [0, N) + # by the same amount modulu $N$. that is add N/2 before addition and then remove it. + half_n = N >> 1 + return {'a': a, 'b': int(a + b + half_n) % N - half_n} def _circuit_diagram_info_(self, _) -> cirq.CircuitDiagramInfo: wire_symbols = ["In(x)"] * int(self.a_dtype.bitsize) @@ -188,6 +199,9 @@ def decompose_from_registers( # reverse the order of qubits for big endian-ness. input_bits = quregs['a'][::-1] output_bits = quregs['b'][::-1] + if self.b_dtype.bitsize == 1: + yield CNOT().on(input_bits[0], output_bits[0]) + return ancillas = context.qubit_manager.qalloc(self.b_dtype.bitsize - 1)[::-1] # Start off the addition by anding into the ancilla yield And().on(input_bits[0], output_bits[0], ancillas[0]) diff --git a/qualtran/bloqs/arithmetic/addition_test.py b/qualtran/bloqs/arithmetic/addition_test.py index 858ea21f2..1337dccc2 100644 --- a/qualtran/bloqs/arithmetic/addition_test.py +++ b/qualtran/bloqs/arithmetic/addition_test.py @@ -316,6 +316,13 @@ def test_classical_add_k_unsigned(bitsize, k, x, cvs, ctrls, result): assert bloq_classical[-1] == result +@pytest.mark.parametrize('bitsize', range(2, 5)) +def test_classical_add_signed_overflow(bitsize): + bloq = Add(QInt(bitsize)) + mx = 2 ** (bitsize - 1) - 1 + assert bloq.call_classically(a=mx, b=mx) == (mx, -2) + + # TODO: write tests for signed integer addition (subtraction) # https://github.com/quantumlib/Qualtran/issues/606 @pytest.mark.parametrize('bitsize,k,x,cvs,ctrls,result', [(5, 2, 0, (1, 0), (1, 0), 2)]) diff --git a/qualtran/bloqs/arithmetic/subtraction.ipynb b/qualtran/bloqs/arithmetic/subtraction.ipynb index d7e5c66f3..584a0dce4 100644 --- a/qualtran/bloqs/arithmetic/subtraction.ipynb +++ b/qualtran/bloqs/arithmetic/subtraction.ipynb @@ -2,7 +2,7 @@ "cells": [ { "cell_type": "markdown", - "id": "5e1ac363", + "id": "3f62c768", "metadata": { "cq.autogen": "title_cell" }, @@ -13,7 +13,7 @@ { "cell_type": "code", "execution_count": null, - "id": "3a2de513", + "id": "1e50eeb1", "metadata": { "cq.autogen": "top_imports" }, @@ -30,7 +30,7 @@ }, { "cell_type": "markdown", - "id": "b91ab938", + "id": "9c974476", "metadata": { "cq.autogen": "Subtract.bloq_doc.md" }, @@ -38,9 +38,11 @@ "## `Subtract`\n", "An n-bit subtraction gate.\n", "\n", - "Implements $U|a\\rangle|b\\rangle \\rightarrow |a\\rangle|a-b\\rangle$ using $4n - 4 T$ gates.\n", + "Implements $U|a\\rangle|b\\rangle \\rightarrow |a\\rangle|a-b\\rangle$ using $4n - 4$ T gates.\n", "\n", - "This first negates $b$ (assuming a two's complement representation), and then adds $a$ into it.\n", + "This construction uses the relation `a - b = ~(~a + b)` to turn the operation into addition. This relation is used in\n", + "[Compilation of Fault-Tolerant Quantum Heuristics for Combinatorial Optimization](https://arxiv.org/pdf/2007.07391)\n", + "to turn addition into subtraction conditioned on a qubit.\n", "\n", "#### Parameters\n", " - `a_dtype`: Quantum datatype used to represent the integer a.\n", @@ -48,13 +50,16 @@ "\n", "#### Registers\n", " - `a`: A a_dtype.bitsize-sized input register (register a above).\n", - " - `b`: A b_dtype.bitsize-sized input/output register (register b above).\n" + " - `b`: A b_dtype.bitsize-sized input/output register (register b above). \n", + "\n", + "#### References\n", + " - [Compilation of Fault-Tolerant Quantum Heuristics for Combinatorial Optimization, page 9](https://arxiv.org/pdf/2007.07391). \n" ] }, { "cell_type": "code", "execution_count": null, - "id": "64dd238f", + "id": "0af034cf", "metadata": { "cq.autogen": "Subtract.bloq_doc.py" }, @@ -65,7 +70,7 @@ }, { "cell_type": "markdown", - "id": "4f77ad26", + "id": "9156f4e7", "metadata": { "cq.autogen": "Subtract.example_instances.md" }, @@ -76,7 +81,7 @@ { "cell_type": "code", "execution_count": null, - "id": "2cd50aa6", + "id": "135553b1", "metadata": { "cq.autogen": "Subtract.sub_symb" }, @@ -89,7 +94,7 @@ { "cell_type": "code", "execution_count": null, - "id": "e2862687", + "id": "84c7207a", "metadata": { "cq.autogen": "Subtract.sub_small" }, @@ -101,7 +106,7 @@ { "cell_type": "code", "execution_count": null, - "id": "8b6584fc", + "id": "e98b3206", "metadata": { "cq.autogen": "Subtract.sub_large" }, @@ -113,7 +118,7 @@ { "cell_type": "code", "execution_count": null, - "id": "27312995", + "id": "e003f5ad", "metadata": { "cq.autogen": "Subtract.sub_diff_size_regs" }, @@ -122,9 +127,22 @@ "sub_diff_size_regs = Subtract(QInt(bitsize=4), QInt(bitsize=16))" ] }, + { + "cell_type": "code", + "execution_count": null, + "id": "4aabfa98", + "metadata": { + "cq.autogen": "Subtract.sub_symp_decomposition" + }, + "outputs": [], + "source": [ + "n = sympy.Symbol('n')\n", + "sub_symp_decomposition = Subtract(QInt(bitsize=n)).decompose_bloq()" + ] + }, { "cell_type": "markdown", - "id": "82c580af", + "id": "9b580c5b", "metadata": { "cq.autogen": "Subtract.graphical_signature.md" }, @@ -135,20 +153,20 @@ { "cell_type": "code", "execution_count": null, - "id": "4e508f72", + "id": "aa412d69", "metadata": { "cq.autogen": "Subtract.graphical_signature.py" }, "outputs": [], "source": [ "from qualtran.drawing import show_bloqs\n", - "show_bloqs([sub_symb, sub_small, sub_large, sub_diff_size_regs],\n", - " ['`sub_symb`', '`sub_small`', '`sub_large`', '`sub_diff_size_regs`'])" + "show_bloqs([sub_symb, sub_small, sub_large, sub_diff_size_regs, sub_symp_decomposition],\n", + " ['`sub_symb`', '`sub_small`', '`sub_large`', '`sub_diff_size_regs`', '`sub_symp_decomposition`'])" ] }, { "cell_type": "markdown", - "id": "a282a369", + "id": "6fad45f3", "metadata": { "cq.autogen": "Subtract.call_graph.md" }, @@ -159,7 +177,7 @@ { "cell_type": "code", "execution_count": null, - "id": "43f1884d", + "id": "6a3d1dbc", "metadata": { "cq.autogen": "Subtract.call_graph.py" }, @@ -173,7 +191,7 @@ }, { "cell_type": "markdown", - "id": "12650fb1", + "id": "c4da55d6", "metadata": { "cq.autogen": "SubtractFrom.bloq_doc.md" }, @@ -200,7 +218,7 @@ { "cell_type": "code", "execution_count": null, - "id": "746c25eb", + "id": "ccc21662", "metadata": { "cq.autogen": "SubtractFrom.bloq_doc.py" }, @@ -211,7 +229,7 @@ }, { "cell_type": "markdown", - "id": "75ea545f", + "id": "b8182c9f", "metadata": { "cq.autogen": "SubtractFrom.example_instances.md" }, @@ -222,7 +240,7 @@ { "cell_type": "code", "execution_count": null, - "id": "e7ee57a8", + "id": "4b56cde9", "metadata": { "cq.autogen": "SubtractFrom.sub_from_symb" }, @@ -235,7 +253,7 @@ { "cell_type": "code", "execution_count": null, - "id": "8123a738", + "id": "cc45ae49", "metadata": { "cq.autogen": "SubtractFrom.sub_from_small" }, @@ -247,7 +265,7 @@ { "cell_type": "code", "execution_count": null, - "id": "e469b17a", + "id": "8e40b0fb", "metadata": { "cq.autogen": "SubtractFrom.sub_from_large" }, @@ -258,7 +276,7 @@ }, { "cell_type": "markdown", - "id": "a32af1e3", + "id": "6d965b01", "metadata": { "cq.autogen": "SubtractFrom.graphical_signature.md" }, @@ -269,7 +287,7 @@ { "cell_type": "code", "execution_count": null, - "id": "0fd2dbd7", + "id": "093d2a0d", "metadata": { "cq.autogen": "SubtractFrom.graphical_signature.py" }, @@ -282,7 +300,7 @@ }, { "cell_type": "markdown", - "id": "47e50576", + "id": "420fd950", "metadata": { "cq.autogen": "SubtractFrom.call_graph.md" }, @@ -293,7 +311,7 @@ { "cell_type": "code", "execution_count": null, - "id": "4564d561", + "id": "7d1b9b0e", "metadata": { "cq.autogen": "SubtractFrom.call_graph.py" }, diff --git a/qualtran/bloqs/arithmetic/subtraction.py b/qualtran/bloqs/arithmetic/subtraction.py index 7091ef0ea..0aada6426 100644 --- a/qualtran/bloqs/arithmetic/subtraction.py +++ b/qualtran/bloqs/arithmetic/subtraction.py @@ -11,9 +11,11 @@ # WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. # See the License for the specific language governing permissions and # limitations under the License. + import math from typing import Dict, Optional, Set, Tuple, TYPE_CHECKING, Union +import numpy as np import sympy from attrs import field, frozen @@ -22,6 +24,7 @@ bloq_example, BloqBuilder, BloqDocSpec, + QAny, QInt, QMontgomeryUInt, QUInt, @@ -30,9 +33,12 @@ Soquet, SoquetT, ) +from qualtran._infra.composite_bloq import CompositeBloq from qualtran.bloqs.arithmetic.addition import Add from qualtran.bloqs.arithmetic.bitwise import BitwiseNot -from qualtran.bloqs.arithmetic.negate import Negate +from qualtran.bloqs.basic_gates import OnEach, XGate +from qualtran.bloqs.bookkeeping import Allocate, Cast, Free +from qualtran.bloqs.mcmt.multi_control_multi_target_pauli import MultiTargetCNOT from qualtran.drawing import Text if TYPE_CHECKING: @@ -45,9 +51,11 @@ class Subtract(Bloq): r"""An n-bit subtraction gate. - Implements $U|a\rangle|b\rangle \rightarrow |a\rangle|a-b\rangle$ using $4n - 4 T$ gates. + Implements $U|a\rangle|b\rangle \rightarrow |a\rangle|a-b\rangle$ using $4n - 4$ T gates. - This first negates $b$ (assuming a two's complement representation), and then adds $a$ into it. + This construction uses the relation `a - b = ~(~a + b)` to turn the operation into addition. This relation is used in + [Compilation of Fault-Tolerant Quantum Heuristics for Combinatorial Optimization](https://arxiv.org/pdf/2007.07391) + to turn addition into subtraction conditioned on a qubit. Args: a_dtype: Quantum datatype used to represent the integer a. @@ -58,6 +66,9 @@ class Subtract(Bloq): Registers: a: A a_dtype.bitsize-sized input register (register a above). b: A b_dtype.bitsize-sized input/output register (register b above). + + References: + [Compilation of Fault-Tolerant Quantum Heuristics for Combinatorial Optimization, page 9](https://arxiv.org/pdf/2007.07391) """ a_dtype: Union[QInt, QUInt, QMontgomeryUInt] = field() @@ -112,8 +123,18 @@ def on_classical_vals( ) -> Dict[str, 'ClassicalValT']: unsigned = isinstance(self.a_dtype, (QUInt, QMontgomeryUInt)) b_bitsize = self.b_dtype.bitsize - N = 2**b_bitsize if unsigned else 2 ** (b_bitsize - 1) - return {'a': a, 'b': int(math.fmod(a - b, N))} + N = 2**b_bitsize + if unsigned: + return {'a': a, 'b': int((a - b) % N)} + # Subtraction of signed integers can result in overflow. In most classical programming languages (e.g. C++) + # what happens when an overflow happens is left as an implementation detail for compiler designers. + # However for quantum subtraction the operation should be unitary and that means that the unitary of + # the bloq should be a permutation matrix. + # If we hold `a` constant then the valid range of values of `b` [-N/2, N/2) gets shifted forward or backwards + # by `a`. to keep the operation unitary overflowing values wrap around. this is the same as moving the range [0, N) + # by the same amount modulu $N$. that is add N/2 before subtraction and then remove it. + half_n = N >> 1 + return {'a': a, 'b': int((a - b + half_n) % N) - half_n} def wire_symbol( self, reg: Optional['Register'], idx: Tuple[int, ...] = tuple() @@ -130,14 +151,54 @@ def wire_symbol( raise ValueError() def build_call_graph(self, ssa: 'SympySymbolAllocator') -> Set['BloqCountT']: - return { - (Negate(self.b_dtype), 1), - (Add(self.a_dtype_as_unsigned, self.b_dtype_as_unsigned), 1), - } + delta = self.b_dtype.bitsize - self.a_dtype.bitsize + return ( + { + (OnEach(self.b_dtype.bitsize, XGate()), 3), + (Add(QUInt(self.b_dtype.bitsize), QUInt(self.b_dtype.bitsize)), 1), + } + .union( + [(MultiTargetCNOT(delta), 2)] if delta and isinstance(self.a_dtype, QInt) else [] + ) + .union([(Allocate(QAny(delta)), 1), (Free(QAny(delta)), 1)] if delta else []) + ) def build_composite_bloq(self, bb: 'BloqBuilder', a: Soquet, b: Soquet) -> Dict[str, 'SoquetT']: - b = bb.add(Negate(self.b_dtype), x=b) - a, b = bb.add(Add(self.a_dtype_as_unsigned, self.b_dtype_as_unsigned), a=a, b=b) # a - b + delta = self.b_dtype.bitsize - self.a_dtype.bitsize + n_bits = self.b_dtype.bitsize + if delta: + if not isinstance(delta, int): + raise ValueError( + 'Symbolic decomposition not supported when a_dtype != b_dtype' + ) # pragma: no cover + a_split = bb.split(a) + prefix = bb.allocate(delta) + if isinstance(self.a_dtype, QInt): + a_split[0], prefix = bb.add( + MultiTargetCNOT(delta), control=a_split[0], targets=prefix + ) + prefix_split = bb.split(prefix) + a_split = np.concatenate([prefix_split, a_split]) + a = bb.join(a_split, QUInt(n_bits)) + else: + a = bb.add(Cast(self.a_dtype, QUInt(n_bits)), reg=a) + b = bb.add(Cast(self.b_dtype, QUInt(n_bits)), reg=b) + a = bb.add(OnEach(n_bits, XGate()), q=a) + a, b = bb.add(Add(QUInt(n_bits)), a=a, b=b) + b = bb.add(OnEach(n_bits, XGate()), q=b) + a = bb.add(OnEach(n_bits, XGate()), q=a) + b = bb.add(Cast(QUInt(n_bits), self.b_dtype), reg=b) + if delta: + a_split = bb.split(a) + prefix = bb.join(a_split[:delta]) + if isinstance(self.a_dtype, QInt): + a_split[delta], prefix = bb.add( + MultiTargetCNOT(delta), control=a_split[delta], targets=prefix + ) + prefix = bb.add(Cast(prefix.reg.dtype, QAny(delta)), reg=prefix) + bb.free(prefix) + a = bb.join(a_split[delta:], QUInt(self.a_dtype.bitsize)) + a = bb.add(Cast(QUInt(self.a_dtype.bitsize), self.a_dtype), reg=a) return {'a': a, 'b': b} @@ -166,8 +227,16 @@ def _sub_diff_size_regs() -> Subtract: return sub_diff_size_regs +@bloq_example +def _sub_symp_decomposition() -> CompositeBloq: + n = sympy.Symbol('n') + sub_symp_decomposition = Subtract(QInt(bitsize=n)).decompose_bloq() + return sub_symp_decomposition + + _SUB_DOC = BloqDocSpec( - bloq_cls=Subtract, examples=[_sub_symb, _sub_small, _sub_large, _sub_diff_size_regs] + bloq_cls=Subtract, + examples=[_sub_symb, _sub_small, _sub_large, _sub_diff_size_regs, _sub_symp_decomposition], ) diff --git a/qualtran/bloqs/arithmetic/subtraction_test.py b/qualtran/bloqs/arithmetic/subtraction_test.py index b6c2575bc..5aa3c488b 100644 --- a/qualtran/bloqs/arithmetic/subtraction_test.py +++ b/qualtran/bloqs/arithmetic/subtraction_test.py @@ -12,6 +12,8 @@ # See the License for the specific language governing permissions and # limitations under the License. +import itertools + import numpy as np import pytest @@ -47,15 +49,45 @@ def test_sub_diff_size_regs(bloq_autotester): bloq_autotester(_sub_diff_size_regs) -def test_subtract_bloq_decomposition(): - gate = Subtract(QInt(3), QInt(5)) +@pytest.mark.parametrize( + ['a_bits', 'b_bits'], [(a, b) for a in range(1, 6) for b in range(a, 6) if a + b <= 10] +) +def test_subtract_bloq_decomposition_unsigned(a_bits, b_bits): + gate = Subtract(QUInt(a_bits), QUInt(b_bits)) qlt_testing.assert_valid_bloq_decomposition(gate) - want = np.zeros((256, 256)) - for a_b in range(256): - a, b = a_b >> 5, a_b & 31 - c = (a - b) % 32 - want[(a << 5) | c][a_b] = 1 + tot = 1 << (a_bits + b_bits) + want = np.zeros((tot, tot)) + max_b = 1 << b_bits + for a_b in range(tot): + a, b = a_b >> b_bits, a_b & (max_b - 1) + c = (a - b) % max_b + want[(a << b_bits) | c][a_b] = 1 + got = gate.tensor_contract() + np.testing.assert_allclose(got, want) + + +def _to_signed_binary(x: int, bits: int): + if x >= 0: + return x + return (~(-x) + 1) % (2 << bits) + + +@pytest.mark.parametrize( + ['a_bits', 'b_bits'], [(a, b) for a in range(1, 6) for b in range(a, 6) if a + b <= 10] +) +def test_subtract_bloq_decomposition_signed(a_bits, b_bits): + gate = Subtract(QInt(a_bits + 1), QInt(b_bits + 1)) + qlt_testing.assert_valid_bloq_decomposition(gate) + + tot = 1 << (a_bits + b_bits + 2) + want = np.zeros((tot, tot)) + for a in range(-(1 << a_bits), (1 << a_bits)): + for b in range(-(1 << b_bits), (1 << b_bits)): + a_binary = _to_signed_binary(a, a_bits) + b_binary = _to_signed_binary(b, b_bits) + c_binary = _to_signed_binary(a - b, b_bits) + want[(a_binary << b_bits << 1) | c_binary, (a_binary << b_bits << 1) | b_binary] = 1 got = gate.tensor_contract() np.testing.assert_allclose(got, want) @@ -73,6 +105,38 @@ def test_subtract_bloq_consistent_counts(): ) +@pytest.mark.parametrize('n_bits', range(1, 10)) +def test_t_complexity(n_bits): + complexity = Subtract(QUInt(n_bits)).t_complexity() + assert complexity.t == 4 * (n_bits - 1) + assert complexity.rotations == 0 + + +@pytest.mark.parametrize('dtype', [QInt, QUInt]) +def test_against_classical_values(dtype): + subtract = Subtract(dtype(3), dtype(5)) + cbloq = subtract.decompose_bloq() + if dtype is QInt: + R1 = range(-4, 4) + R2 = range(-16, 16) + else: + R1 = range(8) + R2 = range(32) + for (a, b) in itertools.product(R1, R2): + ref = subtract.call_classically(a=a, b=b) + comp = cbloq.call_classically(a=a, b=b) + assert ref == comp + + +@pytest.mark.parametrize('bitsize', range(2, 5)) +def test_classical_add_signed_overflow(bitsize): + bloq = Subtract(QInt(bitsize)) + cbloq = bloq.decompose_bloq() + mn = -(2 ** (bitsize - 1)) + assert bloq.call_classically(a=0, b=mn) == (0, mn) + assert cbloq.call_classically(a=0, b=mn) == (0, mn) + + def test_sub_from_symb(bloq_autotester): bloq_autotester(_sub_from_symb) @@ -96,7 +160,3 @@ def test_subtract_from_bloq_decomposition(): want[(a << 4) | c][a_b] = 1 got = gate.tensor_contract() np.testing.assert_allclose(got, want) - - -def test_subtract_from_bloq_consistent_counts(): - qlt_testing.assert_equivalent_bloq_counts(SubtractFrom(QInt(3)), generalizer=ignore_split_join) diff --git a/qualtran/simulation/classical_sim.py b/qualtran/simulation/classical_sim.py index 2c74108ea..d62e016e3 100644 --- a/qualtran/simulation/classical_sim.py +++ b/qualtran/simulation/classical_sim.py @@ -62,13 +62,7 @@ def ints_to_bits( w: The bit width of the returned bitstrings. """ x = np.atleast_1d(x) - if not np.issubdtype(x.dtype, np.uint): - assert np.all(x >= 0) - assert np.iinfo(x.dtype).bits <= 64 - x = x.astype(np.uint64) - assert w <= np.iinfo(x.dtype).bits - mask = 2 ** np.arange(w - 1, 0 - 1, -1, dtype=x.dtype).reshape((w, 1)) - return (x & mask).astype(bool).astype(np.uint8).T + return np.array([list(map(int, np.binary_repr(v, width=w))) for v in x], dtype=np.uint8) def _get_in_vals( diff --git a/qualtran/simulation/classical_sim_test.py b/qualtran/simulation/classical_sim_test.py index b9b422d06..bd56280b7 100644 --- a/qualtran/simulation/classical_sim_test.py +++ b/qualtran/simulation/classical_sim_test.py @@ -55,6 +55,10 @@ def test_int_to_bits(): bitstrings = ints_to_bits(nums, w=23) assert bitstrings.shape == (100, 23) + nums = rs.randint(-(2**22), 2**22, size=(100,), dtype=np.int64) + bitstrings = ints_to_bits(nums, w=23) + assert bitstrings.shape == (100, 23) + for num, bs in zip(nums, bitstrings): ref_bs = cirq.big_endian_int_to_bits(int(num), bit_count=23) np.testing.assert_array_equal(ref_bs, bs) @@ -63,9 +67,8 @@ def test_int_to_bits(): (bitstring,) = ints_to_bits(2, w=8) assert bitstring.tolist() == [0, 0, 0, 0, 0, 0, 1, 0] - # check bounds - with pytest.raises(AssertionError): - ints_to_bits([4, -2], w=8) + bitstring = ints_to_bits([31, -1], w=6) + assert bitstring.tolist() == [[0, 1, 1, 1, 1, 1], [1, 1, 1, 1, 1, 1]] def test_dtype_validation():