diff --git a/qiskit/circuit/library/__init__.py b/qiskit/circuit/library/__init__.py index a121165d453b..5306f63cc651 100644 --- a/qiskit/circuit/library/__init__.py +++ b/qiskit/circuit/library/__init__.py @@ -113,6 +113,7 @@ SdgGate SwapGate iSwapGate + SiSwapGate SXGate SXdgGate TGate diff --git a/qiskit/circuit/library/pauli_evolution.py b/qiskit/circuit/library/pauli_evolution.py index 3fa52c00569f..e949291b0398 100644 --- a/qiskit/circuit/library/pauli_evolution.py +++ b/qiskit/circuit/library/pauli_evolution.py @@ -19,7 +19,7 @@ from qiskit.circuit.gate import Gate from qiskit.circuit.parameterexpression import ParameterExpression -from qiskit.synthesis import EvolutionSynthesis, LieTrotter +from qiskit.synthesis.evolution import EvolutionSynthesis, LieTrotter from qiskit.quantum_info import Pauli, SparsePauliOp diff --git a/qiskit/circuit/library/standard_gates/__init__.py b/qiskit/circuit/library/standard_gates/__init__.py index f085c1f3a5b8..b472a70b93b9 100644 --- a/qiskit/circuit/library/standard_gates/__init__.py +++ b/qiskit/circuit/library/standard_gates/__init__.py @@ -11,7 +11,67 @@ # that they have been altered from the originals. """ -Standard gates +============================================================= +Standard gates (:mod:`qiskit.circuit.library.standard_gates`) +============================================================= + +.. autosummary:: + :toctree: ../stubs/ + + C3XGate + C3SXGate + C4XGate + CCXGate + DCXGate + CHGate + CPhaseGate + CRXGate + CRYGate + CRZGate + CSwapGate + CSXGate + CUGate + CU1Gate + CU3Gate + CXGate + CYGate + CZGate + CCZGate + HGate + IGate + MCPhaseGate + PhaseGate + RCCXGate + RC3XGate + RXGate + RXXGate + RYGate + RYYGate + RZGate + RZZGate + RZXGate + XXMinusYYGate + XXPlusYYGate + ECRGate + SGate + SdgGate + CSGate + CSdgGate + SwapGate + iSwapGate + SiSwapGate + SXGate + SXdgGate + TGate + TdgGate + UGate + U1Gate + U2Gate + U3Gate + XGate + YGate + ZGate + """ from .h import HGate, CHGate @@ -31,6 +91,7 @@ from .s import SGate, SdgGate, CSGate, CSdgGate from .swap import SwapGate, CSwapGate from .iswap import iSwapGate +from .siswap import SiSwapGate from .sx import SXGate, SXdgGate, CSXGate from .dcx import DCXGate from .t import TGate, TdgGate @@ -99,6 +160,7 @@ def get_standard_gate_name_mapping(): CSdgGate(), SwapGate(), iSwapGate(), + SiSwapGate(), SXdgGate(), TGate(), TdgGate(), diff --git a/qiskit/circuit/library/standard_gates/equivalence_library.py b/qiskit/circuit/library/standard_gates/equivalence_library.py index fa17cf1ba3e0..fcc0e6f9c065 100644 --- a/qiskit/circuit/library/standard_gates/equivalence_library.py +++ b/qiskit/circuit/library/standard_gates/equivalence_library.py @@ -41,6 +41,7 @@ SwapGate, CSwapGate, iSwapGate, + SiSwapGate, SXGate, SXdgGate, CSXGate, @@ -762,6 +763,23 @@ def_iswap.append(inst, qargs, cargs) _sel.add_equivalence(iSwapGate(), def_iswap) +# SiSwapGate +# +# ┌──────────┐ ┌────────────┐┌────────────┐ +# q_0: ┤0 ├ q_0: ┤0 ├┤0 ├ +# │ Siswap │ ≡ │ Rxx(-π/4) ││ Ryy(-π/4) │ +# q_1: ┤1 ├ q_1: ┤1 ├┤1 ├ +# └──────────┘ └────────────┘└────────────┘ + +q = QuantumRegister(2, "q") +def_siswap = QuantumCircuit(q) +for inst, qargs, cargs in [ + (RXXGate(-pi / 4), [q[0], q[1]], []), + (RYYGate(-pi / 4), [q[0], q[1]], []), +]: + def_siswap.append(inst, qargs, cargs) +_sel.add_equivalence(SiSwapGate(), def_siswap) + # SXGate # global phase: π/4 # ┌────┐ ┌─────┐┌───┐┌─────┐ diff --git a/qiskit/circuit/library/standard_gates/siswap.py b/qiskit/circuit/library/standard_gates/siswap.py new file mode 100644 index 000000000000..1c0ff5f0d5fb --- /dev/null +++ b/qiskit/circuit/library/standard_gates/siswap.py @@ -0,0 +1,77 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2017, 2023. +# +# 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. + +"""sqrt(iSWAP) gate.""" + +import numpy as np +from qiskit.circuit.gate import Gate + + +class SiSwapGate(Gate): + r"""sqrt(iSWAP) gate. + + A 2-qubit symmetric gate from the iSWAP (or XY) family. + It has Weyl chamber coordinates (π/8, π/8, 0). + + Can be applied to a :class:`~qiskit.circuit.QuantumCircuit` + with the :meth:`~qiskit.circuit.QuantumCircuit.SiSwap` method. + + .. parsed-literal:: + + ┌────────────┐┌────────────┐ + q_0: ┤0 ├┤0 ├ + │ Rxx(-π/4) ││ Ryy(-π/4) │ + q_1: ┤1 ├┤1 ├ + └────────────┘└────────────┘ + + .. math:: + \sqrt{\mathrm{iSwap}}\ q_0, q_1 = + \begin{pmatrix} + 1 & 0 & 0 & 0 \\ + 0 & \frac{1}{\sqrt{2}} & \frac{i}{\sqrt{2}} & 0 \\ + 0 & \frac{i}{\sqrt{2}} & \frac{1}{\sqrt{2}} & 0 \\ + 0 & 0 & 0 & 1 + \end{pmatrix} + """ + + def __init__(self): + """Create new SiSwap gate.""" + super().__init__("siswap", 2, []) + + def _define(self): + """ + gate siswap a, b { rxx(-pi/4) a, b; ryy(-pi/4) a, b; } + """ + # pylint: disable=cyclic-import + from qiskit.circuit.quantumcircuit import QuantumRegister, QuantumCircuit + from .rxx import RXXGate + from .ryy import RYYGate + + q = QuantumRegister(2, "q") + qc = QuantumCircuit(q, name=self.name) + rules = [(RXXGate(-np.pi / 4), [q[0], q[1]], []), (RYYGate(-np.pi / 4), [q[0], q[1]], [])] + for instr, qargs, cargs in rules: + qc._append(instr, qargs, cargs) + + self.definition = qc + + def __array__(self, dtype=None): + """Return a numpy.array for the DCX gate.""" + return np.array( + [ + [1, 0, 0, 0], + [0, 1 / np.sqrt(2), 1j / np.sqrt(2), 0], + [0, 1j / np.sqrt(2), 1 / np.sqrt(2), 0], + [0, 0, 0, 1], + ], + dtype=dtype, + ) diff --git a/qiskit/quantum_info/synthesis/two_qubit_decompose.py b/qiskit/quantum_info/synthesis/two_qubit_decompose.py index 6e580265e944..b232cb7fec34 100644 --- a/qiskit/quantum_info/synthesis/two_qubit_decompose.py +++ b/qiskit/quantum_info/synthesis/two_qubit_decompose.py @@ -1092,17 +1092,17 @@ def __call__( basis_fidelity: Optional[float] = None, approximate: bool = True, *, - _num_basis_uses: int = None, + _num_basis_uses: Optional[int] = None, ) -> QuantumCircuit: """Decompose a two-qubit `unitary` over fixed basis + SU(2) using the best approximation given that each basis application has a finite `basis_fidelity`. Args: - unitary (Operator or ndarray): 4x4 unitary to synthesize. - basis_fidelity (float or None): Fidelity to be assumed for applications of KAK Gate. + unitary: 4x4 unitary to synthesize. + basis_fidelity: Fidelity to be assumed for applications of KAK Gate. If given, overrides basis_fidelity given at init. - approximate (bool): Approximates if basis fidelities are less than 1.0. - _num_basis_uses (int): force a particular approximation by passing a number in [0, 3]. + approximate: Approximates if basis fidelities are less than 1.0. + _num_basis_uses: force a particular approximation by passing a number in [0, 3]. Returns: QuantumCircuit: Synthesized circuit. Raises: diff --git a/qiskit/quantum_info/synthesis/xx_decompose/decomposer.py b/qiskit/quantum_info/synthesis/xx_decompose/decomposer.py index 3d5714b5c2ab..55da6caaf2f3 100644 --- a/qiskit/quantum_info/synthesis/xx_decompose/decomposer.py +++ b/qiskit/quantum_info/synthesis/xx_decompose/decomposer.py @@ -15,7 +15,6 @@ """ from __future__ import annotations import heapq -import math from operator import itemgetter from typing import Callable @@ -27,31 +26,13 @@ from qiskit.quantum_info.operators import Operator from qiskit.quantum_info.synthesis.one_qubit_decompose import ONE_QUBIT_EULER_BASIS_GATES from qiskit.quantum_info.synthesis.two_qubit_decompose import TwoQubitWeylDecomposition +from qiskit.synthesis.su4.utils import average_infidelity from .circuits import apply_reflection, apply_shift, canonical_xx_circuit from .utilities import EPSILON from .polytopes import XXPolytope -def _average_infidelity(p, q): - """ - Computes the infidelity distance between two points p, q expressed in positive canonical - coordinates. - """ - - a0, b0, c0 = p - a1, b1, c1 = q - - return 1 - 1 / 20 * ( - 4 - + 16 - * ( - math.cos(a0 - a1) ** 2 * math.cos(b0 - b1) ** 2 * math.cos(c0 - c1) ** 2 - + math.sin(a0 - a1) ** 2 * math.sin(b0 - b1) ** 2 * math.sin(c0 - c1) ** 2 - ) - ) - - class XXDecomposer: """ A class for optimal decomposition of 2-qubit unitaries into 2-qubit basis gates of XX type @@ -155,7 +136,7 @@ def _best_decomposition(canonical_coordinate, available_strengths): strength_polytope = XXPolytope.from_strengths(*[x / 2 for x in sequence]) candidate_point = strength_polytope.nearest(canonical_coordinate) - candidate_cost = sequence_cost + _average_infidelity( + candidate_cost = sequence_cost + average_infidelity( canonical_coordinate, candidate_point ) diff --git a/qiskit/synthesis/__init__.py b/qiskit/synthesis/__init__.py index eb3adae2955a..73246e48e442 100644 --- a/qiskit/synthesis/__init__.py +++ b/qiskit/synthesis/__init__.py @@ -87,6 +87,13 @@ SolovayKitaevDecomposition generate_basic_approximations +SU(4) (2-qubit unitary) Synthesis +================================= + +.. autosummary:: + :toctree: ../stubs/ + + SiSwapDecomposer """ from .evolution import ( @@ -122,3 +129,4 @@ synth_cnotdihedral_general, ) from .discrete_basis import SolovayKitaevDecomposition, generate_basic_approximations +from .su4 import SiSwapDecomposer diff --git a/qiskit/synthesis/su4/__init__.py b/qiskit/synthesis/su4/__init__.py new file mode 100644 index 000000000000..1203e78389fd --- /dev/null +++ b/qiskit/synthesis/su4/__init__.py @@ -0,0 +1,15 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2023. +# +# 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. + +"""Two-qubit unitary synthesis methods.""" + +from .siswap_decompose import SiSwapDecomposer diff --git a/qiskit/synthesis/su4/siswap_decompose.py b/qiskit/synthesis/su4/siswap_decompose.py new file mode 100644 index 000000000000..6ed5a7e4584a --- /dev/null +++ b/qiskit/synthesis/su4/siswap_decompose.py @@ -0,0 +1,286 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2023 +# +# 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. + +# pylint: disable=invalid-name, non-ascii-name + +"""Synthesis of two-qubit unitaries using at most 3 applications of the sqrt(iSWAP) gate.""" + +import cmath +from typing import Optional, List, Union +import numpy as np + +from qiskit.circuit.quantumcircuit import QuantumCircuit +from qiskit.quantum_info.operators import Operator +from qiskit.quantum_info.synthesis.two_qubit_decompose import TwoQubitWeylDecomposition +from qiskit.circuit.library import SiSwapGate, XGate, YGate, ZGate +from qiskit.synthesis.su4.utils import find_min_point, average_infidelity + + +_EPS = 1e-12 + +# the polytope accessible by two applications of the SiSwapGate +# (red region in paper) +_ID = np.array([0, 0, 0]) +_CNOT = np.array([np.pi / 4, 0, 0]) +_ISWAP = np.array([np.pi / 4, np.pi / 4, 0]) +_MID = np.array([np.pi / 4, np.pi / 8, np.pi / 8]) +_MIDDG = np.array([np.pi / 4, np.pi / 8, -np.pi / 8]) +_SISWAP = np.array([np.pi / 8, np.pi / 8, 0]) +_POLYTOPE = np.array([_ID, _CNOT, _ISWAP, _MID, _MIDDG]) + + +class SiSwapDecomposer: + """ + A class for decomposing 2-qubit unitaries into at most 3 uses of the sqrt(iSWAP) gate. + + This basis is attractive in that it is shorter/cheaper than iSWAP and generates a large + volume (~79%) of Haar-random SU(4) unitaries with only two applications, and an even + larger volume when approximation is enabled. + + Args: + euler_basis (list(str)): single-qubit gates basis in the decomposition. + + Reference: + 1. C. Huang et al, + Towards ultra-high fidelity quantum operations: + SQiSW gate as a native two-qubit gate (2021) + `arXiv:2105.06074 `_ + """ + + def __init__(self, euler_basis: Optional[List[str]]): + # decomposer for the local single-qubit gates + from qiskit.transpiler.passes.optimization.optimize_1q_decomposition import ( + Optimize1qGatesDecomposition, # pylint: disable=cyclic-import + ) + + euler_basis = euler_basis or ["u"] + self._decomposer1q = Optimize1qGatesDecomposition(euler_basis) + + def __call__( + self, + unitary: Union[Operator, np.ndarray], + basis_fidelity: Optional[float] = 1.0, + approximate: bool = True, + *, + _num_basis_uses: Optional[int] = None, + ) -> QuantumCircuit: + """Decompose a two-qubit unitary into using the sqrt(iSWAP) gate. + + Args: + unitary: A 4x4 unitary to synthesize. + basis_fidelity: Fidelity of the iSWAP gate. + approximate: Approximates if basis fidelities are less than 1.0. + _num_basis_uses: Force a particular approximation by passing a number in [0, 3]. + + Returns: + QuantumCircuit: Synthesized circuit. + """ + if not approximate: + basis_fidelity = 1.0 + + u_decomp = TwoQubitWeylDecomposition(Operator(unitary)) + + x = u_decomp.a + y = u_decomp.b + z = u_decomp.c + + A1 = Operator(u_decomp.K1r) + A2 = Operator(u_decomp.K1l) + B1 = Operator(u_decomp.K2r) + B2 = Operator(u_decomp.K2l) + + p = np.array([x, y, z]) + + if abs(z) <= x - y + _EPS: + polytope_projection = p + else: + polytope_projection = find_min_point([v - p for v in _POLYTOPE]) + p + + candidate_points = [ + _ID, # 0 applications + _SISWAP, # 1 application + polytope_projection, # 2 applications + p, # 3 applications + ] + + if _num_basis_uses is None: + expected_fidelities = [ + (1 - average_infidelity(p, q)) * basis_fidelity**i + for i, q in enumerate(candidate_points) + ] + best_nbasis = int(np.argmax(expected_fidelities)) # tiebreaks with smallest + else: + best_nbasis = _num_basis_uses + + p = candidate_points[best_nbasis] + + circuit = _can_circuit(*p) + circuit = circuit.compose(B1.to_instruction(), [0], front=True) + circuit = circuit.compose(B2.to_instruction(), [1], front=True) + circuit = circuit.compose(A1.to_instruction(), [0]) + circuit = circuit.compose(A2.to_instruction(), [1]) + + # FIXME: cleaner to track global phase during construction + i = np.where(~np.isclose(np.ravel(Operator(circuit).data), 0.0))[0][0] + phase_diff = cmath.phase(Operator(unitary).data.flat[i] / Operator(circuit).data.flat[i]) + circuit.global_phase += phase_diff + + return self._decomposer1q(circuit) + + +def _can_circuit(x, y, z): + """ + Find a circuit that implements the canonical gate + CAN(x, y, z) := RXX(x) . RYY(y) . RZZ(z) + """ + circuit = QuantumCircuit(2) + if np.allclose([x, y, z], _ID): + pass + elif np.allclose([x, y, z], _SISWAP): + circuit.append(SiSwapGate(), [0, 1]) + elif abs(z) <= x - y + _EPS: # red region + V = _interleaving_single_qubit_gates(x, y, z) + v_decomp = TwoQubitWeylDecomposition(Operator(V)) + can = _remove_pre_post(V, v_decomp) + circuit = circuit.compose(can, [0, 1]) + + # Due to the symmetry of Weyl chamber CAN(pi/4, y, z) ~ CAN(pi/4, y, -z) + # we might get a V operator that implements CAN(pi/4, y, -z) instead + # we catch this case and fix it by local gates. + if not np.isclose(v_decomp.c, z): + circuit = circuit.compose(XGate(), [0], front=True) + circuit = circuit.compose(ZGate(), [1], front=True) + circuit = circuit.compose(YGate(), [1]) + else: + if z < 0: # blue region + # CAN(x, y, z) ~ CAN(x, y, -z)† + # so we decompose the adjoint, replace SiSwap with a template in + # terms of SiSwap†, then invert the whole thing. + inverse_circuit = _can_circuit(x, y, -z) + inverse_circuit_with_siswap_dg = QuantumCircuit(2) + for instruction in inverse_circuit: + if isinstance(instruction.operation, SiSwapGate): + inverse_circuit_with_siswap_dg.z(0) + inverse_circuit_with_siswap_dg.append(SiSwapGate().inverse(), [0, 1]) + inverse_circuit_with_siswap_dg.z(0) + else: + inverse_circuit_with_siswap_dg.append(instruction) + + V = inverse_circuit_with_siswap_dg.inverse() + v_decomp = TwoQubitWeylDecomposition(Operator(V)) + can = _remove_pre_post(V, v_decomp) + circuit = circuit.compose(can, [0, 1]) + else: + # make a circuit using 1 SiSwap that is able to bring a red point to here + follow = QuantumCircuit(2) + + # x -> x + pi/8 + # y -> y + pi/8 + follow = follow.compose(SiSwapGate(), [0, 1]) + + eigenphase_crossing = False + if x > np.pi / 8: # green region + # y -> y + pi/8 + # z -> z + pi/8 + follow = follow.compose(YGate().power(1 / 2), [0], front=True) + follow = follow.compose(YGate().power(1 / 2), [1], front=True) + follow = follow.compose(YGate().power(-1 / 2), [0]) + follow = follow.compose(YGate().power(-1 / 2), [1]) + if y + z < np.pi / 4: # eigenphase crossing condition: a_2 - pi/4 < a_3 + pi/4 + eigenphase_crossing = True + # corresponding red coordinates + y -= np.pi / 8 + z -= np.pi / 8 + else: # purple region + # x -> x - pi/8 + # z -> z + pi/8 + follow = follow.compose(XGate().power(1 / 2), [0], front=True) + follow = follow.compose(XGate().power(1 / 2), [1], front=True) + follow = follow.compose(XGate().power(-1 / 2), [0]) + follow = follow.compose(XGate().power(-1 / 2), [1]) + follow = follow.compose(ZGate(), [0], front=True) + follow = follow.compose(ZGate(), [0]) + if y + z < np.pi / 8: # eigenphase crossing condition: a_2 - pi/4 < a_3 + eigenphase_crossing = True + # corresponding red coordinates + x += np.pi / 8 + z -= np.pi / 8 + + if eigenphase_crossing: + y, z = -z, -y + + # final 3xSiSwap circuit: red --> fix crossing --> green or purple + red_decomp = _can_circuit(x, y, z) + circuit = circuit.compose(red_decomp, [0, 1]) + if eigenphase_crossing: + # y, z -> -z, -y + circuit = circuit.compose(XGate().power(1 / 2), [0], front=True) + circuit = circuit.compose(XGate().power(1 / 2), [1], front=True) + circuit = circuit.compose(XGate().power(-1 / 2), [0]) + circuit = circuit.compose(XGate().power(-1 / 2), [1]) + circuit = circuit.compose(XGate(), [0], front=True) + circuit = circuit.compose(XGate(), [0]) + circuit = circuit.compose(follow, [0, 1]) + + return circuit + + +def _interleaving_single_qubit_gates(x, y, z): + """ + Find the single-qubit gates given the interaction coefficients + (x, y, z) ∈ W′ when sandwiched by two SiSwap gates. + Return the SiSwap sandwich. + """ + C = np.sin(x + y - z) * np.sin(x - y + z) * np.sin(-x - y - z) * np.sin(-x + y + z) + if abs(C) < _EPS: + C = 0.0 + + α = np.arccos(np.cos(2 * x) - np.cos(2 * y) + np.cos(2 * z) + 2 * np.sqrt(C)) + + β = np.arccos(np.cos(2 * x) - np.cos(2 * y) + np.cos(2 * z) - 2 * np.sqrt(C)) + + s = 4 * (np.cos(x) ** 2) * (np.cos(z) ** 2) * (np.sin(y) ** 2) + t = np.cos(2 * x) * np.cos(2 * y) * np.cos(2 * z) + sign_z = 1 if z >= 0 else -1 + γ = np.arccos(sign_z * np.sqrt(s / (s + t))) + + # create V (sandwich) operator + V = QuantumCircuit(2) + V.append(SiSwapGate(), [0, 1]) + V.rz(γ, 0) + V.rx(α, 0) + V.rz(γ, 0) + V.rx(β, 1) + V.append(SiSwapGate(), [0, 1]) + + # the returned circuit is the SiSwap sandwich + return V + + +def _remove_pre_post(circuit, decomp): + """ + Given a circuit and its Weyl decomposition, multiply local gates before and after + it to get a new circuit which is equivalent to RXX.RYY.RZZ (up to global phase). + """ + D1 = Operator(decomp.K1r) + D2 = Operator(decomp.K1l) + E1 = Operator(decomp.K2r) + E2 = Operator(decomp.K2l) + + new_circuit = QuantumCircuit(2) + new_circuit.append(E1.adjoint(), [0]) + new_circuit.append(E2.adjoint(), [1]) + new_circuit.compose(circuit, inplace=True) + new_circuit.append(D1.adjoint(), [0]) + new_circuit.append(D2.adjoint(), [1]) + + return new_circuit diff --git a/qiskit/synthesis/su4/utils.py b/qiskit/synthesis/su4/utils.py new file mode 100644 index 000000000000..3890057e0e9b --- /dev/null +++ b/qiskit/synthesis/su4/utils.py @@ -0,0 +1,88 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2023 +# +# 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. + +# pylint: disable=invalid-name + +"""Utilities for SU(4) synthesis: Weyl chamber, polytopes, distances, etc.""" + +import numpy as np + + +_EPS = 1e-12 + + +def find_min_point(P): + """Find the closest point on convex polytope P to the origin (in Euclidean distance). + + P is given as a list of points forming its vertices. + + Function from + https://scipy-cookbook.readthedocs.io/items/Finding_Convex_Hull_Minimum_Point.html + """ + if len(P) == 1: + return P[0] + + P = [np.array(p) for p in P] + + # Step 0. Choose a point from C(P) + x = P[np.array([np.dot(p, p) for p in P]).argmin()] + + while True: + # Step 1. \alpha_k := min{x_{k-1}^T p | p \in P} + p_alpha = P[np.array([np.dot(x, p) for p in P]).argmin()] + + if np.dot(x, x - p_alpha) < _EPS: + return np.array([float(i) for i in x]) + + Pk = [p for p in P if abs(np.dot(x, p - p_alpha)) < _EPS] + + # Step 2. P_k := { p | p \in P and x_{k-1}^T p = \alpha_k} + P_Pk = [p for p in P if not np.array([(p == q).all() for q in Pk]).any()] + + if len(Pk) == len(P): + return np.array([float(i) for i in x]) + + y = find_min_point(Pk) + + p_beta = P_Pk[np.array([np.dot(y, p) for p in P_Pk]).argmin()] + + if np.dot(y, y - p_beta) < _EPS: + return np.array([float(i) for i in y]) + + # Step 4. + P_aux = [p for p in P_Pk if (np.dot(y - x, y - p) > _EPS) and (np.dot(x, y - p) != 0)] + p_lambda = P_aux[np.array([np.dot(y, y - p) / np.dot(x, y - p) for p in P_aux]).argmin()] + lam = np.dot(x, p_lambda - y) / np.dot(y - x, y - p_lambda) + + x += lam * (y - x) + + +def average_infidelity(p, q): + """Computes the infidelity distance between two points p, q expressed in + positive canonical coordinates. + + Average gate fidelity is :math:`Fbar = (d + |Tr (Utarget \\cdot U^dag)|^2) / d(d+1)` + + [1] M. Horodecki, P. Horodecki and R. Horodecki, PRA 60, 1888 (1999) + """ + + a0, b0, c0 = p + a1, b1, c1 = q + + return 1 - 1 / 20 * ( + 4 + + 16 + * ( + np.cos(a0 - a1) ** 2 * np.cos(b0 - b1) ** 2 * np.cos(c0 - c1) ** 2 + + np.sin(a0 - a1) ** 2 * np.sin(b0 - b1) ** 2 * np.sin(c0 - c1) ** 2 + ) + ) diff --git a/releasenotes/notes/sqisw-decomposition-2478b9f0b7c54044.yaml b/releasenotes/notes/sqisw-decomposition-2478b9f0b7c54044.yaml new file mode 100644 index 000000000000..8f8a23ad6ef9 --- /dev/null +++ b/releasenotes/notes/sqisw-decomposition-2478b9f0b7c54044.yaml @@ -0,0 +1,12 @@ +--- +features: + - | + Added a new gate to the standard_gates library :class:`.SiSwapGate` + to represent the square-root of an iSWAP. This gate has the desirable + properties of being faster/higher-fidelity than a full iSWAP, but also + synthesizing the volume of two-qubit unitaries with fewer applications. + + An associated decomposition method that synthesizes arbitrary two-qubit + unitaries into two or three applications of SiSwap has been added in + :class:`.SiSwapDecomposer`. + diff --git a/test/python/synthesis/test_siswap_synthesis.py b/test/python/synthesis/test_siswap_synthesis.py new file mode 100644 index 000000000000..8dc87075ecd2 --- /dev/null +++ b/test/python/synthesis/test_siswap_synthesis.py @@ -0,0 +1,137 @@ +# This code is part of Qiskit. +# +# (C) Copyright IBM 2017, 2019. +# +# 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. + +# pylint: disable=invalid-name + +"""Test synthesis of two-qubit unitaries into SiSwap gates.""" + +import unittest +from test import combine +from ddt import ddt +import numpy as np + +from qiskit.synthesis.su4 import SiSwapDecomposer +from qiskit.quantum_info import random_unitary, Operator, average_gate_fidelity +from qiskit.test import QiskitTestCase +from qiskit.circuit.library import ( + SwapGate, + iSwapGate, + CXGate, + IGate, + XGate, + RXXGate, + RYYGate, + RZZGate, +) + +_EPS = 1e-12 + + +@ddt +class TestSiSwapSynth(QiskitTestCase): + """Test synthesis of SU(4)s over SiSwap basis.""" + + @combine(seed=range(50)) + def test_siswap_random(self, seed): + """Test synthesis of 50 random SU(4)s.""" + u = random_unitary(4, seed=seed) + decomposer = SiSwapDecomposer(euler_basis=["rz", "ry"]) + circuit = decomposer(u) + self.assertLessEqual(circuit.count_ops().get("siswap", 0), 3) + self.assertEqual(Operator(circuit), Operator(u)) + + @combine(seed=range(50, 100), basis_fidelity=[0.98, 0.99, 1.0]) + def test_siswap_random_approx(self, seed, basis_fidelity): + """Test synthesis of 50 random SU(4)s with approximation.""" + u = random_unitary(4, seed=seed) + decomposer = SiSwapDecomposer(euler_basis=["rz", "ry"]) + circuit = decomposer(u, basis_fidelity=basis_fidelity) + self.assertLessEqual(circuit.count_ops().get("siswap", 0), 3) + self.assertGreaterEqual(average_gate_fidelity(Operator(circuit), u), basis_fidelity - _EPS) + + @combine(corner=[SwapGate(), SwapGate().power(1 / 2), SwapGate().power(1 / 32)]) + def test_siswap_corners_3_use(self, corner): + """Test synthesis of some 3-use special corner cases.""" + u = Operator(corner) + decomposer = SiSwapDecomposer(euler_basis=["rz", "ry"]) + circuit = decomposer(u) + self.assertEqual(circuit.count_ops().get("siswap", 0), 3) + self.assertEqual(Operator(circuit), Operator(u)) + + @combine( + corner=[ + iSwapGate(), + CXGate(), + CXGate().power(-1 / 2), + ] + ) + def test_siswap_corners_2_uses(self, corner): + """Test synthesis of some 2-use special corner cases.""" + u = Operator(corner) + decomposer = SiSwapDecomposer(euler_basis=["u"]) + circuit = decomposer(u) + self.assertEqual(circuit.count_ops().get("siswap", 0), 2) + self.assertEqual(Operator(circuit), Operator(u)) + + @combine( + corner=[ + iSwapGate().power(1 / 2), + iSwapGate().power(-1 / 2), + ] + ) + def test_siswap_corners_1_use(self, corner): + """Test synthesis of some 1-use special corner cases.""" + u = Operator(corner) + decomposer = SiSwapDecomposer(euler_basis=["u"]) + circuit = decomposer(u) + self.assertEqual(circuit.count_ops().get("siswap", 0), 1) + self.assertEqual(Operator(circuit), Operator(u)) + + @combine(corner=[Operator(IGate()) ^ Operator(IGate()), Operator(IGate()) ^ Operator(XGate())]) + def test_siswap_corners_0_use(self, corner): + """Test synthesis of some 0-use special corner cases.""" + u = Operator(corner) + decomposer = SiSwapDecomposer(euler_basis=["u"]) + circuit = decomposer(u) + self.assertEqual(circuit.count_ops().get("siswap", 0), 0) + self.assertEqual(Operator(circuit), Operator(u)) + + @combine( + p=[ + [np.pi / 4, np.pi / 6, -np.pi / 12], # red and blue intersection (bottom face) + [np.pi / 4, np.pi / 6, np.pi / 12], # inverse of above (and locally equivalent) + [np.pi / 4, np.pi / 8, np.pi / 8], # half-way between CX and SWAP (peak red polytope) + [np.pi / 4, np.pi / 8, -np.pi / 8], # inverse of above (and locally equivalent) + [np.pi / 4, np.pi / 16, np.pi / 16], # quarter-way between CX and SWAP + [np.pi / 4, np.pi / 16, np.pi / 16], # inverse of above (and locally equivalent) + [np.pi / 6, np.pi / 8, np.pi / 24], # red and blue and green intersection + [np.pi / 6, np.pi / 8, -np.pi / 24], # inverse of above (not locally equivalent) + [np.pi / 16, np.pi / 24, np.pi / 48], # red and blue and purple intersection + [np.pi / 16, np.pi / 24, -np.pi / 48], # inverse of above (not locally equivalent) + [np.pi / 6, np.pi / 12, -np.pi / 16], # inside red polytope + ] + ) + def test_siswap_special_points(self, p): + """Test special points in the Weyl chamber related to SiSwap polytopes.""" + u = ( + Operator(RXXGate(-2 * p[0])) + & Operator(RYYGate(-2 * p[1])) + & Operator(RZZGate(-2 * p[2])) + ) + decomposer = SiSwapDecomposer(euler_basis=["u"]) + circuit = decomposer(u) + self.assertEqual(circuit.count_ops().get("siswap", 0), 2) + self.assertEqual(Operator(circuit), Operator(u)) + + +if __name__ == "__main__": + unittest.main()