Use singleton matrices for unparametrised standard gates (Qiskit#10296)

* Use singleton matrices for unparametrised standard gates

This makes the array form of standard gates with zero parameters
singleton class attributes that reject modification. The class-level
`__array__` methods are updated to return exactly the same instance,
except in very unusual circumstances, which means that
`Gate.to_matrix()` and `numpy.asarray()` calls on the objects will
return the same instance. This avoids a decent amount of construction
time, and avoids several Python-space list allocations and array
allocations.

The dtypes of the static arrays are all standardised to by complex128.
Gate matrices are in general unitary, `Gate.to_matrix()` already
enforces a cast to `complex128`.  For gates that allowed their dtypes to
be inferred, there were several cases where native ints and floats would
be used, meaning that `Gate.to_matrix()` would also involve an extra
matrix allocation to hold the cast, which just wasted time.

For standard controlled gates, we store both the closed- and
open-controlled matrices singly controlled gates.  For gates with more
than one control, we only store the "all ones" controlled case, as a
memory/speed trade-off; open controls are much less common than closed
controls.

For the most part this won't have an effect on peak memory usage, since
all the allocated matrices in standard Qiskit usage would be freed by
the garbage collector almost immediately.  This will, however, reduce
construction costs and garbage-collector pressure, since fewer
allocations+frees will occur, and no calculations will need to be done.

* Store only all-ones controls for large matrices

* Fix lint

* Use metaprogramming decorator to make `__array__` methods

Instead of defining the array functions manually for each class, this
adds a small amount of metaprogramming that adds them in with the
correct `ndarray` properties set, including for controlled gates.

Author: jakelishman

qiskit/circuit/_utils.py

@@ -95,3 +95,63 @@ def _ctrl_state_to_int(ctrl_state, num_ctrl_qubits):
     else:
         raise CircuitError(f"invalid control state specification: {repr(ctrl_state)}")
     return ctrl_state_std
+
+
+def with_gate_array(base_array):
+    """Class decorator that adds an ``__array__`` method to a :class:`.Gate` instance that returns a
+    singleton nonwritable view onto the complex matrix described by ``base_array``."""
+    nonwritable = numpy.array(base_array, dtype=numpy.complex128)
+    nonwritable.setflags(write=False)
+
+    def __array__(_self, dtype=None):
+        return numpy.asarray(nonwritable, dtype=dtype)
+
+    def decorator(cls):
+        if hasattr(cls, "__array__"):
+            raise RuntimeError("Refusing to decorate a class that already has '__array__' defined.")
+        cls.__array__ = __array__
+        return cls
+
+    return decorator
+
+
+def with_controlled_gate_array(base_array, num_ctrl_qubits, cached_states=None):
+    """Class decorator that adds an ``__array__`` method to a :class:`.ControlledGate` instance that
+    returns singleton nonwritable views onto a relevant precomputed complex matrix for the given
+    control state.
+
+    If ``cached_states`` is not given, then all possible control states are precomputed.  If it is
+    given, it should be an iterable of integers, and only these control states will be cached."""
+    base = numpy.asarray(base_array, dtype=numpy.complex128)
+
+    def matrix_for_control_state(state):
+        out = numpy.asarray(
+            _compute_control_matrix(base, num_ctrl_qubits, state),
+            dtype=numpy.complex128,
+        )
+        out.setflags(write=False)
+        return out
+
+    if cached_states is None:
+        nonwritables = [matrix_for_control_state(state) for state in range(2**num_ctrl_qubits)]
+
+        def __array__(self, dtype=None):
+            return numpy.asarray(nonwritables[self.ctrl_state], dtype=dtype)
+
+    else:
+        nonwritables = {state: matrix_for_control_state(state) for state in cached_states}
+
+        def __array__(self, dtype=None):
+            if (out := nonwritables.get(self.ctrl_state)) is not None:
+                return numpy.asarray(out, dtype=dtype)
+            return numpy.asarray(
+                _compute_control_matrix(base, num_ctrl_qubits, self.ctrl_state), dtype=dtype
+            )
+
+    def decorator(cls):
+        if hasattr(cls, "__array__"):
+            raise RuntimeError("Refusing to decorate a class that already has '__array__' defined.")
+        cls.__array__ = __array__
+        return cls
+
+    return decorator

qiskit/circuit/library/standard_gates/dcx.py

@@ -12,11 +12,12 @@
 
 """Double-CNOT gate."""
 
-import numpy as np
 from qiskit.circuit.gate import Gate
 from qiskit.circuit.quantumregister import QuantumRegister
+from qiskit.circuit._utils import with_gate_array
 
 
+@with_gate_array([[1, 0, 0, 0], [0, 0, 0, 1], [0, 1, 0, 0], [0, 0, 1, 0]])
 class DCXGate(Gate):
     r"""Double-CNOT gate.
 
@@ -66,7 +67,3 @@ def _define(self):
             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, 0, 0, 1], [0, 1, 0, 0], [0, 0, 1, 0]], dtype=dtype)

qiskit/circuit/library/standard_gates/ecr.py

@@ -14,12 +14,16 @@
 from math import sqrt
 import numpy as np
 
+from qiskit.circuit._utils import with_gate_array
 from qiskit.circuit.gate import Gate
 from qiskit.circuit.quantumregister import QuantumRegister
 from .rzx import RZXGate
 from .x import XGate
 
 
+@with_gate_array(
+    sqrt(0.5) * np.array([[0, 1, 0, 1.0j], [1, 0, -1.0j, 0], [0, 1.0j, 0, 1], [-1.0j, 0, 1, 0]])
+)
 class ECRGate(Gate):
     r"""An echoed cross-resonance gate.
 
@@ -106,14 +110,3 @@ def _define(self):
     def inverse(self):
         """Return inverse ECR gate (itself)."""
         return ECRGate()  # self-inverse
-
-    def to_matrix(self):
-        """Return a numpy.array for the ECR gate."""
-        return (
-            1
-            / sqrt(2)
-            * np.array(
-                [[0, 1, 0, 1.0j], [1, 0, -1.0j, 0], [0, 1.0j, 0, 1], [-1.0j, 0, 1, 0]],
-                dtype=complex,
-            )
-        )

qiskit/circuit/library/standard_gates/h.py

@@ -17,10 +17,14 @@
 from qiskit.circuit.controlledgate import ControlledGate
 from qiskit.circuit.gate import Gate
 from qiskit.circuit.quantumregister import QuantumRegister
+from qiskit.circuit._utils import with_gate_array, with_controlled_gate_array
 from .t import TGate, TdgGate
 from .s import SGate, SdgGate
 
+_H_ARRAY = 1 / sqrt(2) * numpy.array([[1, 1], [1, -1]], dtype=numpy.complex128)
 
+
+@with_gate_array(_H_ARRAY)
 class HGate(Gate):
     r"""Single-qubit Hadamard gate.
 
@@ -99,11 +103,8 @@ def inverse(self):
         r"""Return inverted H gate (itself)."""
         return HGate()  # self-inverse
 
-    def __array__(self, dtype=None):
-        """Return a Numpy.array for the H gate."""
-        return numpy.array([[1, 1], [1, -1]], dtype=dtype) / numpy.sqrt(2)
-
 
+@with_controlled_gate_array(_H_ARRAY, num_ctrl_qubits=1)
 class CHGate(ControlledGate):
     r"""Controlled-Hadamard gate.
 
@@ -160,16 +161,6 @@ class CHGate(ControlledGate):
                     0 & 0 & \frac{1}{\sqrt{2}} & -\frac{1}{\sqrt{2}}
                 \end{pmatrix}
     """
-    # Define class constants. This saves future allocation time.
-    _sqrt2o2 = 1 / sqrt(2)
-    _matrix1 = numpy.array(
-        [[1, 0, 0, 0], [0, _sqrt2o2, 0, _sqrt2o2], [0, 0, 1, 0], [0, _sqrt2o2, 0, -_sqrt2o2]],
-        dtype=complex,
-    )
-    _matrix0 = numpy.array(
-        [[_sqrt2o2, 0, _sqrt2o2, 0], [0, 1, 0, 0], [_sqrt2o2, 0, -_sqrt2o2, 0], [0, 0, 0, 1]],
-        dtype=complex,
-    )
 
     def __init__(self, label: Optional[str] = None, ctrl_state: Optional[Union[int, str]] = None):
         """Create new CH gate."""
@@ -212,10 +203,3 @@ def _define(self):
     def inverse(self):
         """Return inverted CH gate (itself)."""
         return CHGate(ctrl_state=self.ctrl_state)  # self-inverse
-
-    def __array__(self, dtype=None):
-        """Return a numpy.array for the CH gate."""
-        mat = self._matrix1 if self.ctrl_state else self._matrix0
-        if dtype:
-            return numpy.asarray(mat, dtype=dtype)
-        return mat

qiskit/circuit/library/standard_gates/i.py

@@ -13,10 +13,11 @@
 """Identity gate."""
 
 from typing import Optional
-import numpy
 from qiskit.circuit.gate import Gate
+from qiskit.circuit._utils import with_gate_array
 
 
+@with_gate_array([[1, 0], [0, 1]])
 class IGate(Gate):
     r"""Identity gate.
 
@@ -52,10 +53,6 @@ def inverse(self):
         """Invert this gate."""
         return IGate()  # self-inverse
 
-    def __array__(self, dtype=None):
-        """Return a numpy.array for the identity gate."""
-        return numpy.array([[1, 0], [0, 1]], dtype=dtype)
-
     def power(self, exponent: float):
         """Raise gate to a power."""
         return IGate()

qiskit/circuit/library/standard_gates/iswap.py

@@ -18,10 +18,12 @@
 
 from qiskit.circuit.gate import Gate
 from qiskit.circuit.quantumregister import QuantumRegister
+from qiskit.circuit._utils import with_gate_array
 
 from .xx_plus_yy import XXPlusYYGate
 
 
+@with_gate_array([[1, 0, 0, 0], [0, 0, 1j, 0], [0, 1j, 0, 0], [0, 0, 0, 1]])
 class iSwapGate(Gate):
     r"""iSWAP gate.
 
@@ -120,10 +122,6 @@ def _define(self):
 
         self.definition = qc
 
-    def __array__(self, dtype=None):
-        """Return a numpy.array for the iSWAP gate."""
-        return np.array([[1, 0, 0, 0], [0, 0, 1j, 0], [0, 1j, 0, 0], [0, 0, 0, 1]], dtype=dtype)
-
     def power(self, exponent: float):
         """Raise gate to a power."""
         return XXPlusYYGate(-np.pi * exponent)

qiskit/circuit/library/standard_gates/s.py

@@ -21,8 +21,14 @@
 from qiskit.circuit.gate import Gate
 from qiskit.circuit.library.standard_gates.p import CPhaseGate, PhaseGate
 from qiskit.circuit.quantumregister import QuantumRegister
+from qiskit.circuit._utils import with_gate_array, with_controlled_gate_array
 
 
+_S_ARRAY = numpy.array([[1, 0], [0, 1j]])
+_SDG_ARRAY = numpy.array([[1, 0], [0, -1j]])
+
+
+@with_gate_array(_S_ARRAY)
 class SGate(Gate):
     r"""Single qubit S gate (Z**0.5).
 
@@ -78,15 +84,12 @@ def inverse(self):
         """Return inverse of S (SdgGate)."""
         return SdgGate()
 
-    def __array__(self, dtype=None):
-        """Return a numpy.array for the S gate."""
-        return numpy.array([[1, 0], [0, 1j]], dtype=dtype)
-
     def power(self, exponent: float):
         """Raise gate to a power."""
         return PhaseGate(0.5 * numpy.pi * exponent)
 
 
+@with_gate_array(_SDG_ARRAY)
 class SdgGate(Gate):
     r"""Single qubit S-adjoint gate (~Z**0.5).
 
@@ -142,15 +145,12 @@ def inverse(self):
         """Return inverse of Sdg (SGate)."""
         return SGate()
 
-    def __array__(self, dtype=None):
-        """Return a numpy.array for the Sdg gate."""
-        return numpy.array([[1, 0], [0, -1j]], dtype=dtype)
-
     def power(self, exponent: float):
         """Raise gate to a power."""
         return PhaseGate(-0.5 * numpy.pi * exponent)
 
 
+@with_controlled_gate_array(_S_ARRAY, num_ctrl_qubits=1)
 class CSGate(ControlledGate):
     r"""Controlled-S gate.
 
@@ -179,23 +179,6 @@ class CSGate(ControlledGate):
                 0 & 0 & 0 & i
             \end{pmatrix}
     """
-    # Define class constants. This saves future allocation time.
-    _matrix1 = numpy.array(
-        [
-            [1, 0, 0, 0],
-            [0, 1, 0, 0],
-            [0, 0, 1, 0],
-            [0, 0, 0, 1j],
-        ]
-    )
-    _matrix0 = numpy.array(
-        [
-            [1, 0, 0, 0],
-            [0, 1, 0, 0],
-            [0, 0, 1j, 0],
-            [0, 0, 0, 1],
-        ]
-    )
 
     def __init__(self, label: Optional[str] = None, ctrl_state: Optional[Union[str, int]] = None):
         """Create new CS gate."""
@@ -213,18 +196,12 @@ def inverse(self):
         """Return inverse of CSGate (CSdgGate)."""
         return CSdgGate(ctrl_state=self.ctrl_state)
 
-    def __array__(self, dtype=None):
-        """Return a numpy.array for the CS gate."""
-        mat = self._matrix1 if self.ctrl_state == 1 else self._matrix0
-        if dtype is not None:
-            return numpy.asarray(mat, dtype=dtype)
-        return mat
-
     def power(self, exponent: float):
         """Raise gate to a power."""
         return CPhaseGate(0.5 * numpy.pi * exponent)
 
 
+@with_controlled_gate_array(_SDG_ARRAY, num_ctrl_qubits=1)
 class CSdgGate(ControlledGate):
     r"""Controlled-S^\dagger gate.
 
@@ -253,23 +230,6 @@ class CSdgGate(ControlledGate):
                 0 & 0 & 0 & -i
             \end{pmatrix}
     """
-    # Define class constants. This saves future allocation time.
-    _matrix1 = numpy.array(
-        [
-            [1, 0, 0, 0],
-            [0, 1, 0, 0],
-            [0, 0, 1, 0],
-            [0, 0, 0, -1j],
-        ]
-    )
-    _matrix0 = numpy.array(
-        [
-            [1, 0, 0, 0],
-            [0, 1, 0, 0],
-            [0, 0, -1j, 0],
-            [0, 0, 0, 1],
-        ]
-    )
 
     def __init__(self, label: Optional[str] = None, ctrl_state: Optional[Union[str, int]] = None):
         """Create new CSdg gate."""
@@ -293,13 +253,6 @@ def inverse(self):
         """Return inverse of CSdgGate (CSGate)."""
         return CSGate(ctrl_state=self.ctrl_state)
 
-    def __array__(self, dtype=None):
-        """Return a numpy.array for the CSdg gate."""
-        mat = self._matrix1 if self.ctrl_state == 1 else self._matrix0
-        if dtype is not None:
-            return numpy.asarray(mat, dtype=dtype)
-        return mat
-
     def power(self, exponent: float):
         """Raise gate to a power."""
         return CPhaseGate(-0.5 * numpy.pi * exponent)