Add ABADecomposer class

opensquirrel/decomposer/aba_decomposer.py

@@ -0,0 +1,102 @@
+from __future__ import annotations
+
+import math
+from collections.abc import Callable, Iterable
+
+from opensquirrel.common import ATOL
+from opensquirrel.decomposer.general_decomposer import Decomposer
+from opensquirrel.ir import BlochSphereRotation, Float, Gate
+from opensquirrel.utils.identity_filter import filter_out_identities
+
+
+class ABADecomposer(Decomposer):
+    def __init__(self, ra: Callable[..., BlochSphereRotation], rb: Callable[..., BlochSphereRotation]):
+        self.ra = ra
+        self.rb = rb
+
+    @staticmethod
+    def get_decomposition_angles(alpha: float, axis: Iterable[float]) -> tuple[float, float, float]:
+        """
+        Gives the angles used in the A-B-A decomposition of the Bloch sphere rotation
+        characterized by a rotation around `axis` of angle `alpha`.
+
+        Parameters:
+            alpha: angle of the Bloch sphere rotation
+            axis: _normalized_ axis of the Bloch sphere rotation
+
+        Returns:
+            A triple (theta1, theta2, theta3) corresponding to the decomposition of the
+            arbitrary Bloch sphere rotation into U = Ra(theta3) Rb(theta2) Ra(theta1)
+
+        """
+        nx, ny, nz = axis
+
+        assert abs(nx**2 + ny**2 + nz**2 - 1) < ATOL, "Axis needs to be normalized"
+
+        assert -math.pi + ATOL < alpha <= math.pi + ATOL, "Angle needs to be normalized"
+
+        if abs(alpha - math.pi) < ATOL:
+            # alpha == pi, math.tan(alpha / 2) is not defined.
+
+            p: float
+            if abs(nz) < ATOL:
+                theta2 = math.pi
+                p = 0
+                m = 2 * math.acos(ny)
+
+            else:
+                p = math.pi
+                theta2 = 2 * math.acos(nz)
+
+                if abs(nz - 1) < ATOL or abs(nz + 1) < ATOL:
+                    m = p  # This can be anything, but setting m = p means theta3 == 0, which is better for gate count.
+                else:
+                    m = 2 * math.acos(ny / math.sqrt(1 - nz**2))
+
+        else:
+            p = 2 * math.atan2(nz * math.sin(alpha / 2), math.cos(alpha / 2))
+
+            acos_argument = math.cos(alpha / 2) * math.sqrt(1 + (nz * math.tan(alpha / 2)) ** 2)
+
+            # This fixes float approximations like 1.0000000000002 which acos doesn't like.
+            acos_argument = max(min(acos_argument, 1.0), -1.0)
+
+            theta2 = 2 * math.acos(acos_argument)
+            theta2 = math.copysign(theta2, alpha)
+
+            if abs(math.sin(theta2 / 2)) < ATOL:
+                m = p  # This can be anything, but setting m = p means theta3 == 0, which is better for gate count.
+            else:
+                acos_argument = ny * math.sin(alpha / 2) / math.sin(theta2 / 2)
+
+                # This fixes float approximations like 1.0000000000002 which acos doesn't like.
+                acos_argument = max(min(acos_argument, 1.0), -1.0)
+
+                m = 2 * math.acos(acos_argument)
+
+        theta1 = (p + m) / 2
+
+        theta3 = p - theta1
+
+        return theta1, theta2, theta3
+
+    def decompose(self, g: Gate) -> list[Gate]:
+        if not isinstance(g, BlochSphereRotation):
+            # Only decomposer single-qubit gates.
+            return [g]
+
+        theta1, theta2, theta3 = self.get_decomposition_angles(g.angle, g.axis)
+
+        a1 = self.ra(g.qubit, Float(theta1))
+        b = self.rb(g.qubit, Float(theta2))
+        a2 = self.ra(g.qubit, Float(theta3))
+
+        # Note: written like this, the decomposition doesn't preserve the global phase, which is fine
+        # since the global phase is a physically irrelevant artifact of the mathematical
+        # model we use to describe the quantum system.
+
+        # Should we want to preserve it, we would need to use a raw BlochSphereRotation, which would then
+        # be an anonymous gate in the resulting decomposed circuit:
+        # z2 = BlochSphereRotation(qubit=g.qubit, angle=theta3, axis=(0, 0, 1), phase = g.phase)
+
+        return filter_out_identities([a1, b, a2])

opensquirrel/decomposer/cnot_decomposer.py

@@ -2,7 +2,7 @@
 
 from opensquirrel.common import ATOL
 from opensquirrel.decomposer.general_decomposer import Decomposer
-from opensquirrel.decomposer.zyz_decomposer import get_zyz_decomposition_angles
+from opensquirrel.decomposer.zyz_decomposer import ZYZDecomposer
 from opensquirrel.default_gates import CNOT, Ry, Rz, X
 from opensquirrel.ir import BlochSphereRotation, ControlledGate, Float, Gate
 from opensquirrel.merger import general_merger
@@ -38,7 +38,7 @@ def decompose(g: Gate) -> [Gate]:
 
         # Try special case first, see https://arxiv.org/pdf/quant-ph/9503016.pdf lemma 5.5
         controlled_rotation_times_x = general_merger.compose_bloch_sphere_rotations(X(target_qubit), g.target_gate)
-        theta0_with_x, theta1_with_x, theta2_with_x = get_zyz_decomposition_angles(
+        theta0_with_x, theta1_with_x, theta2_with_x = ZYZDecomposer().get_decomposition_angles(
             controlled_rotation_times_x.angle, controlled_rotation_times_x.axis
         )
         if abs((theta0_with_x - theta2_with_x) % (2 * math.pi)) < ATOL:
@@ -58,7 +58,7 @@ def decompose(g: Gate) -> [Gate]:
                 + [Rz(q=g.control_qubit, theta=Float(g.target_gate.phase - math.pi / 2))]
             )
 
-        theta0, theta1, theta2 = get_zyz_decomposition_angles(g.target_gate.angle, g.target_gate.axis)
+        theta0, theta1, theta2 = ZYZDecomposer().get_decomposition_angles(g.target_gate.angle, g.target_gate.axis)
 
         A = [Ry(q=target_qubit, theta=Float(theta1 / 2)), Rz(q=target_qubit, theta=Float(theta2))]
 

opensquirrel/decomposer/xyx_decomposer.py

@@ -1,97 +1,14 @@
-import math
-from typing import Tuple
+from __future__ import annotations
 
-from opensquirrel.common import ATOL
-from opensquirrel.decomposer.general_decomposer import Decomposer
-from opensquirrel.default_gates import Rx, Ry
-from opensquirrel.ir import BlochSphereRotation, Float, Gate
-from opensquirrel.utils.identity_filter import filter_out_identities
-
-
-def get_xyx_decomposition_angles(alpha: float, axis: Tuple[float, float, float]) -> Tuple[float, float, float]:
-    """
-    Gives the angles used in the X-Y-X decomposition of the Bloch sphere rotation
-    characterized by a rotation around `axis` of angle `alpha`.
-
-    Parameters:
-        alpha: angle of the Bloch sphere rotation
-        axis: _normalized_ axis of the Bloch sphere rotation
-
-    Returns:
-        A triple (theta1, theta2, theta3) corresponding to the decomposition of the
-        arbitrary Bloch sphere rotation into U = Rx(theta3) Ry(theta2) Rx(theta1)
-
-    """
-    nx, ny, nz = axis
-
-    assert abs(nx**2 + ny**2 + nz**2 - 1) < ATOL, "Axis needs to be normalized"
-
-    assert -math.pi + ATOL < alpha <= math.pi + ATOL, "Angle needs to be normalized"
-
-    if abs(alpha - math.pi) < ATOL:
-        # alpha == pi, math.tan(alpha / 2) is not defined.
-
-        if abs(nx) < ATOL:
-            theta2 = math.pi
-            p = 0
-            m = 2 * math.acos(ny)
-
-        else:
-            p = math.pi
-            theta2 = 2 * math.acos(nx)
-
-            if abs(nx - 1) < ATOL or abs(nx + 1) < ATOL:
-                m = p  # This can be anything, but setting m = p means theta3 == 0, which is better for gate count.
-            else:
-                m = 2 * math.acos(ny / math.sqrt(1 - nx**2))
-
-    else:
-        p = 2 * math.atan2(nx * math.sin(alpha / 2), math.cos(alpha / 2))
+from collections.abc import Iterable
 
-        acos_argument = math.cos(alpha / 2) * math.sqrt(1 + (nx * math.tan(alpha / 2)) ** 2)
-
-        # This fixes float approximations like 1.0000000000002 which acos doesn't like.
-        acos_argument = max(min(acos_argument, 1.0), -1.0)
-
-        theta2 = 2 * math.acos(acos_argument)
-        theta2 = math.copysign(theta2, alpha)
-
-        if abs(math.sin(theta2 / 2)) < ATOL:
-            m = p  # This can be anything, but setting m = p means theta3 == 0, which is better for gate count.
-        else:
-            acos_argument = ny * math.sin(alpha / 2) / math.sin(theta2 / 2)
-
-            # This fixes float approximations like 1.0000000000002 which acos doesn't like.
-            acos_argument = max(min(acos_argument, 1.0), -1.0)
-
-            m = 2 * math.acos(acos_argument)
-
-    theta1 = (p + m) / 2
-
-    theta3 = p - theta1
-
-    return theta1, theta2, theta3
-
-
-class XYXDecomposer(Decomposer):
-    @staticmethod
-    def decompose(g: Gate) -> [Gate]:
-        if not isinstance(g, BlochSphereRotation):
-            # Only decomposer single-qubit gates.
-            return [g]
-
-        theta1, theta2, theta3 = get_xyx_decomposition_angles(g.angle, g.axis)
-
-        x1 = Rx(g.qubit, Float(theta1))
-        y = Ry(g.qubit, Float(theta2))
-        x2 = Rx(g.qubit, Float(theta3))
+from opensquirrel.decomposer.aba_decomposer import ABADecomposer
+from opensquirrel.default_gates import Rx, Ry
 
-        # Note: written like this, the decomposition doesn't preserve the global phase, which is fine
-        # since the global phase is a physically irrelevant artifact of the mathematical
-        # model we use to describe the quantum system.
 
-        # Should we want to preserve it, we would need to use a raw BlochSphereRotation, which would then
-        # be an anonymous gate in the resulting decomposed circuit:
-        # z2 = BlochSphereRotation(qubit=g.qubit, angle=theta3, axis=(0, 0, 1), phase = g.phase)
+class XYXDecomposer(ABADecomposer):
+    def __init__(self):
+        ABADecomposer.__init__(self, Rx, Ry)
 
-        return filter_out_identities([x1, y, x2])
+    def get_decomposition_angles(self, alpha: float, axis: Iterable[float]) -> tuple[float, float, float]:
+        return ABADecomposer.get_decomposition_angles(alpha, axis)

opensquirrel/decomposer/zyz_decomposer.py

@@ -1,97 +1,14 @@
-import math
-from typing import Tuple
+from __future__ import annotations
 
-from opensquirrel.common import ATOL
-from opensquirrel.decomposer.general_decomposer import Decomposer
-from opensquirrel.default_gates import Ry, Rz
-from opensquirrel.ir import BlochSphereRotation, Float, Gate
-from opensquirrel.utils.identity_filter import filter_out_identities
-
-
-def get_zyz_decomposition_angles(alpha: float, axis: Tuple[float, float, float]) -> Tuple[float, float, float]:
-    """
-    Gives the angles used in the Z-Y-Z decomposition of the Bloch sphere rotation
-    characterized by a rotation around `axis` of angle `alpha`.
-
-    Parameters:
-        alpha: angle of the Bloch sphere rotation
-        axis: _normalized_ axis of the Bloch sphere rotation
-
-    Returns:
-        A triple (theta1, theta2, theta3) corresponding to the decomposition of the
-        arbitrary Bloch sphere rotation into U = Rz(theta3) Ry(theta2) Rz(theta1)
-
-    """
-    nx, ny, nz = axis
-
-    assert abs(nx**2 + ny**2 + nz**2 - 1) < ATOL, "Axis needs to be normalized"
-
-    assert -math.pi + ATOL < alpha <= math.pi + ATOL, "Angle needs to be normalized"
-
-    if abs(alpha - math.pi) < ATOL:
-        # alpha == pi, math.tan(alpha / 2) is not defined.
-
-        if abs(nz) < ATOL:
-            theta2 = math.pi
-            p = 0
-            m = 2 * math.acos(ny)
-
-        else:
-            p = math.pi
-            theta2 = 2 * math.acos(nz)
-
-            if abs(nz - 1) < ATOL or abs(nz + 1) < ATOL:
-                m = p  # This can be anything, but setting m = p means theta3 == 0, which is better for gate count.
-            else:
-                m = 2 * math.acos(ny / math.sqrt(1 - nz**2))
-
-    else:
-        p = 2 * math.atan2(nz * math.sin(alpha / 2), math.cos(alpha / 2))
+from collections.abc import Iterable
 
-        acos_argument = math.cos(alpha / 2) * math.sqrt(1 + (nz * math.tan(alpha / 2)) ** 2)
-
-        # This fixes float approximations like 1.0000000000002 which acos doesn't like.
-        acos_argument = max(min(acos_argument, 1.0), -1.0)
-
-        theta2 = 2 * math.acos(acos_argument)
-        theta2 = math.copysign(theta2, alpha)
-
-        if abs(math.sin(theta2 / 2)) < ATOL:
-            m = p  # This can be anything, but setting m = p means theta3 == 0, which is better for gate count.
-        else:
-            acos_argument = ny * math.sin(alpha / 2) / math.sin(theta2 / 2)
-
-            # This fixes float approximations like 1.0000000000002 which acos doesn't like.
-            acos_argument = max(min(acos_argument, 1.0), -1.0)
-
-            m = 2 * math.acos(acos_argument)
-
-    theta1 = (p + m) / 2
-
-    theta3 = p - theta1
-
-    return theta1, theta2, theta3
-
-
-class ZYZDecomposer(Decomposer):
-    @staticmethod
-    def decompose(g: Gate) -> [Gate]:
-        if not isinstance(g, BlochSphereRotation):
-            # Only decomposer single-qubit gates.
-            return [g]
-
-        theta1, theta2, theta3 = get_zyz_decomposition_angles(g.angle, g.axis)
-
-        z1 = Rz(g.qubit, Float(theta1))
-        y = Ry(g.qubit, Float(theta2))
-        z2 = Rz(g.qubit, Float(theta3))
+from opensquirrel.decomposer.aba_decomposer import ABADecomposer
+from opensquirrel.default_gates import Ry, Rz
 
-        # Note: written like this, the decomposition doesn't preserve the global phase, which is fine
-        # since the global phase is a physically irrelevant artifact of the mathematical
-        # model we use to describe the quantum system.
 
-        # Should we want to preserve it, we would need to use a raw BlochSphereRotation, which would then
-        # be an anonymous gate in the resulting decomposed circuit:
-        # z2 = BlochSphereRotation(qubit=g.qubit, angle=theta3, axis=(0, 0, 1), phase = g.phase)
+class ZYZDecomposer(ABADecomposer):
+    def __init__(self):
+        ABADecomposer.__init__(self, Rz, Ry)
 
-        return filter_out_identities([z1, y, z2])
+    def get_decomposition_angles(self, alpha: float, axis: Iterable[float]) -> tuple[float, float, float]:
+        return ABADecomposer.get_decomposition_angles(alpha, axis)

test/test_integration.py

@@ -361,7 +361,7 @@ def test_export_quantify_scheduler(self):
             """
         )
 
-        circuit.decompose(decomposer=CNOTDecomposer)
+        circuit.decompose(decomposer=CNOTDecomposer())
 
         # Quantify-scheduler prefers CZ.
         circuit.replace(
@@ -378,7 +378,7 @@ def test_export_quantify_scheduler(self):
 
         # FIXME: for best gate count we need a Z-XY decomposer.
         # See https://github.com/QuTech-Delft/OpenSquirrel/issues/98
-        circuit.decompose(decomposer=ZYZDecomposer)
+        circuit.decompose(decomposer=ZYZDecomposer())
 
         if importlib.util.find_spec("quantify_scheduler") is None:
             with self.assertRaisesRegex(