diff --git a/qiskit/circuit/instruction.py b/qiskit/circuit/instruction.py index d6930fbe859e..934a237dfcc6 100644 --- a/qiskit/circuit/instruction.py +++ b/qiskit/circuit/instruction.py @@ -117,14 +117,14 @@ def __eq__(self, other): try: if numpy.shape(self_param) == numpy.shape(other_param) \ and numpy.allclose(self_param, other_param, - atol=_CUTOFF_PRECISION): + atol=_CUTOFF_PRECISION, rtol=0): continue except TypeError: pass try: if numpy.isclose(float(self_param), float(other_param), - atol=_CUTOFF_PRECISION): + atol=_CUTOFF_PRECISION, rtol=0): continue except TypeError: pass diff --git a/qiskit/circuit/quantumcircuit.py b/qiskit/circuit/quantumcircuit.py index 2a274b90dd69..7bfa5af82922 100644 --- a/qiskit/circuit/quantumcircuit.py +++ b/qiskit/circuit/quantumcircuit.py @@ -1926,12 +1926,10 @@ def global_phase(self, angle): if isinstance(angle, ParameterExpression) and angle.parameters: self._global_phase = angle else: - # Set the phase to the [-2 * pi, 2 * pi] interval + # Set the phase to the [0, 2π) interval angle = float(angle) if not angle: self._global_phase = 0 - elif angle < 0: - self._global_phase = angle % (-2 * np.pi) else: self._global_phase = angle % (2 * np.pi) diff --git a/qiskit/dagcircuit/dagcircuit.py b/qiskit/dagcircuit/dagcircuit.py index a17f4cdb7f36..7ee8766f7110 100644 --- a/qiskit/dagcircuit/dagcircuit.py +++ b/qiskit/dagcircuit/dagcircuit.py @@ -171,12 +171,10 @@ def global_phase(self, angle): if isinstance(angle, ParameterExpression): self._global_phase = angle else: - # Set the phase to the [-2 * pi, 2 * pi] interval + # Set the phase to the [0, 2π) interval angle = float(angle) if not angle: self._global_phase = 0 - elif angle < 0: - self._global_phase = angle % (-2 * math.pi) else: self._global_phase = angle % (2 * math.pi) @@ -888,10 +886,11 @@ def __eq__(self, other): # Try to convert to float, but in case of unbound ParameterExpressions # a TypeError will be raise, fallback to normal equality in those # cases + try: self_phase = float(self.global_phase) other_phase = float(other.global_phase) - if not np.isclose(self_phase, other_phase): + if abs((self_phase - other_phase + np.pi) % (2*np.pi) - np.pi) > 1.E-10: # TODO: atol? return False except TypeError: if self.global_phase != other.global_phase: @@ -915,7 +914,6 @@ def __eq__(self, other): for regname, reg in other.qregs.items()] other_creg_indices = [(regname, [other_bit_indices[bit] for bit in reg]) for regname, reg in other.cregs.items()] - if ( self_qreg_indices != other_qreg_indices or self_creg_indices != other_creg_indices diff --git a/qiskit/dagcircuit/dagdependency.py b/qiskit/dagcircuit/dagdependency.py index 407bb669112c..636af3f1ca05 100644 --- a/qiskit/dagcircuit/dagdependency.py +++ b/qiskit/dagcircuit/dagdependency.py @@ -107,12 +107,10 @@ def global_phase(self, angle): if isinstance(angle, ParameterExpression): self._global_phase = angle else: - # Set the phase to the [-2 * pi, 2 * pi] interval + # Set the phase to the [0, 2π) interval angle = float(angle) if not angle: self._global_phase = 0 - elif angle < 0: - self._global_phase = angle % (-2 * math.pi) else: self._global_phase = angle % (2 * math.pi) diff --git a/qiskit/quantum_info/synthesis/one_qubit_decompose.py b/qiskit/quantum_info/synthesis/one_qubit_decompose.py index 5aaf8a4ba0c4..6759217363fc 100644 --- a/qiskit/quantum_info/synthesis/one_qubit_decompose.py +++ b/qiskit/quantum_info/synthesis/one_qubit_decompose.py @@ -15,6 +15,7 @@ """ import math +import cmath import numpy as np import scipy.linalg as la @@ -22,7 +23,8 @@ from qiskit.circuit.quantumregister import QuantumRegister from qiskit.circuit.library.standard_gates import (UGate, PhaseGate, U3Gate, U2Gate, U1Gate, RXGate, RYGate, - RZGate, RGate, SXGate) + RZGate, RGate, SXGate, + XGate) from qiskit.exceptions import QiskitError from qiskit.quantum_info.operators.predicates import is_unitary_matrix @@ -38,6 +40,7 @@ 'ZYZ': ['rz', 'ry'], 'ZXZ': ['rz', 'rx'], 'XYX': ['rx', 'ry'], + 'ZSXX': ['rz', 'sx', 'x'], 'ZSX': ['rz', 'sx'], } @@ -83,8 +86,13 @@ class OneQubitEulerDecomposer: :math:`U_1(\theta+\pi).R_X\left(\frac{\pi}{2}\right).U_1(\lambda)` * - 'ZSX' - :math:`Z(\phi) Y(\theta) Z(\lambda)` - - :math:`e^{i\gamma} U_1(\phi+\pi).R_X\left(\frac{\pi}{2}\right).` - :math:`R_Z(\theta+\pi).S_X\left(\frac{\pi}{2}\right).U_1(\lambda)` + - :math:`e^{i\gamma} R_Z(\phi+\pi).\sqrt{X}.` + :math:`R_Z(\theta+\pi).\sqrt{X}.R_Z(\lambda)` + * - 'ZSXX' + - :math:`Z(\phi) Y(\theta) Z(\lambda)` + - :math:`e^{i\gamma} R_Z(\phi+\pi).\sqrt{X}.R_Z(\theta+\pi).\sqrt{X}.R_Z(\lambda)` + or + :math:`e^{i\gamma} R_Z(\phi+\pi).X.R_Z(\lambda)` * - 'U1X' - :math:`Z(\phi) Y(\theta) Z(\lambda)` - :math:`e^{i\gamma} U_1(\phi+\pi).R_X\left(\frac{\pi}{2}\right).` @@ -98,7 +106,8 @@ class OneQubitEulerDecomposer: def __init__(self, basis='U3'): """Initialize decomposer - Supported bases are: 'U', 'PSX', 'ZSX', 'U321', 'U3', 'U1X', 'RR', 'ZYZ', 'ZXZ', 'XYX'. + Supported bases are: 'U', 'PSX', 'ZSXX', 'ZSX', 'U321', 'U3', 'U1X', 'RR', 'ZYZ', 'ZXZ', + 'XYX'. Args: basis (str): the decomposition basis [Default: 'U3'] @@ -167,6 +176,7 @@ def basis(self, basis): 'U': (self._params_u3, self._circuit_u), 'PSX': (self._params_u1x, self._circuit_psx), 'ZSX': (self._params_u1x, self._circuit_zsx), + 'ZSXX': (self._params_u1x, self._circuit_zsxx), 'U1X': (self._params_u1x, self._circuit_u1x), 'RR': (self._params_zyz, self._circuit_rr), 'ZYZ': (self._params_zyz, self._circuit_zyz), @@ -207,7 +217,7 @@ def _params_zyz(mat): # We rescale the input matrix to be special unitary (det(U) = 1) # This ensures that the quaternion representation is real coeff = la.det(mat)**(-0.5) - phase = -np.angle(coeff) + phase = -cmath.phase(coeff) su_mat = coeff * mat # U in SU(2) # OpenQASM SU(2) parameterization: # U[0, 0] = exp(-i(phi+lambda)/2) * cos(theta/2) @@ -215,10 +225,10 @@ def _params_zyz(mat): # U[1, 0] = exp(i(phi-lambda)/2) * sin(theta/2) # U[1, 1] = exp(i(phi+lambda)/2) * cos(theta/2) theta = 2 * math.atan2(abs(su_mat[1, 0]), abs(su_mat[0, 0])) - phiplambda = 2 * np.angle(su_mat[1, 1]) - phimlambda = 2 * np.angle(su_mat[1, 0]) - phi = (phiplambda + phimlambda) / 2.0 - lam = (phiplambda - phimlambda) / 2.0 + phiplambda2 = cmath.phase(su_mat[1, 1]) + phimlambda2 = cmath.phase(su_mat[1, 0]) + phi = (phiplambda2 + phimlambda2) + lam = (phiplambda2 - phimlambda2) return theta, phi, lam, phase @staticmethod @@ -243,7 +253,8 @@ def _params_xyx(mat): ]], dtype=complex) theta, phi, lam, phase = OneQubitEulerDecomposer._params_zyz(mat_zyz) - return -theta, phi, lam, phase + newphi, newlam = _mod_2pi(phi+np.pi), _mod_2pi(lam+np.pi) + return theta, newphi, newlam, phase + (newphi + newlam - phi - lam)/2 @staticmethod def _params_u3(mat): @@ -271,17 +282,31 @@ def _circuit_zyz(theta, phase, simplify=True, atol=DEFAULT_ATOL): + gphase = phase - (phi+lam)/2 qr = QuantumRegister(1, 'qr') - circuit = QuantumCircuit(qr, global_phase=phase) - if simplify and math.isclose(theta, 0.0, abs_tol=atol): - circuit._append(RZGate(phi + lam), [qr[0]], []) + circuit = QuantumCircuit(qr) + if not simplify: + atol = -1.0 + if abs(theta) < atol: + tot = _mod_2pi(phi + lam, atol) + if abs(tot) > atol: + circuit._append(RZGate(tot), [qr[0]], []) + gphase += tot/2 + circuit.global_phase = gphase return circuit - if not simplify or not math.isclose(lam, 0.0, abs_tol=atol): + if abs(theta - np.pi) < atol: + gphase += phi + lam, phi = lam-phi, 0 + lam = _mod_2pi(lam, atol) + if abs(lam) > atol: + gphase += lam/2 circuit._append(RZGate(lam), [qr[0]], []) - if not simplify or not math.isclose(theta, 0.0, abs_tol=atol): - circuit._append(RYGate(theta), [qr[0]], []) - if not simplify or not math.isclose(phi, 0.0, abs_tol=atol): + circuit._append(RYGate(theta), [qr[0]], []) + phi = _mod_2pi(phi, atol) + if abs(phi) > atol: + gphase += phi/2 circuit._append(RZGate(phi), [qr[0]], []) + circuit.global_phase = gphase return circuit @staticmethod @@ -291,17 +316,31 @@ def _circuit_zxz(theta, phase, simplify=True, atol=DEFAULT_ATOL): + gphase = phase - (phi+lam)/2 qr = QuantumRegister(1, 'qr') - circuit = QuantumCircuit(qr, global_phase=phase) - if simplify and math.isclose(theta, 0.0, abs_tol=atol): - circuit._append(RZGate(phi + lam), [qr[0]], []) + circuit = QuantumCircuit(qr) + if not simplify: + atol = -1.0 + if abs(theta) < atol: + tot = _mod_2pi(phi + lam) + if abs(tot) > atol: + circuit._append(RZGate(tot), [qr[0]], []) + gphase += tot/2 + circuit.global_phase = gphase return circuit - if not simplify or not math.isclose(lam, 0.0, abs_tol=atol): + if abs(theta - np.pi) < atol: + gphase += phi + lam, phi = lam - phi, 0 + lam = _mod_2pi(lam, atol) + if abs(lam) > atol: + gphase += lam/2 circuit._append(RZGate(lam), [qr[0]], []) - if not simplify or not math.isclose(theta, 0.0, abs_tol=atol): - circuit._append(RXGate(theta), [qr[0]], []) - if not simplify or not math.isclose(phi, 0.0, abs_tol=atol): + circuit._append(RXGate(theta), [qr[0]], []) + phi = _mod_2pi(phi, atol) + if abs(phi) > atol: + gphase += phi/2 circuit._append(RZGate(phi), [qr[0]], []) + circuit.global_phase = gphase return circuit @staticmethod @@ -311,17 +350,31 @@ def _circuit_xyx(theta, phase, simplify=True, atol=DEFAULT_ATOL): + gphase = phase - (phi+lam)/2 qr = QuantumRegister(1, 'qr') - circuit = QuantumCircuit(qr, global_phase=phase) - if simplify and math.isclose(theta, 0.0, abs_tol=atol): - circuit._append(RXGate(phi + lam), [qr[0]], []) + circuit = QuantumCircuit(qr) + if not simplify: + atol = -1.0 + if abs(theta) < atol: + tot = _mod_2pi(phi + lam, atol) + if abs(tot) > atol: + circuit._append(RXGate(tot), [qr[0]], []) + gphase += tot/2 + circuit.global_phase = gphase return circuit - if not simplify or not math.isclose(lam, 0.0, abs_tol=atol): + if abs(theta - np.pi) < atol: + gphase += phi + lam, phi = lam-phi, 0 + lam = _mod_2pi(lam, atol) + if abs(lam) > atol: + gphase += lam/2 circuit._append(RXGate(lam), [qr[0]], []) - if not simplify or not math.isclose(theta, 0.0, abs_tol=atol): - circuit._append(RYGate(theta), [qr[0]], []) - if not simplify or not math.isclose(phi, 0.0, abs_tol=atol): + circuit._append(RYGate(theta), [qr[0]], []) + phi = _mod_2pi(phi, atol) + if abs(phi) > atol: + gphase += phi/2 circuit._append(RXGate(phi), [qr[0]], []) + circuit.global_phase = gphase return circuit @staticmethod @@ -331,10 +384,12 @@ def _circuit_u3(theta, phase, simplify=True, atol=DEFAULT_ATOL): - # pylint: disable=unused-argument qr = QuantumRegister(1, 'qr') circuit = QuantumCircuit(qr, global_phase=phase) - circuit._append(U3Gate(theta, phi, lam), [qr[0]], []) + phi = _mod_2pi(phi, atol) + lam = _mod_2pi(lam, atol) + if not simplify or abs(theta) > atol or abs(phi) > atol or abs(lam) > atol: + circuit._append(U3Gate(theta, phi, lam), [qr[0]], []) return circuit @staticmethod @@ -344,18 +399,18 @@ def _circuit_u321(theta, phase, simplify=True, atol=DEFAULT_ATOL): - rtol = 1e-9 # default is 1e-5, too far from atol=1e-12 qr = QuantumRegister(1, 'qr') circuit = QuantumCircuit(qr, global_phase=phase) - if simplify and (math.isclose(theta, 0.0, abs_tol=atol, rel_tol=rtol)): - phi_lam = phi + lam - if not (math.isclose(phi_lam, 0.0, abs_tol=atol, rel_tol=rtol) or - math.isclose(phi_lam, 2*np.pi, abs_tol=atol, rel_tol=rtol)): - circuit._append(U1Gate(_mod2pi(phi+lam)), [qr[0]], []) - elif simplify and math.isclose(theta, np.pi/2, abs_tol=atol, rel_tol=rtol): - circuit._append(U2Gate(phi, lam), [qr[0]], []) + if not simplify: + atol = -1.0 + if abs(theta) < atol: + tot = _mod_2pi(phi + lam, atol) + if abs(tot) > atol: + circuit._append(U1Gate(tot), [qr[0]], []) + elif abs(theta - np.pi/2) < atol: + circuit._append(U2Gate(_mod_2pi(phi, atol), _mod_2pi(lam, atol)), [qr[0]], []) else: - circuit._append(U3Gate(theta, phi, lam), [qr[0]], []) + circuit._append(U3Gate(theta, _mod_2pi(phi, atol), _mod_2pi(lam, atol)), [qr[0]], []) return circuit @staticmethod @@ -365,10 +420,48 @@ def _circuit_u(theta, phase, simplify=True, atol=DEFAULT_ATOL): - # pylint: disable=unused-argument qr = QuantumRegister(1, 'qr') circuit = QuantumCircuit(qr, global_phase=phase) - circuit._append(UGate(theta, phi, lam), [qr[0]], []) + if not simplify: + atol = -1.0 + phi = _mod_2pi(phi, atol) + lam = _mod_2pi(lam, atol) + if abs(theta) > atol or abs(phi) > atol or abs(lam) > atol: + circuit._append(UGate(theta, phi, lam), [qr[0]], []) + return circuit + + @staticmethod + def _circuit_psx_gen(theta, phi, lam, phase, atol, pfun, xfun, xpifun=None): + """Generic X90, phase decomposition""" + qr = QuantumRegister(1, 'qr') + circuit = QuantumCircuit(qr, global_phase=phase) + # Check for decomposition into minimimal number required SX pulses + if np.abs(theta) < atol: + # Zero SX gate decomposition + pfun(circuit, qr, lam + phi) + return circuit + if abs(theta - np.pi/2) < atol: + # Single SX gate decomposition + pfun(circuit, qr, lam - np.pi/2) + xfun(circuit, qr) + pfun(circuit, qr, phi + np.pi/2) + return circuit + # General two-SX gate decomposition + # Shift theta and phi so decomposition is + # P(phi).SX.P(theta).SX.P(lam) + if abs(theta-np.pi) < atol: + circuit.global_phase += lam + phi, lam = phi-lam, 0 + circuit.global_phase -= np.pi/2 + pfun(circuit, qr, lam) + if xpifun and abs(_mod_2pi(theta + np.pi)) < atol: + xpifun(circuit, qr) + else: + xfun(circuit, qr) + pfun(circuit, qr, theta + np.pi) + xfun(circuit, qr) + pfun(circuit, qr, phi + np.pi) + return circuit @staticmethod @@ -378,52 +471,18 @@ def _circuit_psx(theta, phase, simplify=True, atol=DEFAULT_ATOL): - # Shift theta and phi so decomposition is - # Phase(phi+pi).SX.Phase(theta+pi).SX.Phase(lam) - theta = _mod2pi(theta + np.pi) - phi = _mod2pi(phi + np.pi) - qr = QuantumRegister(1, 'qr') - circuit = QuantumCircuit(qr, global_phase=phase - np.pi / 2) - # Check for decomposition into minimimal number required SX gates - abs_theta = abs(theta) - if simplify and math.isclose(abs_theta, np.pi, abs_tol=atol): - lam_phi_theta = _mod2pi(lam + phi + theta) - abs_lam_phi_theta = _mod2pi(abs(lam + phi + theta)) - if not (math.isclose(abs_lam_phi_theta, 0., abs_tol=atol) or - math.isclose(abs_lam_phi_theta, 2*np.pi, abs_tol=atol)): - circuit._append(PhaseGate(lam_phi_theta), [qr[0]], []) - circuit.global_phase += np.pi / 2 - elif simplify and (math.isclose(abs_theta, np.pi/2, abs_tol=atol) or - math.isclose(abs_theta, 3*np.pi/2, abs_tol=atol)): - lam_theta = _mod2pi(lam + theta) - abs_lam_theta = _mod2pi(abs(lam + theta)) - if not (math.isclose(abs_lam_theta, 0, abs_tol=atol) or - math.isclose(abs_lam_theta, 2*np.pi, abs_tol=atol)): - circuit._append(PhaseGate(lam_theta), [qr[0]], []) - circuit._append(SXGate(), [qr[0]], []) - phi_theta = _mod2pi(phi + theta) - abs_phi_theta = _mod2pi(abs(phi_theta)) - if not (math.isclose(abs_phi_theta, 0, abs_tol=atol) or - math.isclose(abs_phi_theta, 2*np.pi, abs_tol=atol)): - circuit._append(PhaseGate(_mod2pi(phi + theta)), [qr[0]], []) - if (math.isclose(theta, -np.pi / 2, abs_tol=atol) or math.isclose( - theta, 3 * np.pi / 2, abs_tol=atol)): - circuit.global_phase += np.pi / 2 - else: - abs_lam = abs(lam) - if not (math.isclose(abs_lam, 0., abs_tol=atol) or - math.isclose(abs_lam, 2*np.pi, abs_tol=atol)): - circuit._append(PhaseGate(lam), [qr[0]], []) - circuit._append(SXGate(), [qr[0]], []) - if not (math.isclose(abs_theta, 0., abs_tol=atol) or - math.isclose(abs_theta, 2*np.pi, abs_tol=atol)): - circuit._append(PhaseGate(theta), [qr[0]], []) - circuit._append(SXGate(), [qr[0]], []) - abs_phi = abs(phi) - if not (math.isclose(abs_phi, 0., abs_tol=atol) or - math.isclose(abs_phi, 2*np.pi, abs_tol=atol)): + if not simplify: + atol = -1.0 + + def fnz(circuit, qr, phi): + phi = _mod_2pi(phi, atol) + if abs(phi) > atol: circuit._append(PhaseGate(phi), [qr[0]], []) - return circuit + + def fnx(circuit, qr): + circuit._append(SXGate(), [qr[0]], []) + + return OneQubitEulerDecomposer._circuit_psx_gen(theta, phi, lam, phase, atol, fnz, fnx) @staticmethod def _circuit_zsx(theta, @@ -432,58 +491,19 @@ def _circuit_zsx(theta, phase, simplify=True, atol=DEFAULT_ATOL): - # Shift theta and phi so decomposition is - # RZ(phi+pi).SX.RZ(theta+pi).SX.RZ(lam) - theta = _mod2pi(theta + np.pi) - phi = _mod2pi(phi + np.pi) - qr = QuantumRegister(1, 'qr') - circuit = QuantumCircuit(qr, global_phase=phase - np.pi / 2) - # Check for decomposition into minimimal number required SX gates - abs_theta = abs(theta) - if simplify and math.isclose(abs_theta, np.pi, abs_tol=atol): - lam_phi_theta = _mod2pi(lam + phi + theta) - abs_lam_phi_theta = _mod2pi(abs(lam + phi + theta)) - if not (math.isclose(abs_lam_phi_theta, 0., abs_tol=atol) or - math.isclose(abs_lam_phi_theta, 2*np.pi, abs_tol=atol)): - circuit._append(RZGate(lam_phi_theta), [qr[0]], []) - circuit.global_phase += 0.5 * lam_phi_theta - circuit.global_phase += np.pi / 2 - elif simplify and (math.isclose(abs_theta, np.pi/2, abs_tol=atol) or - math.isclose(abs_theta, 3*np.pi/2, abs_tol=atol)): - lam_theta = _mod2pi(lam + theta) - abs_lam_theta = _mod2pi(abs(lam + theta)) - if not (math.isclose(abs_lam_theta, 0, abs_tol=atol) or - math.isclose(abs_lam_theta, 2*np.pi, abs_tol=atol)): - circuit._append(RZGate(lam_theta), [qr[0]], []) - circuit.global_phase += 0.5 * lam_theta - circuit._append(SXGate(), [qr[0]], []) - phi_theta = _mod2pi(phi + theta) - abs_phi_theta = _mod2pi(abs(phi_theta)) - if not (math.isclose(abs_phi_theta, 0, abs_tol=atol) or - math.isclose(abs_phi_theta, 2*np.pi, abs_tol=atol)): - circuit._append(RZGate(phi_theta), [qr[0]], []) - circuit.global_phase += 0.5 * phi_theta - if (math.isclose(theta, -np.pi / 2, abs_tol=atol) or - math.isclose(theta, 3 * np.pi / 2, abs_tol=atol)): - circuit.global_phase += np.pi / 2 - else: - abs_lam = abs(lam) - if not (math.isclose(abs_lam, 0., abs_tol=atol) or - math.isclose(abs_lam, 2*np.pi, abs_tol=atol)): - circuit._append(RZGate(lam), [qr[0]], []) - circuit.global_phase += 0.5 * lam - circuit._append(SXGate(), [qr[0]], []) - if not (math.isclose(abs_theta, 0., abs_tol=atol) or - math.isclose(abs_theta, 2*np.pi, abs_tol=atol)): - circuit._append(RZGate(theta), [qr[0]], []) - circuit.global_phase += 0.5 * theta - circuit._append(SXGate(), [qr[0]], []) - abs_phi = abs(phi) - if not (math.isclose(abs_phi, 0., abs_tol=atol) or - math.isclose(abs_phi, 2*np.pi, abs_tol=atol)): + if not simplify: + atol = -1.0 + + def fnz(circuit, qr, phi): + phi = _mod_2pi(phi, atol) + if abs(phi) > atol: circuit._append(RZGate(phi), [qr[0]], []) - circuit.global_phase += 0.5 * phi - return circuit + circuit.global_phase += phi/2 + + def fnx(circuit, qr): + circuit._append(SXGate(), [qr[0]], []) + + return OneQubitEulerDecomposer._circuit_psx_gen(theta, phi, lam, phase, atol, fnz, fnx) @staticmethod def _circuit_u1x(theta, @@ -492,34 +512,44 @@ def _circuit_u1x(theta, phase, simplify=True, atol=DEFAULT_ATOL): - # Shift theta and phi so decomposition is - # U1(phi).X90.U1(theta).X90.U1(lam) - theta += np.pi - phi += np.pi - # Check for decomposition into minimimal number required X90 pulses - if simplify and math.isclose(abs(theta), np.pi, abs_tol=atol): - # Zero X90 gate decomposition - qr = QuantumRegister(1, 'qr') - circuit = QuantumCircuit(qr, global_phase=phase) - circuit._append(U1Gate(lam + phi + theta), [qr[0]], []) - return circuit - if simplify and math.isclose(abs(theta), np.pi/2, abs_tol=atol): - # Single X90 gate decomposition - qr = QuantumRegister(1, 'qr') - circuit = QuantumCircuit(qr, global_phase=phase) - circuit._append(U1Gate(lam + theta), [qr[0]], []) + if not simplify: + atol = -1.0 + + def fnz(circuit, qr, phi): + phi = _mod_2pi(phi, atol) + if abs(phi) > atol: + circuit._append(U1Gate(phi), [qr[0]], []) + + def fnx(circuit, qr): + circuit.global_phase += np.pi/4 circuit._append(RXGate(np.pi / 2), [qr[0]], []) - circuit._append(U1Gate(phi + theta), [qr[0]], []) - return circuit - # General two-X90 gate decomposition - qr = QuantumRegister(1, 'qr') - circuit = QuantumCircuit(qr, global_phase=phase) - circuit._append(U1Gate(lam), [qr[0]], []) - circuit._append(RXGate(np.pi / 2), [qr[0]], []) - circuit._append(U1Gate(theta), [qr[0]], []) - circuit._append(RXGate(np.pi / 2), [qr[0]], []) - circuit._append(U1Gate(phi), [qr[0]], []) - return circuit + + return OneQubitEulerDecomposer._circuit_psx_gen(theta, phi, lam, phase, atol, fnz, fnx) + + @staticmethod + def _circuit_zsxx(theta, + phi, + lam, + phase, + simplify=True, + atol=DEFAULT_ATOL): + if not simplify: + atol = -1.0 + + def fnz(circuit, qr, phi): + phi = _mod_2pi(phi, atol) + if abs(phi) > atol: + circuit._append(RZGate(phi), [qr[0]], []) + circuit.global_phase += phi/2 + + def fnx(circuit, qr): + circuit._append(SXGate(), [qr[0]], []) + + def fnxpi(circuit, qr): + circuit._append(XGate(), [qr[0]], []) + + return OneQubitEulerDecomposer._circuit_psx_gen(theta, phi, lam, phase, atol, + fnz, fnx, fnxpi) @staticmethod def _circuit_rr(theta, @@ -530,14 +560,19 @@ def _circuit_rr(theta, atol=DEFAULT_ATOL): qr = QuantumRegister(1, 'qr') circuit = QuantumCircuit(qr, global_phase=phase) - if not simplify or not math.isclose(theta, -np.pi, abs_tol=atol): - circuit._append(RGate(theta + np.pi, np.pi / 2 - lam), [qr[0]], []) - circuit._append(RGate(-np.pi, 0.5 * (phi - lam + np.pi)), [qr[0]], []) + if not simplify: + atol = -1.0 + if abs(theta) < atol and abs(phi) < atol and abs(lam) < atol: + return circuit + if abs(theta - np.pi) > atol: + circuit._append(RGate(theta - np.pi, _mod_2pi(np.pi / 2 - lam, atol)), [qr[0]], []) + circuit._append(RGate(np.pi, _mod_2pi(0.5 * (phi - lam + np.pi), atol)), [qr[0]], []) return circuit -def _mod2pi(angle): - if angle >= 0: - return math.fmod(angle, 2*np.pi) - else: - return math.fmod(angle, -2*np.pi) +def _mod_2pi(angle: float, atol: float = 0): + """Wrap angle into interval [-π,π). If within atol of the endpoint, clamp to -π""" + wrapped = (angle+np.pi) % (2*np.pi) - np.pi + if abs(wrapped - np.pi) < atol: + wrapped = -np.pi + return wrapped diff --git a/qiskit/quantum_info/synthesis/two_qubit_decompose.py b/qiskit/quantum_info/synthesis/two_qubit_decompose.py index 94acbbe25af1..6da1548b4967 100644 --- a/qiskit/quantum_info/synthesis/two_qubit_decompose.py +++ b/qiskit/quantum_info/synthesis/two_qubit_decompose.py @@ -23,29 +23,34 @@ Gambetta, J. M. Validating quantum computers using randomized model circuits. arXiv:1811.12926 [quant-ph] (2018). """ +import cmath import math +import io +import base64 import warnings +from typing import ClassVar, Optional + +import logging import numpy as np import scipy.linalg as la from qiskit.circuit.quantumregister import QuantumRegister from qiskit.circuit.quantumcircuit import QuantumCircuit -from qiskit.circuit.library.standard_gates.x import CXGate -from qiskit.circuit.tools import pi_check +from qiskit.circuit.library.standard_gates import CXGate, RXGate, RYGate, RZGate from qiskit.exceptions import QiskitError from qiskit.quantum_info.operators import Operator -from qiskit.quantum_info.operators.predicates import is_unitary_matrix from qiskit.quantum_info.synthesis.weyl import weyl_coordinates -from qiskit.quantum_info.synthesis.one_qubit_decompose import OneQubitEulerDecomposer +from qiskit.quantum_info.synthesis.one_qubit_decompose import OneQubitEulerDecomposer, DEFAULT_ATOL -_CUTOFF_PRECISION = 1e-12 +logger = logging.getLogger(__name__) def decompose_two_qubit_product_gate(special_unitary_matrix): """Decompose U = Ul⊗Ur where U in SU(4), and Ul, Ur in SU(2). Throws QiskitError if this isn't possible. """ + special_unitary_matrix = np.asarray(special_unitary_matrix, dtype=complex) # extract the right component R = special_unitary_matrix[:2, :2].copy() detR = R[0, 0]*R[1, 1] - R[0, 1]*R[1, 0] @@ -64,10 +69,10 @@ def decompose_two_qubit_product_gate(special_unitary_matrix): if abs(detL) < 0.9: raise QiskitError("decompose_two_qubit_product_gate: unable to decompose: detL < 0.9") L /= np.sqrt(detL) - phase = np.angle(detL) / 2 + phase = cmath.phase(detL) / 2 temp = np.kron(L, R) - deviation = np.abs(np.abs(temp.conj(temp).T.dot(special_unitary_matrix).trace()) - 4) + deviation = abs(abs(temp.conj().T.dot(special_unitary_matrix).trace()) - 4) if deviation > 1.E-13: raise QiskitError("decompose_two_qubit_product_gate: decomposition failed: " "deviation too large: {}".format(deviation)) @@ -86,21 +91,60 @@ def decompose_two_qubit_product_gate(special_unitary_matrix): [-1, 0]], dtype=complex) _ipz = np.array([[1j, 0], [0, -1j]], dtype=complex) +_id = np.array([[1, 0], + [0, 1]], dtype=complex) + +class TwoQubitWeylDecomposition(): + """Decompose two-qubit unitary U = (K1l⊗K1r).Exp(i a xx + i b yy + i c zz).(K2l⊗K2r) , where U ∈ + U(4), (K1l|K1r|K2l|K2r) ∈ SU(2), and we stay in the "Weyl Chamber" 𝜋/4 ≥ a ≥ b ≥ |c| -class TwoQubitWeylDecomposition: - """ Decompose two-qubit unitary U = (K1l⊗K1r).Exp(i a xx + i b yy + i c zz).(K2l⊗K2r) , - where U ∈ U(4), (K1l|K1r|K2l|K2r) ∈ SU(2), and we stay in the "Weyl Chamber" - 𝜋/4 ≥ a ≥ b ≥ |c| + This is an abstract factory class that instantiates itself as specialized subclasses based on + the fidelity, such that the approximation error from specialization has an average gate fidelity + at least as high as requested. The specialized subclasses have unique canonical representations + thus avoiding problems of numerical stability. + + Passing non-None fidelity to specializations is treated as an assertion, raising QiskitError if + forcing the specialization is more approximate than asserted. """ - def __init__(self, unitary_matrix, eps=1e-15): - """The flip into the Weyl Chamber is described in B. Kraus and J. I. Cirac, - Phys. Rev. A 63, 062309 (2001). + # The parameters of the decomposition: + a: float + b: float + c: float + global_phase: float + K1l: np.ndarray + K2l: np.ndarray + K1r: np.ndarray + K2r: np.ndarray + + unitary_matrix: np.ndarray # The unitary that was input + requested_fidelity: Optional[float] # None means no automatic specialization + calculated_fidelity: float # Fidelity after specialization + + _original_decomposition: "TwoQubitWeylDecomposition" + _is_flipped_from_original: bool # The approx is closest to a Weyl reflection of the original? + + _default_1q_basis: ClassVar[str] = 'ZYZ' # Default one qubit basis (explicit parameterization) + + def __init_subclass__(cls, **kwargs): + """Subclasses should be concrete, not factories. + + Make explicitly-instantiated subclass __new__ call base __new__ with fidelity=None""" + super().__init_subclass__(**kwargs) + cls.__new__ = (lambda cls, *a, fidelity=None, **k: + TwoQubitWeylDecomposition.__new__(cls, *a, fidelity=None, **k)) - FIXME: There's a cleaner-seeming method based on choosing branch cuts carefully, in - Andrew M. Childs, Henry L. Haselgrove, and Michael A. Nielsen, Phys. Rev. A 68, 052311, - but I wasn't able to get that to work. + @staticmethod + def __new__(cls, unitary_matrix, *, fidelity=(1.-1.E-9)): + """Perform the Weyl chamber decomposition, and optionally choose a specialized subclass. + + The flip into the Weyl Chamber is described in B. Kraus and J. I. Cirac, Phys. Rev. A 63, + 062309 (2001). + + FIXME: There's a cleaner-seeming method based on choosing branch cuts carefully, in Andrew + M. Childs, Henry L. Haselgrove, and Michael A. Nielsen, Phys. Rev. A 68, 052311, but I + wasn't able to get that to work. The overall decomposition scheme is taken from Drury and Love, arXiv:0806.4015 [quant-ph]. """ @@ -109,25 +153,23 @@ def __init__(self, unitary_matrix, eps=1e-15): pi4 = np.pi/4 # Make U be in SU(4) - U = unitary_matrix.copy() + U = np.array(unitary_matrix, dtype=complex, copy=True) detU = la.det(U) U *= detU**(-0.25) - global_phase = np.angle(detU) / 4 + global_phase = cmath.phase(detU) / 4 Up = _Bd.dot(U).dot(_B) M2 = Up.T.dot(Up) - M2.real[abs(M2.real) < eps] = 0.0 - M2.imag[abs(M2.imag) < eps] = 0.0 # M2 is a symmetric complex matrix. We need to decompose it as M2 = P D P^T where # P ∈ SO(4), D is diagonal with unit-magnitude elements. # D, P = la.eig(M2) # this can fail for certain kinds of degeneracy - for i in range(100): # FIXME: this randomized algorithm is horrendous - state = np.random.default_rng(i) + state = np.random.default_rng(2020) + for _ in range(100): # FIXME: this randomized algorithm is horrendous M2real = state.normal()*M2.real + state.normal()*M2.imag _, P = np.linalg.eigh(M2real) D = P.T.dot(M2).dot(P).diagonal() - if np.allclose(P.dot(np.diag(D)).dot(P.T), M2, rtol=1.0e-13, atol=1.0e-13): + if np.allclose(P.dot(np.diag(D)).dot(P.T), M2, rtol=0, atol=1.0E-13): break else: raise QiskitError("TwoQubitWeylDecomposition: failed to diagonalize M2") @@ -150,11 +192,7 @@ def __init__(self, unitary_matrix, eps=1e-15): # Find K1, K2 so that U = K1.A.K2, with K being product of single-qubit unitaries K1 = _B.dot(Up).dot(P).dot(np.diag(np.exp(1j*d))).dot(_Bd) - K1.real[abs(K1.real) < eps] = 0.0 - K1.imag[abs(K1.imag) < eps] = 0.0 K2 = _B.dot(P.T).dot(_Bd) - K2.real[abs(K2.real) < eps] = 0.0 - K2.imag[abs(K2.imag) < eps] = 0.0 K1l, K1r, phase_l = decompose_two_qubit_product_gate(K1) K2l, K2r, phase_r = decompose_two_qubit_product_gate(K2) @@ -205,39 +243,384 @@ def __init__(self, unitary_matrix, eps=1e-15): K1l = K1l.dot(_ipz) K1r = K1r.dot(_ipz) global_phase -= pi2 - self.a = cs[1] - self.b = cs[0] - self.c = cs[2] - self.K1l = K1l - self.K1r = K1r - self.K2l = K2l - self.K2r = K2r - self.global_phase = global_phase + + a, b, c = cs[1], cs[0], cs[2] + + # Save the non-specialized decomposition for later comparison + od = super().__new__(TwoQubitWeylDecomposition) + od.a = a + od.b = b + od.c = c + od.K1l = K1l + od.K1r = K1r + od.K2l = K2l + od.K2r = K2r + od.global_phase = global_phase + od.requested_fidelity = fidelity + od.calculated_fidelity = 1.0 + od.unitary_matrix = np.array(unitary_matrix, dtype=complex, copy=True) + od.unitary_matrix.setflags(write=False) + od._original_decomposition = None + od._is_flipped_from_original = False + + def is_close(ap, bp, cp): + da, db, dc = a-ap, b-bp, c-cp + tr = 4*complex(math.cos(da)*math.cos(db)*math.cos(dc), + math.sin(da)*math.sin(db)*math.sin(dc)) + fid = trace_to_fid(tr) + return fid >= fidelity + + if fidelity is None: # Don't specialize if None + instance = super().__new__(TwoQubitWeylGeneral + if cls is TwoQubitWeylDecomposition else cls) + elif is_close(0, 0, 0): + instance = super().__new__(TwoQubitWeylIdEquiv) + elif is_close(pi4, pi4, pi4) or is_close(pi4, pi4, -pi4): + instance = super().__new__(TwoQubitWeylSWAPEquiv) + elif (lambda x: is_close(x, x, x))(_closest_partial_swap(a, b, c)): + instance = super().__new__(TwoQubitWeylPartialSWAPEquiv) + elif (lambda x: is_close(x, x, -x))(_closest_partial_swap(a, b, -c)): + instance = super().__new__(TwoQubitWeylPartialSWAPFlipEquiv) + elif is_close(a, 0, 0): + instance = super().__new__(TwoQubitWeylControlledEquiv) + elif is_close(pi4, pi4, c): + instance = super().__new__(TwoQubitWeylMirrorControlledEquiv) + elif is_close((a+b)/2, (a+b)/2, c): + instance = super().__new__(TwoQubitWeylfSimaabEquiv) + elif is_close(a, (b+c)/2, (b+c)/2): + instance = super().__new__(TwoQubitWeylfSimabbEquiv) + elif is_close(a, (b-c)/2, (c-b)/2): + instance = super().__new__(TwoQubitWeylfSimabmbEquiv) + else: + instance = super().__new__(TwoQubitWeylGeneral) + + instance._original_decomposition = od + return instance + + def __init__(self, unitary_matrix, fidelity=None): + del unitary_matrix # unused in __init__ (used in new) + od = self._original_decomposition + self.a, self.b, self.c = od.a, od.b, od.c + self.K1l, self.K1r = od.K1l, od.K1r + self.K2l, self.K2r = od.K2l, od.K2r + self.global_phase = od.global_phase + self.unitary_matrix = od.unitary_matrix + self.requested_fidelity = fidelity + self._is_flipped_from_original = False + self.specialize() + + # Update the phase after specialization: + if self._is_flipped_from_original: + da, db, dc = (np.pi/2-od.a)-self.a, od.b-self.b, -od.c-self.c + tr = 4 * complex(math.cos(da)*math.cos(db)*math.cos(dc), + math.sin(da)*math.sin(db)*math.sin(dc)) + else: + da, db, dc = od.a-self.a, od.b-self.b, od.c-self.c + tr = 4 * complex(math.cos(da)*math.cos(db)*math.cos(dc), + math.sin(da)*math.sin(db)*math.sin(dc)) + self.global_phase += cmath.phase(tr) + self.calculated_fidelity = trace_to_fid(tr) + if logger.isEnabledFor(logging.DEBUG): + actual_fidelity = self.actual_fidelity() + logger.debug("Requested fidelity: %s calculated fidelity: %s actual fidelity %s", + self.requested_fidelity, self.calculated_fidelity, actual_fidelity) + if abs(self.calculated_fidelity - actual_fidelity) > 1.E-12: + logger.warning("Requested fidelity different from actual by %s", + self.calculated_fidelity - actual_fidelity) + if self.requested_fidelity and self.calculated_fidelity + 1.E-13 < self.requested_fidelity: + raise QiskitError(f"{self.__class__.__name__}: " + f"calculated fidelity: {self.calculated_fidelity} " + f"is worse than requested fidelity: {self.requested_fidelity}.") + + def specialize(self): + """Make changes to the decomposition to comply with any specialization. + + Do update a, b, c, k1l, k1r, k2l, k2r, _is_flipped_from_original to round to the + specialization. Do not update the global phase, since this gets done in generic + __init__()""" + raise NotImplementedError + + def circuit(self, *, euler_basis: Optional[str] = None, + simplify=False, atol=DEFAULT_ATOL) -> QuantumCircuit: + """Returns Weyl decomposition in circuit form. + + simplify, atol arguments are passed to OneQubitEulerDecomposer""" + if euler_basis is None: + euler_basis = self._default_1q_basis + oneq_decompose = OneQubitEulerDecomposer(euler_basis) + c1l, c1r, c2l, c2r = (oneq_decompose(k, simplify=simplify, atol=atol) + for k in (self.K1l, self.K1r, self.K2l, self.K2r)) + circ = QuantumCircuit(2, global_phase=self.global_phase) + circ.compose(c2r, [0], inplace=True) + circ.compose(c2l, [1], inplace=True) + self._weyl_gate(simplify, circ, atol) + circ.compose(c1r, [0], inplace=True) + circ.compose(c1l, [1], inplace=True) + return circ + + def _weyl_gate(self, simplify, circ: QuantumCircuit, atol): + """Appends Ud(a, b, c) to the circuit. + + Can be overriden in subclasses for special cases""" + if not simplify or abs(self.a) > atol: + circ.rxx(-self.a*2, 0, 1) + if not simplify or abs(self.b) > atol: + circ.ryy(-self.b*2, 0, 1) + if not simplify or abs(self.c) > atol: + circ.rzz(-self.c*2, 0, 1) + + def actual_fidelity(self, **kwargs) -> float: + """Calculates the actual fidelity of the decomposed circuit to the input unitary""" + circ = self.circuit(**kwargs) + trace = np.trace(Operator(circ).data.T.conj() @ self.unitary_matrix) + return trace_to_fid(trace) def __repr__(self): - # FIXME: this is worth making prettier since it's very useful for debugging - return ("{}\n{}\n{}\nUd({}, {}, {})\n{}\n{}\n".format( - pi_check(self.global_phase), - np.array_str(self.K1l), - np.array_str(self.K1r), - self.a, self.b, self.c, - np.array_str(self.K2l), - np.array_str(self.K2r))) + """Represent with enough precision to allow copy-paste debugging of all corner cases + """ + prefix = f"{type(self).__qualname__}.from_bytes(" + with io.BytesIO() as f: + np.save(f, self.unitary_matrix, allow_pickle=False) + b64 = base64.encodebytes(f.getvalue()).splitlines() + b64ascii = [repr(x) for x in b64] + b64ascii[-1] += "," + pretty = [f'# {x.rstrip()}' for x in str(self).splitlines()] + indent = '\n' + ' '*4 + lines = ([prefix] + pretty + b64ascii + + [f"requested_fidelity={self.requested_fidelity},", + f"calculated_fidelity={self.calculated_fidelity},", + f"actual_fidelity={self.actual_fidelity()},", + f"abc={(self.a, self.b, self.c)})"]) + return indent.join(lines) + + @classmethod + def from_bytes(cls, bytes_in: bytes, *, requested_fidelity: float, **kwargs + ) -> "TwoQubitWeylDecomposition": + """Decode bytes into TwoQubitWeylDecomposition. Used by __repr__""" + del kwargs # Unused (just for display) + b64 = base64.decodebytes(bytes_in) + with io.BytesIO(b64) as f: + arr = np.load(f, allow_pickle=False) + return cls(arr, fidelity=requested_fidelity) + + def __str__(self): + pre = f"{self.__class__.__name__}(\n\t" + circ_indent = "\n\t".join(self.circuit(simplify=True).draw("text").lines(-1)) + return f"{pre}{circ_indent}\n)" + + +class TwoQubitWeylIdEquiv(TwoQubitWeylDecomposition): + """U ~ Ud(0,0,0) ~ Id + + This gate binds 0 parameters, we make it canonical by setting + K2l = Id , + K2r = Id . + """ + + def specialize(self): + self.a = self.b = self.c = 0. + self.K1l = self.K1l @ self.K2l + self.K1r = self.K1r @ self.K2r + self.K2l = _id.copy() + self.K2r = _id.copy() + + +class TwoQubitWeylSWAPEquiv(TwoQubitWeylDecomposition): + """U ~ Ud(𝜋/4, 𝜋/4, 𝜋/4) ~ U(𝜋/4, 𝜋/4, -𝜋/4) ~ SWAP + + This gate binds 0 parameters, we make it canonical by setting + K2l = Id , + K2r = Id . + """ + def specialize(self): + if self.c > 0: + self.K1l = self.K1l @ self.K2r + self.K1r = self.K1r @ self.K2l + else: + self._is_flipped_from_original = True + self.K1l = self.K1l @ _ipz @ self.K2r + self.K1r = self.K1r @ _ipz @ self.K2l + self.global_phase = self.global_phase + np.pi/2 + self.a = self.b = self.c = np.pi/4 + self.K2l = _id.copy() + self.K2r = _id.copy() + + def _weyl_gate(self, simplify, circ: QuantumCircuit, atol): + del self, simplify, atol # unused + circ.swap(0, 1) + circ.global_phase -= 3*np.pi/4 + + +def _closest_partial_swap(a, b, c) -> float: + """A good approximation to the best value x to get the minimum + trace distance for Ud(x, x, x) from Ud(a, b, c) + """ + m = (a + b + c) / 3 + am, bm, cm = a-m, b-m, c-m + ab, bc, ca = a-b, b-c, c-a + + return m + am * bm * cm * (6 + ab*ab + bc*bc * ca*ca) / 18 + + +class TwoQubitWeylPartialSWAPEquiv(TwoQubitWeylDecomposition): + """U ~ Ud(α𝜋/4, α𝜋/4, α𝜋/4) ~ SWAP**α + + This gate binds 3 parameters, we make it canonical by setting: + K2l = Id . + """ + def specialize(self): + self.a = self.b = self.c = _closest_partial_swap(self.a, self.b, self.c) + self.K1l = self.K1l @ self.K2l + self.K1r = self.K1r @ self.K2l + self.K2r = self.K2l.T.conj() @ self.K2r + self.K2l = _id.copy() + + +class TwoQubitWeylPartialSWAPFlipEquiv(TwoQubitWeylDecomposition): + """U ~ Ud(α𝜋/4, α𝜋/4, -α𝜋/4) ~ SWAP**α + + (a non-equivalent root of SWAP from the TwoQubitWeylPartialSWAPEquiv + similar to how x = (±sqrt(x))**2 ) + + This gate binds 3 parameters, we make it canonical by setting: + K2l = Id . + """ + + def specialize(self): + self.a = self.b = _closest_partial_swap(self.a, self.b, -self.c) + self.c = -self.a + self.K1l = self.K1l @ self.K2l + self.K1r = self.K1r @ _ipz @ self.K2l @ _ipz + self.K2r = _ipz @ self.K2l.T.conj() @ _ipz @ self.K2r + self.K2l = _id.copy() + + +_oneq_xyx = OneQubitEulerDecomposer('XYX') +_oneq_zyz = OneQubitEulerDecomposer('ZYZ') + + +class TwoQubitWeylControlledEquiv(TwoQubitWeylDecomposition): + """U ~ Ud(α, 0, 0) ~ Ctrl-U + + This gate binds 4 parameters, we make it canonical by setting: + K2l = Ry(θl).Rx(λl) , + K2r = Ry(θr).Rx(λr) . + """ + _default_1q_basis = 'XYX' + + def specialize(self): + self.b = self.c = 0 + k2ltheta, k2lphi, k2llambda, k2lphase = _oneq_xyx.angles_and_phase(self.K2l) + k2rtheta, k2rphi, k2rlambda, k2rphase = _oneq_xyx.angles_and_phase(self.K2r) + self.global_phase += k2lphase + k2rphase + self.K1l = self.K1l @ np.asarray(RXGate(k2lphi)) + self.K1r = self.K1r @ np.asarray(RXGate(k2rphi)) + self.K2l = np.asarray(RYGate(k2ltheta)) @ np.asarray(RXGate(k2llambda)) + self.K2r = np.asarray(RYGate(k2rtheta)) @ np.asarray(RXGate(k2rlambda)) + + +class TwoQubitWeylMirrorControlledEquiv(TwoQubitWeylDecomposition): + """U ~ Ud(𝜋/4, 𝜋/4, α) ~ SWAP . Ctrl-U + + This gate binds 4 parameters, we make it canonical by setting: + K2l = Ry(θl).Rz(λl) , + K2r = Ry(θr).Rz(λr) . + """ + def specialize(self): + self.a = self.b = np.pi/4 + k2ltheta, k2lphi, k2llambda, k2lphase = _oneq_zyz.angles_and_phase(self.K2l) + k2rtheta, k2rphi, k2rlambda, k2rphase = _oneq_zyz.angles_and_phase(self.K2r) + self.global_phase += k2lphase + k2rphase + self.K1r = self.K1r @ np.asarray(RZGate(k2lphi)) + self.K1l = self.K1l @ np.asarray(RZGate(k2rphi)) + self.K2l = np.asarray(RYGate(k2ltheta)) @ np.asarray(RZGate(k2llambda)) + self.K2r = np.asarray(RYGate(k2rtheta)) @ np.asarray(RZGate(k2rlambda)) + + def _weyl_gate(self, simplify, circ: QuantumCircuit, atol): + circ.swap(0, 1) + circ.rzz((np.pi/4 - self.c) * 2, 0, 1) + circ.global_phase += np.pi/4 + + +# These next 3 gates use the definition of fSim from https://arxiv.org/pdf/2001.08343.pdf eq (1) +class TwoQubitWeylfSimaabEquiv(TwoQubitWeylDecomposition): + """U ~ Ud(α, α, β), α ≥ |β| + + This gate binds 5 parameters, we make it canonical by setting: + K2l = Ry(θl).Rz(λl) . + """ + def specialize(self): + self.a = self.b = (self.a + self.b)/2 + k2ltheta, k2lphi, k2llambda, k2lphase = _oneq_zyz.angles_and_phase(self.K2l) + self.global_phase += k2lphase + self.K1r = self.K1r @ np.asarray(RZGate(k2lphi)) + self.K1l = self.K1l @ np.asarray(RZGate(k2lphi)) + self.K2l = np.asarray(RYGate(k2ltheta)) @ np.asarray(RZGate(k2llambda)) + self.K2r = np.asarray(RZGate(-k2lphi)) @ self.K2r + + +class TwoQubitWeylfSimabbEquiv(TwoQubitWeylDecomposition): + """U ~ Ud(α, β, β), α ≥ β + + This gate binds 5 parameters, we make it canonical by setting: + K2l = Ry(θl).Rx(λl) . + """ + _default_1q_basis = 'XYX' + + def specialize(self): + self.b = self.c = (self.b + self.c)/2 + k2ltheta, k2lphi, k2llambda, k2lphase = _oneq_xyx.angles_and_phase(self.K2l) + self.global_phase += k2lphase + self.K1r = self.K1r @ np.asarray(RXGate(k2lphi)) + self.K1l = self.K1l @ np.asarray(RXGate(k2lphi)) + self.K2l = np.asarray(RYGate(k2ltheta)) @ np.asarray(RXGate(k2llambda)) + self.K2r = np.asarray(RXGate(-k2lphi)) @ self.K2r + + +class TwoQubitWeylfSimabmbEquiv(TwoQubitWeylDecomposition): + """U ~ Ud(α, β, -β), α ≥ β ≥ 0 + + This gate binds 5 parameters, we make it canonical by setting: + K2l = Ry(θl).Rx(λl) . + """ + + _default_1q_basis = 'XYX' + + def specialize(self): + self.b = (self.b - self.c)/2 + self.c = -self.b + k2ltheta, k2lphi, k2llambda, k2lphase = _oneq_xyx.angles_and_phase(self.K2l) + self.global_phase += k2lphase + self.K1r = self.K1r @ _ipz @ np.asarray(RXGate(k2lphi)) @ _ipz + self.K1l = self.K1l @ np.asarray(RXGate(k2lphi)) + self.K2l = np.asarray(RYGate(k2ltheta)) @ np.asarray(RXGate(k2llambda)) + self.K2r = _ipz @ np.asarray(RXGate(-k2lphi)) @ _ipz @ self.K2r + + +class TwoQubitWeylGeneral(TwoQubitWeylDecomposition): + """U has no special symmetry. + + This gate binds all 6 possible parameters, so there is no need to make the single-qubit + pre-/post-gates canonical. + """ + def specialize(self): + pass # Nothing to do def Ud(a, b, c): """Generates the array Exp(i(a xx + b yy + c zz)) """ - return np.array([[np.exp(1j*c)*np.cos(a-b), 0, 0, 1j*np.exp(1j*c)*np.sin(a-b)], - [0, np.exp(-1j*c)*np.cos(a+b), 1j*np.exp(-1j*c)*np.sin(a+b), 0], - [0, 1j*np.exp(-1j*c)*np.sin(a+b), np.exp(-1j*c)*np.cos(a+b), 0], - [1j*np.exp(1j*c)*np.sin(a-b), 0, 0, np.exp(1j*c)*np.cos(a-b)]], dtype=complex) + return np.array([[cmath.exp(1j*c)*math.cos(a-b), 0, 0, 1j*cmath.exp(1j*c)*math.sin(a-b)], + [0, cmath.exp(-1j*c)*math.cos(a+b), 1j*cmath.exp(-1j*c)*math.sin(a+b), 0], + [0, 1j*cmath.exp(-1j*c)*math.sin(a+b), cmath.exp(-1j*c)*math.cos(a+b), 0], + [1j*cmath.exp(1j*c)*math.sin(a-b), 0, 0, cmath.exp(1j*c)*math.cos(a-b)]], + dtype=complex) def trace_to_fid(trace): """Average gate fidelity is :math:`Fbar = (d + |Tr (Utarget \\cdot U^dag)|^2) / d(d+1)` M. Horodecki, P. Horodecki and R. Horodecki, PRA 60, 1888 (1999)""" - return (4 + np.abs(trace)**2)/20 + return (4 + abs(trace)**2)/20 def rz_array(theta): @@ -245,8 +628,8 @@ def rz_array(theta): Rz(theta) = diag(exp(-i*theta/2),exp(i*theta/2)) """ - return np.array([[np.exp(-1j*theta/2.0), 0], - [0, np.exp(1j*theta/2.0)]], dtype=complex) + return np.array([[cmath.exp(-1j*theta/2.0), 0], + [0, cmath.exp(1j*theta/2.0)]], dtype=complex) class TwoQubitBasisDecomposer(): @@ -272,32 +655,32 @@ def __init__(self, gate, basis_fidelity=1.0, euler_basis=None): self._decomposer1q = OneQubitEulerDecomposer('U3') # FIXME: find good tolerances - self.is_supercontrolled = np.isclose(basis.a, np.pi/4) and np.isclose(basis.c, 0.) + self.is_supercontrolled = math.isclose(basis.a, np.pi/4) and math.isclose(basis.c, 0.) # Create some useful matrices U1, U2, U3 are equivalent to the basis, # expand as Ui = Ki1.Ubasis.Ki2 b = basis.b - K11l = 1/(1+1j) * np.array([[-1j*np.exp(-1j*b), np.exp(-1j*b)], - [-1j*np.exp(1j*b), -np.exp(1j*b)]], dtype=complex) - K11r = 1/np.sqrt(2) * np.array([[1j*np.exp(-1j*b), -np.exp(-1j*b)], - [np.exp(1j*b), -1j*np.exp(1j*b)]], dtype=complex) + K11l = 1/(1+1j) * np.array([[-1j*cmath.exp(-1j*b), cmath.exp(-1j*b)], + [-1j*cmath.exp(1j*b), -cmath.exp(1j*b)]], dtype=complex) + K11r = 1/math.sqrt(2) * np.array([[1j*cmath.exp(-1j*b), -cmath.exp(-1j*b)], + [cmath.exp(1j*b), -1j*cmath.exp(1j*b)]], dtype=complex) K12l = 1/(1+1j) * np.array([[1j, 1j], [-1, 1]], dtype=complex) - K12r = 1/np.sqrt(2) * np.array([[1j, 1], - [-1, -1j]], dtype=complex) - K32lK21l = 1/np.sqrt(2) * np.array([[1+1j*np.cos(2*b), 1j*np.sin(2*b)], - [1j*np.sin(2*b), 1-1j*np.cos(2*b)]], dtype=complex) - K21r = 1/(1-1j) * np.array([[-1j*np.exp(-2j*b), np.exp(-2j*b)], - [1j*np.exp(2j*b), np.exp(2j*b)]], dtype=complex) - K22l = 1/np.sqrt(2) * np.array([[1, -1], - [1, 1]], dtype=complex) + K12r = 1/math.sqrt(2) * np.array([[1j, 1], + [-1, -1j]], dtype=complex) + K32lK21l = 1/math.sqrt(2) * np.array([[1+1j*np.cos(2*b), 1j*np.sin(2*b)], + [1j*np.sin(2*b), 1-1j*np.cos(2*b)]], dtype=complex) + K21r = 1/(1-1j) * np.array([[-1j*cmath.exp(-2j*b), cmath.exp(-2j*b)], + [1j*cmath.exp(2j*b), cmath.exp(2j*b)]], dtype=complex) + K22l = 1/math.sqrt(2) * np.array([[1, -1], + [1, 1]], dtype=complex) K22r = np.array([[0, 1], [-1, 0]], dtype=complex) - K31l = 1/np.sqrt(2) * np.array([[np.exp(-1j*b), np.exp(-1j*b)], - [-np.exp(1j*b), np.exp(1j*b)]], dtype=complex) - K31r = 1j * np.array([[np.exp(1j*b), 0], - [0, -np.exp(-1j*b)]], dtype=complex) - K32r = 1/(1-1j) * np.array([[np.exp(1j*b), -np.exp(-1j*b)], - [-1j*np.exp(1j*b), -1j*np.exp(-1j*b)]], dtype=complex) + K31l = 1/math.sqrt(2) * np.array([[cmath.exp(-1j*b), cmath.exp(-1j*b)], + [-cmath.exp(1j*b), cmath.exp(1j*b)]], dtype=complex) + K31r = 1j * np.array([[cmath.exp(1j*b), 0], + [0, -cmath.exp(-1j*b)]], dtype=complex) + K32r = 1/(1-1j) * np.array([[cmath.exp(1j*b), -cmath.exp(-1j*b)], + [-1j*cmath.exp(1j*b), -1j*cmath.exp(-1j*b)]], dtype=complex) k1ld = basis.K1l.T.conj() k1rd = basis.K1r.T.conj() k2ld = basis.K2l.T.conj() @@ -330,7 +713,8 @@ def __init__(self, gate, basis_fidelity=1.0, euler_basis=None): # In the future could use different decomposition functions for different basis classes, etc if not self.is_supercontrolled: warnings.warn("Only know how to decompose properly for supercontrolled basis gate. " - "This gate is ~Ud({}, {}, {})".format(basis.a, basis.b, basis.c)) + "This gate is ~Ud({}, {}, {})".format(basis.a, basis.b, basis.c), + stacklevel=2) self.decomposition_fns = [self.decomp0, self.decomp1, self.decomp2_supercontrolled, @@ -342,16 +726,17 @@ def traces(self, target): # Future gotcha: extending this to non-supercontrolled basis. # Careful: closest distance between a1,b1,c1 and a2,b2,c2 may be between reflections. # This doesn't come up if either c1==0 or c2==0 but otherwise be careful. - - return [4*(np.cos(target.a)*np.cos(target.b)*np.cos(target.c) + - 1j*np.sin(target.a)*np.sin(target.b)*np.sin(target.c)), - 4*(np.cos(np.pi/4-target.a)*np.cos(self.basis.b-target.b)*np.cos(target.c) + - 1j*np.sin(np.pi/4-target.a)*np.sin(self.basis.b-target.b)*np.sin(target.c)), - 4*np.cos(target.c), + ta, tb, tc = target.a, target.b, target.c + bb = self.basis.b + return [4*complex(math.cos(ta)*math.cos(tb)*math.cos(tc), + math.sin(ta)*math.sin(tb)*math.sin(tc)), + 4*complex(math.cos(math.pi/4-ta)*math.cos(bb-tb)*math.cos(tc), + math.sin(math.pi/4-ta)*math.sin(bb-tb)*math.sin(tc)), + 4*math.cos(tc), 4] @staticmethod - def decomp0(target, eps=1e-15): + def decomp0(target): """Decompose target ~Ud(x, y, z) with 0 uses of the basis gate. Result Ur has trace: :math:`|Tr(Ur.Utarget^dag)| = 4|(cos(x)cos(y)cos(z)+ j sin(x)sin(y)sin(z)|`, @@ -359,10 +744,6 @@ def decomp0(target, eps=1e-15): U0l = target.K1l.dot(target.K2l) U0r = target.K1r.dot(target.K2r) - U0l.real[abs(U0l.real) < eps] = 0.0 - U0l.imag[abs(U0l.imag) < eps] = 0.0 - U0r.real[abs(U0r.real) < eps] = 0.0 - U0r.imag[abs(U0r.imag) < eps] = 0.0 return U0r, U0l def decomp1(self, target): @@ -422,33 +803,22 @@ def decomp3_supercontrolled(self, target): return U3r, U3l, U2r, U2l, U1r, U1l, U0r, U0l - def __call__(self, target, basis_fidelity=None): + def __call__(self, target, basis_fidelity=None, *, _num_basis_uses=None) -> QuantumCircuit: """Decompose a two-qubit unitary over fixed basis + SU(2) using the best approximation given that each basis application has a finite fidelity. + + You can force a particular approximation by passing _num_basis_uses. """ basis_fidelity = basis_fidelity or self.basis_fidelity - if hasattr(target, 'to_operator'): - # If input is a BaseOperator subclass this attempts to convert - # the object to an Operator so that we can extract the underlying - # numpy matrix from `Operator.data`. - target = target.to_operator().data - if hasattr(target, 'to_matrix'): - # If input is Gate subclass or some other class object that has - # a to_matrix method this will call that method. - target = target.to_matrix() - # Convert to numpy array incase not already an array target = np.asarray(target, dtype=complex) - # Check input is a 2-qubit unitary - if target.shape != (4, 4): - raise QiskitError("TwoQubitBasisDecomposer: expected 4x4 matrix for target") - if not is_unitary_matrix(target): - raise QiskitError("TwoQubitBasisDecomposer: target matrix is not unitary.") target_decomposed = TwoQubitWeylDecomposition(target) traces = self.traces(target_decomposed) expected_fidelities = [trace_to_fid(traces[i]) * basis_fidelity**i for i in range(4)] - best_nbasis = np.argmax(expected_fidelities) + best_nbasis = int(np.argmax(expected_fidelities)) + if _num_basis_uses is not None: + best_nbasis = _num_basis_uses decomposition = self.decomposition_fns[best_nbasis](target_decomposed) decomposition_euler = [self._decomposer1q._decompose(x) for x in decomposition] @@ -471,16 +841,12 @@ def num_basis_gates(self, unitary): """ Computes the number of basis gates needed in a decomposition of input unitary """ - if hasattr(unitary, 'to_operator'): - unitary = unitary.to_operator().data - if hasattr(unitary, 'to_matrix'): - unitary = unitary.to_matrix() unitary = np.asarray(unitary, dtype=complex) a, b, c = weyl_coordinates(unitary)[:] - traces = [4*(np.cos(a)*np.cos(b)*np.cos(c)+1j*np.sin(a)*np.sin(b)*np.sin(c)), - 4*(np.cos(np.pi/4-a)*np.cos(self.basis.b-b)*np.cos(c) + - 1j*np.sin(np.pi/4-a)*np.sin(self.basis.b-b)*np.sin(c)), - 4*np.cos(c), + traces = [4*(math.cos(a)*math.cos(b)*math.cos(c)+1j*math.sin(a)*math.sin(b)*math.sin(c)), + 4*(math.cos(np.pi/4-a)*math.cos(self.basis.b-b)*math.cos(c) + + 1j*math.sin(np.pi/4-a)*math.sin(self.basis.b-b)*math.sin(c)), + 4*math.cos(c), 4] return np.argmax([trace_to_fid(traces[i]) * self.basis_fidelity**i for i in range(4)]) diff --git a/releasenotes/notes/ZSXX-decomposition-90fead72778e410e.yaml b/releasenotes/notes/ZSXX-decomposition-90fead72778e410e.yaml new file mode 100644 index 000000000000..d5c878aafb37 --- /dev/null +++ b/releasenotes/notes/ZSXX-decomposition-90fead72778e410e.yaml @@ -0,0 +1,10 @@ +--- +features: + - | + The decomposition methods for single-qubit gates in + :class:`~qiskit.quantum_info.OneQubitEulerDecomposer` has been expanded to now + also include the ``'ZSXX'`` basis, for making use of direct :math:`X` gate as well as :math:`\sqrt{X}`. +fixes: + - | + Tightened up tolerances, improved repeatability and simplification, fixed several global-phase-tracking bugs + for one- and two-qubit gate synthesis. diff --git a/test/python/circuit/test_controlled_gate.py b/test/python/circuit/test_controlled_gate.py index 3a0382379dd7..a0048397fa1c 100644 --- a/test/python/circuit/test_controlled_gate.py +++ b/test/python/circuit/test_controlled_gate.py @@ -1190,11 +1190,11 @@ def test_single_controlled_rotation_gates(self, gate, cgate): self.log.info('\n%s', str(uqc)) # these limits could be changed if gate.name == 'ry': - self.assertLessEqual(uqc.size(), 32) + self.assertLessEqual(uqc.size(), 32, f"\n{uqc}") elif gate.name == 'rz': - self.assertLessEqual(uqc.size(), 42) + self.assertLessEqual(uqc.size(), 43, f"\n{uqc}") else: - self.assertLessEqual(uqc.size(), 20) + self.assertLessEqual(uqc.size(), 20, f"\n{uqc}") def test_composite(self): """Test composite gate count.""" @@ -1209,7 +1209,7 @@ def test_composite(self): unroller = Unroller(['u', 'cx']) uqc = dag_to_circuit(unroller.run(dag)) self.log.info('%s gate count: %d', uqc.name, uqc.size()) - self.assertLessEqual(uqc.size(), 95) # this limit could be changed + self.assertLessEqual(uqc.size(), 96, f"\n{uqc}") # this limit could be changed @ddt diff --git a/test/python/circuit/test_unitary.py b/test/python/circuit/test_unitary.py index 9c0f0dd03b3d..f6e2597afd68 100644 --- a/test/python/circuit/test_unitary.py +++ b/test/python/circuit/test_unitary.py @@ -258,7 +258,10 @@ def test_qasm_2q_unitary(self): qr = QuantumRegister(2, 'q0') cr = ClassicalRegister(1, 'c0') qc = QuantumCircuit(qr, cr) - matrix = numpy.eye(4) + matrix = numpy.asarray([[0, 0, 0, 1], + [0, 0, 1, 0], + [0, 1, 0, 0], + [1, 0, 0, 0]]) unitary_gate = UnitaryGate(matrix, label="custom_gate") qc.x(qr[0]) @@ -271,8 +274,8 @@ def test_qasm_2q_unitary(self): "creg c0[1];\n" \ "x q0[0];\n" \ "gate custom_gate p0,p1 {\n" \ - "\tu3(0,0,0) p0;\n" \ - "\tu3(0,0,0) p1;\n" \ + "\tu3(pi,-pi/2,pi/2) p0;\n" \ + "\tu3(pi,pi/2,-pi/2) p1;\n" \ "}\n" \ "custom_gate q0[0],q0[1];\n" \ "custom_gate q0[1],q0[0];\n" diff --git a/test/python/compiler/test_transpiler.py b/test/python/compiler/test_transpiler.py index 466c75198d1d..f8e4ca5523fd 100644 --- a/test/python/compiler/test_transpiler.py +++ b/test/python/compiler/test_transpiler.py @@ -557,7 +557,8 @@ def test_optimize_to_nothing(self): basis_gates=['u3', 'u2', 'u1', 'cx']) expected = QuantumCircuit(QuantumRegister(2, 'q'), global_phase=-np.pi/2) - self.assertEqual(after, expected) + msg = f"after:\n{after}\nexpected:\n{expected}" + self.assertEqual(after, expected, msg=msg) def test_pass_manager_empty(self): """Test passing an empty PassManager() to the transpiler. @@ -608,12 +609,12 @@ def test_initialize_reset_should_be_removed(self): qc.initialize([1.0 / math.sqrt(2), -1.0 / math.sqrt(2)], [qr[0]]) expected = QuantumCircuit(qr) - expected.append(U3Gate(1.5708, 0, 0), [qr[0]]) + expected.append(U3Gate(np.pi/2, 0, 0), [qr[0]]) expected.reset(qr[0]) - expected.append(U3Gate(1.5708, 3.1416, 0), [qr[0]]) + expected.append(U3Gate(np.pi/2, -np.pi, 0), [qr[0]]) after = transpile(qc, basis_gates=['reset', 'u3'], optimization_level=1) - self.assertEqual(after, expected) + self.assertEqual(after, expected, msg=f"after:\n{after}\nexpected:\n{expected}") def test_initialize_FakeMelbourne(self): """Test that the zero-state resets are remove in a device not supporting them. @@ -1014,12 +1015,12 @@ def test_no_infinite_loop(self, optimization_level): # for the second and third RZ gates in the U3 decomposition. expected = QuantumCircuit(1, global_phase=-np.pi/2 - 0.5 * (0.2 + np.pi) - 0.5 * 3 * np.pi) expected.sx(0) - expected.p(np.pi + 0.2, 0) + expected.p(-np.pi + 0.2, 0) expected.sx(0) - expected.p(np.pi, 0) + expected.p(-np.pi, 0) - error_message = "\nOutput circuit:\n%s\nExpected circuit:\n%s" % ( - str(out), str(expected)) + error_message = (f"\nOutput circuit:\n{out!s}\n{Operator(out).data}\n" + f"Expected circuit:\n{expected!s}\n{Operator(expected).data}") self.assertEqual(out, expected, error_message) @data(0, 1, 2, 3) diff --git a/test/python/quantum_info/states/test_statevector.py b/test/python/quantum_info/states/test_statevector.py index 491dede2c774..9dfe2bd217ef 100644 --- a/test/python/quantum_info/states/test_statevector.py +++ b/test/python/quantum_info/states/test_statevector.py @@ -944,7 +944,7 @@ def test_global_phase(self): qc2 = transpile(qc, basis_gates=['p']) sv = Statevector.from_instruction(qc2) expected = np.array([0.96891242-0.24740396j, 0]) - self.assertEqual(float(qc2.global_phase), -1/4) + self.assertEqual(float(qc2.global_phase), 2*np.pi - 0.25) self.assertEqual(sv, Statevector(expected)) def test_reverse_qargs(self): diff --git a/test/python/quantum_info/test_synthesis.py b/test/python/quantum_info/test_synthesis.py index 3fd909845297..3d5c40c2c326 100644 --- a/test/python/quantum_info/test_synthesis.py +++ b/test/python/quantum_info/test_synthesis.py @@ -13,26 +13,39 @@ """Tests for quantum synthesis methods.""" import unittest +import contextlib +import logging from test import combine from ddt import ddt import numpy as np -import scipy.linalg as la -from qiskit import execute +from qiskit import execute, QiskitError from qiskit.circuit import QuantumCircuit, QuantumRegister from qiskit.extensions import UnitaryGate from qiskit.circuit.library import (HGate, IGate, SdgGate, SGate, U3Gate, UGate, XGate, YGate, ZGate, CXGate, CZGate, - iSwapGate, RXXGate) + iSwapGate, RXXGate, RXGate, RYGate, RZGate) from qiskit.providers.basicaer import UnitarySimulatorPy from qiskit.quantum_info.operators import Operator from qiskit.quantum_info.random import random_unitary from qiskit.quantum_info.synthesis.one_qubit_decompose import OneQubitEulerDecomposer from qiskit.quantum_info.synthesis.two_qubit_decompose import (TwoQubitWeylDecomposition, + TwoQubitWeylIdEquiv, + TwoQubitWeylSWAPEquiv, + TwoQubitWeylPartialSWAPEquiv, + TwoQubitWeylPartialSWAPFlipEquiv, + TwoQubitWeylfSimaabEquiv, + TwoQubitWeylfSimabbEquiv, + TwoQubitWeylfSimabmbEquiv, + TwoQubitWeylControlledEquiv, + TwoQubitWeylMirrorControlledEquiv, + TwoQubitWeylGeneral, two_qubit_cnot_decompose, TwoQubitBasisDecomposer, - Ud) + Ud, + decompose_two_qubit_product_gate) + from qiskit.quantum_info.synthesis.ion_decompose import cnot_rxx_decompose from qiskit.test import QiskitTestCase @@ -78,8 +91,7 @@ def make_hard_thetas_oneq(smallest=1e-18, factor=3.2, steps=22, phi=0.7, lam=0.9 class CheckDecompositions(QiskitTestCase): """Implements decomposition checkers.""" - def check_one_qubit_euler_angles(self, operator, basis='U3', tolerance=1e-12, - phase_equal=True): + def check_one_qubit_euler_angles(self, operator, basis='U3', tolerance=1e-14, simplify=False): """Check OneQubitEulerDecomposer works for the given unitary""" target_unitary = operator.data if basis is None: @@ -87,21 +99,49 @@ def check_one_qubit_euler_angles(self, operator, basis='U3', tolerance=1e-12, decomp_unitary = U3Gate(*angles).to_matrix() else: decomposer = OneQubitEulerDecomposer(basis) - decomp_unitary = Operator(decomposer(target_unitary)).data - # Add global phase to make special unitary - target_unitary *= la.det(target_unitary) ** (-0.5) - decomp_unitary *= la.det(decomp_unitary) ** (-0.5) + decomp_unitary = Operator(decomposer(target_unitary, simplify=simplify)).data maxdist = np.max(np.abs(target_unitary - decomp_unitary)) - if not phase_equal and maxdist > 0.1: - maxdist = np.max(np.abs(target_unitary + decomp_unitary)) self.assertTrue(np.abs(maxdist) < tolerance, "Operator {}: Worst distance {}".format(operator, maxdist)) - # FIXME: should be possible to set this tolerance tighter after improving the function - def check_two_qubit_weyl_decomposition(self, target_unitary, tolerance=1.e-7): + @contextlib.contextmanager + def assertDebugOnly(self): # FIXME: when at python 3.10+ replace with assertNoLogs + """Context manager, asserts log is emitted at level DEBUG but no higher""" + with self.assertLogs("qiskit.quantum_info.synthesis", "DEBUG") as ctx: + yield + for i in range(len(ctx.records)): + self.assertLessEqual(ctx.records[i].levelno, logging.DEBUG, + msg=f"Unexpected logging entry: {ctx.output[i]}") + self.assertIn("Requested fidelity:", ctx.records[i].getMessage()) + + def assertRoundTrip(self, weyl1: TwoQubitWeylDecomposition): + """Fail if eval(repr(weyl1)) not equal to weyl1""" + repr1 = repr(weyl1) + with self.assertDebugOnly(): + weyl2: TwoQubitWeylDecomposition = eval(repr1) # pylint: disable=eval-used + msg_base = f"weyl1:\n{repr1}\nweyl2:\n{repr(weyl2)}" + self.assertEqual(type(weyl1), type(weyl2), msg_base) + maxdiff = np.max(abs(weyl1.unitary_matrix-weyl2.unitary_matrix)) + self.assertEqual(maxdiff, 0, msg=f"Unitary matrix differs by {maxdiff}\n" + msg_base) + self.assertEqual(weyl1.a, weyl2.a, msg=msg_base) + self.assertEqual(weyl1.b, weyl2.b, msg=msg_base) + self.assertEqual(weyl1.c, weyl2.c, msg=msg_base) + maxdiff = np.max(np.abs(weyl1.K1l - weyl2.K1l)) + self.assertEqual(maxdiff, 0, msg=f"K1l matrix differs by {maxdiff}" + msg_base) + maxdiff = np.max(np.abs(weyl1.K1r - weyl2.K1r)) + self.assertEqual(maxdiff, 0, msg=f"K1r matrix differs by {maxdiff}" + msg_base) + maxdiff = np.max(np.abs(weyl1.K2l - weyl2.K2l)) + self.assertEqual(maxdiff, 0, msg=f"K2l matrix differs by {maxdiff}" + msg_base) + maxdiff = np.max(np.abs(weyl1.K2r - weyl2.K2r)) + self.assertEqual(maxdiff, 0, msg=f"K2r matrix differs by {maxdiff}" + msg_base) + self.assertEqual(weyl1.requested_fidelity, weyl2.requested_fidelity, msg_base) + + def check_two_qubit_weyl_decomposition(self, target_unitary, tolerance=1.e-12): """Check TwoQubitWeylDecomposition() works for a given operator""" # pylint: disable=invalid-name - decomp = TwoQubitWeylDecomposition(target_unitary) + with self.assertDebugOnly(): + decomp = TwoQubitWeylDecomposition(target_unitary, fidelity=None) + # self.assertRoundTrip(decomp) # Too slow op = np.exp(1j * decomp.global_phase) * Operator(np.eye(4)) for u, qs in ( (decomp.K2r, [0]), @@ -113,12 +153,48 @@ def check_two_qubit_weyl_decomposition(self, target_unitary, tolerance=1.e-7): op = op.compose(u, qs) decomp_unitary = op.data maxdist = np.max(np.abs(target_unitary - decomp_unitary)) - self.assertTrue(np.abs(maxdist) < tolerance, - "Unitary {}: Worst distance {}".format(target_unitary, maxdist)) - - def check_exact_decomposition(self, target_unitary, decomposer, tolerance=1.e-7): + self.assertLess(np.abs(maxdist), tolerance, + f"{decomp}\nactual fid: {decomp.actual_fidelity()}\n" + f"Unitary {target_unitary}:\nWorst distance {maxdist}") + + def check_two_qubit_weyl_specialization(self, target_unitary, fidelity, + expected_specialization, expected_gates): + """Check that the two qubit Weyl decomposition gets specialized as expected""" + + # Loop to check both for implicit and explicity specialization + for decomposer in (TwoQubitWeylDecomposition, expected_specialization): + with self.assertDebugOnly(): + decomp = decomposer(target_unitary, fidelity=fidelity) + self.assertRoundTrip(decomp) + self.assertEqual(np.max(np.abs(decomp.unitary_matrix - target_unitary)), 0, + "Incorrect saved unitary in the decomposition.") + self.assertIsInstance(decomp, expected_specialization, + "Incorrect Weyl specialization.") + circ = decomp.circuit(simplify=True) + self.assertDictEqual(dict(circ.count_ops()), expected_gates, + f"Gate counts of {decomposer.__name__}") + actual_fid = decomp.actual_fidelity() + self.assertAlmostEqual(decomp.calculated_fidelity, actual_fid, places=13) + self.assertGreaterEqual(actual_fid, fidelity, + f"fidelity of {decomposer.__name__}") + actual_unitary = Operator(circ).data + trace = np.trace(actual_unitary.T.conj() @ target_unitary) + self.assertAlmostEqual(trace.imag, 0, places=13, + msg=f"Real trace for {decomposer.__name__}") + with self.assertDebugOnly(): + decomp2 = expected_specialization(target_unitary, fidelity=None) # Shouldn't raise + self.assertRoundTrip(decomp2) + if expected_specialization is not TwoQubitWeylGeneral: + with self.assertRaises(QiskitError) as exc: + _ = expected_specialization(target_unitary, fidelity=1.) + self.assertIn("worse than requested", exc.exception.message) + + def check_exact_decomposition(self, target_unitary, decomposer, + tolerance=1.e-12, num_basis_uses=None): """Check exact decomposition for a particular target""" - decomp_circuit = decomposer(target_unitary) + decomp_circuit = decomposer(target_unitary, _num_basis_uses=num_basis_uses) + if num_basis_uses is not None: + self.assertEqual(num_basis_uses, decomp_circuit.count_ops().get('unitary', 0)) result = execute(decomp_circuit, UnitarySimulatorPy(), optimization_level=0).result() decomp_unitary = result.get_unitary() maxdist = np.max(np.abs(target_unitary - decomp_unitary)) @@ -148,62 +224,200 @@ def test_euler_angles_1q_random(self, seed): self.check_one_qubit_euler_angles(unitary) +ANGEXP_ZYZ = [ # Special cases for ZYZ type expansions + [(1.E-13, 0.1, -0.1, 0), (0, 0)], + [(1.E-13, 0.2, -0.1, 0), (1, 0)], + [(1.E-13, np.pi, np.pi, 0), (0, 0)], + [(1.E-13, np.pi, np.pi, np.pi), (0, 0)], + [(np.pi, np.pi, np.pi, 0), (0, 1)], + [(np.pi-1.E-13, np.pi, np.pi, np.pi), (0, 1)], + [(np.pi, 0.1, 0.2, 0), (1, 1)], + [(np.pi, 0.2, 0.2, 0), (0, 1)], + [(1.E-13, 0.1, 0.2, 0), (1, 0)], + [(0.1, 0.2, 1.E-13, 0), (1, 1)], + [(0.1, 0., 0., 0), (0, 1)], + [(0.1, 1.E-13, 0.2, 0), (1, 1)], + [(0.1, 0.2, 0.3, 0), (2, 1)] +] +ANGEXP_PSX = [ # Special cases for Z.X90.Z.X90.Z type expansions + [(0.0, 0.1, -0.1), (0, 0)], + [(0.0, 0.1, 0.2), (1, 0)], + [(-np.pi/2, 0.2, 0.0), (2, 1)], + [(np.pi/2, 0.0, 0.21), (2, 1)], + [(np.pi/2, 0.12, 0.2), (2, 1)], + [(np.pi/2, -np.pi/2, 0.21), (1, 1)], + [(np.pi, np.pi, 0), (0, 2)], + [(np.pi, np.pi+0.1, 0.1), (0, 2)], + [(np.pi, np.pi+0.2, -0.1), (1, 2)], + [(0.1, 0.2, 0.3), (3, 2)], +] + + +@ddt +class TestOneQubitEulerSpecial(CheckDecompositions): + """Test special cases for OneQubitEulerDecomposer. + + FIXME: Currently these are more like smoke tests that exercise each of the code paths + and shapes of decompositions that can be made, but they don't check all the corner cases + where a wrap by 2*pi might happen, etc + """ + + def check_oneq_special_cases(self, target, basis, expected_gates=None, tolerance=1.E-12,): + """Check OneQubitEulerDecomposer produces the expected gates""" + decomposer = OneQubitEulerDecomposer(basis) + circ = decomposer(target, simplify=True) + data = Operator(circ).data + maxdist = np.max(np.abs(target.data - data)) + trace = np.trace(data.T.conj() @ target) + self.assertLess(np.abs(maxdist), tolerance, + f"Worst case distance: {maxdist}, trace: {trace}\n" + f"Target:\n{target}\nActual:\n{data}\n{circ}") + if expected_gates is not None: + self.assertDictEqual(dict(circ.count_ops()), expected_gates, + f"Circuit:\n{circ}") + + @combine(angexp=ANGEXP_ZYZ) + def test_special_ZYZ(self, angexp): + """Special cases of ZYZ. {angexp[0]}""" + a, b, c, d = angexp[0] + exp = {('rz', 'ry')[g]: angexp[1][g] for g in (0, 1) if angexp[1][g]} + tgt = np.exp(1j*d)*RZGate(b).to_matrix() @ RYGate(a).to_matrix() @ RZGate(c).to_matrix() + self.check_oneq_special_cases(tgt, 'ZYZ', exp) + + @combine(angexp=ANGEXP_ZYZ) + def test_special_ZXZ(self, angexp): + """Special cases of ZXZ. {angexp[0]}""" + a, b, c, d = angexp[0] + exp = {('rz', 'rx')[g]: angexp[1][g] for g in (0, 1) if angexp[1][g]} + tgt = np.exp(1j*d)*RZGate(b).to_matrix() @ RXGate(a).to_matrix() @ RZGate(c).to_matrix() + self.check_oneq_special_cases(tgt, 'ZXZ', exp) + + @combine(angexp=ANGEXP_ZYZ) + def test_special_XYX(self, angexp): + """Special cases of XYX. {angexp[0]}""" + a, b, c, d = angexp[0] + exp = {('rx', 'ry')[g]: angexp[1][g] for g in (0, 1) if angexp[1][g]} + tgt = np.exp(1j*d)*RXGate(b).to_matrix() @ RYGate(a).to_matrix() @ RXGate(c).to_matrix() + self.check_oneq_special_cases(tgt, 'XYX', exp) + + def test_special_U321(self): + """Special cases of U321""" + self.check_oneq_special_cases(U3Gate(0.0, 0.1, -0.1).to_matrix(), 'U321', {}) + self.check_oneq_special_cases(U3Gate(0.0, 0.11, 0.2).to_matrix(), 'U321', {'u1': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.2, 0.0).to_matrix(), 'U321', {'u2': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.0, 0.2).to_matrix(), 'U321', {'u2': 1}) + self.check_oneq_special_cases(U3Gate(0.11, 0.27, 0.3).to_matrix(), 'U321', {'u3': 1}) + + def test_special_U3(self): + """Special cases of U3""" + self.check_oneq_special_cases(U3Gate(0.0, 0.1, -0.1).to_matrix(), 'U3', {}) + self.check_oneq_special_cases(U3Gate(0.0, 0.1, 0.2).to_matrix(), 'U3', {'u3': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.2, 0.0).to_matrix(), 'U3', {'u3': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.0, 0.2).to_matrix(), 'U3', {'u3': 1}) + self.check_oneq_special_cases(U3Gate(0.11, 0.27, 0.3).to_matrix(), 'U3', {'u3': 1}) + + def test_special_U(self): + """Special cases of U""" + self.check_oneq_special_cases(U3Gate(0.0, 0.1, -0.1).to_matrix(), 'U', {}) + self.check_oneq_special_cases(U3Gate(0.0, 0.1, 0.2).to_matrix(), 'U', {'u': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.2, 0.0).to_matrix(), 'U', {'u': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.0, 0.2).to_matrix(), 'U', {'u': 1}) + self.check_oneq_special_cases(U3Gate(0.1, 0.2, 0.3).to_matrix(), 'U', {'u': 1}) + + def test_special_RR(self): + """Special cases of RR""" + self.check_oneq_special_cases(U3Gate(0.0, 0.1, -0.1).to_matrix(), 'RR', {}) + self.check_oneq_special_cases(U3Gate(0.0, 0.1, 0.2).to_matrix(), 'RR', {'r': 2}) + self.check_oneq_special_cases(U3Gate(-np.pi, 0.2, 0.0).to_matrix(), 'RR', {'r': 1}) + self.check_oneq_special_cases(U3Gate(np.pi, 0.0, 0.2).to_matrix(), 'RR', {'r': 1}) + self.check_oneq_special_cases(U3Gate(0.1, 0.2, 0.3).to_matrix(), 'RR', {'r': 2}) + + def test_special_U1X(self): + """Special cases of U1X""" + self.check_oneq_special_cases(U3Gate(0.0, 0.1, -0.1).to_matrix(), 'U1X', {}) + self.check_oneq_special_cases(U3Gate(0.0, 0.1, 0.2).to_matrix(), 'U1X', {'u1': 1}) + self.check_oneq_special_cases(U3Gate(-np.pi/2, 0.2, 0.0).to_matrix(), 'U1X', + {'u1': 2, 'rx': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.0, 0.21).to_matrix(), 'U1X', + {'u1': 2, 'rx': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.12, 0.2).to_matrix(), 'U1X', + {'u1': 2, 'rx': 1}) + self.check_oneq_special_cases(U3Gate(0.1, 0.2, 0.3).to_matrix(), 'U1X', {'u1': 3, 'rx': 2}) + + @combine(angexp=ANGEXP_PSX) + def test_special_PSX(self, angexp): + """Special cases of PSX. {angexp[0]}""" + a, b, c = angexp[0] + tgt = U3Gate(a, b, c).to_matrix() + exp = {('p', 'sx')[g]: angexp[1][g] for g in (0, 1) if angexp[1][g]} + self.check_oneq_special_cases(tgt, 'PSX', exp) + + def test_special_ZSX(self): + """Special cases of ZSX""" + self.check_oneq_special_cases(U3Gate(0.0, 0.1, -0.1).to_matrix(), 'ZSX', {}) + self.check_oneq_special_cases(U3Gate(0.0, 0.1, 0.2).to_matrix(), 'ZSX', {'rz': 1}) + self.check_oneq_special_cases(U3Gate(-np.pi/2, 0.2, 0.0).to_matrix(), 'ZSX', + {'rz': 2, 'sx': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.0, 0.21).to_matrix(), 'ZSX', + {'rz': 2, 'sx': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.12, 0.2).to_matrix(), 'ZSX', + {'rz': 2, 'sx': 1}) + self.check_oneq_special_cases(U3Gate(0.1, 0.2, 0.3).to_matrix(), 'ZSX', {'rz': 3, 'sx': 2}) + + def test_special_ZSXX(self): + """Special cases of ZSXX""" + self.check_oneq_special_cases(U3Gate(0.0, 0.1, -0.1).to_matrix(), 'ZSXX', {}) + self.check_oneq_special_cases(U3Gate(0.0, 0.1, 0.2).to_matrix(), 'ZSXX', {'rz': 1}) + self.check_oneq_special_cases(U3Gate(-np.pi/2, 0.2, 0.0).to_matrix(), 'ZSXX', + {'rz': 2, 'sx': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.0, 0.21).to_matrix(), 'ZSXX', + {'rz': 2, 'sx': 1}) + self.check_oneq_special_cases(U3Gate(np.pi/2, 0.12, 0.2).to_matrix(), 'ZSXX', + {'rz': 2, 'sx': 1}) + self.check_oneq_special_cases(U3Gate(0.1, 0.2, 0.3).to_matrix(), 'ZSXX', {'rz': 3, 'sx': 2}) + self.check_oneq_special_cases(U3Gate(np.pi, 0.2, 0.3).to_matrix(), 'ZSXX', + {'rz': 1, 'x': 1}) + self.check_oneq_special_cases(U3Gate(np.pi, -np.pi/2, np.pi/2).to_matrix(), 'ZSXX', + {'x': 1}) + + +ONEQ_BASES = ['U3', "U321", 'U', 'U1X', 'PSX', 'ZSX', 'ZSXX', 'ZYZ', 'ZXZ', 'XYX', 'RR'] +SIMP_TOL = [(False, 1.e-14), (True, 1.E-12)] # Please don't broaden the tolerance (fix the decomp) + + @ddt class TestOneQubitEulerDecomposer(CheckDecompositions): """Test OneQubitEulerDecomposer""" - def check_one_qubit_euler_angles(self, operator, basis='U3', - tolerance=1e-12, - phase_equal=True): - """Check euler_angles_1q works for the given unitary""" - decomposer = OneQubitEulerDecomposer(basis) - with self.subTest(operator=operator): - target_unitary = operator.data - decomp_unitary = Operator(decomposer(target_unitary)).data - # Add global phase to make special unitary - target_unitary *= la.det(target_unitary) ** (-0.5) - decomp_unitary *= la.det(decomp_unitary) ** (-0.5) - maxdist = np.max(np.abs(target_unitary - decomp_unitary)) - if not phase_equal and maxdist > 0.1: - maxdist = np.max(np.abs(target_unitary + decomp_unitary)) - self.assertTrue(np.abs(maxdist) < tolerance, "Worst distance {}".format(maxdist)) - - @combine(basis=['U3', 'U1X', 'PSX', 'ZSX', 'ZYZ', 'ZXZ', 'XYX', 'RR'], - name='test_one_qubit_clifford_{basis}_basis') - def test_one_qubit_clifford_all_basis(self, basis): + @combine(basis=ONEQ_BASES, simp_tol=SIMP_TOL, + name='test_one_qubit_clifford_{basis}_basis_simplify_{simp_tol[0]}') + def test_one_qubit_clifford_all_basis(self, basis, simp_tol): """Verify for {basis} basis and all Cliffords.""" for clifford in ONEQ_CLIFFORDS: - self.check_one_qubit_euler_angles(clifford, basis) - - @combine(basis_tolerance=[('U3', 1e-12), - ('XYX', 1e-12), - ('ZXZ', 1e-12), - ('ZYZ', 1e-12), - ('U1X', 1e-7), - ('PSX', 1e-7), - ('ZSX', 1e-7), - ('RR', 1e-12)], - name='test_one_qubit_hard_thetas_{basis_tolerance[0]}_basis') - # Lower tolerance for U1X test since decomposition since it is - # less numerically accurate. This is due to it having 5 matrix - # multiplications and the X90 gates - def test_one_qubit_hard_thetas_all_basis(self, basis_tolerance): - """Verify for {basis_tolerance[0]} basis and close-to-degenerate theta.""" + self.check_one_qubit_euler_angles( + clifford, basis, simplify=simp_tol[0], tolerance=simp_tol[1]) + + @combine(basis=ONEQ_BASES, simp_tol=SIMP_TOL, + name='test_one_qubit_hard_thetas_{basis}_basis_simplify_{simp_tol[0]}') + def test_one_qubit_hard_thetas_all_basis(self, basis, simp_tol): + """Verify for {basis} basis and close-to-degenerate theta.""" for gate in HARD_THETA_ONEQS: - self.check_one_qubit_euler_angles(Operator(gate), basis_tolerance[0], - basis_tolerance[1]) + self.check_one_qubit_euler_angles( + Operator(gate), basis, simplify=simp_tol[0], tolerance=simp_tol[1]) - @combine(basis=['U3', 'U1X', 'PSX', 'ZSX', 'ZYZ', 'ZXZ', 'XYX', 'RR'], seed=range(50), - name='test_one_qubit_random_{basis}_basis_{seed}') - def test_one_qubit_random_all_basis(self, basis, seed): + @combine(basis=ONEQ_BASES, simp_tol=SIMP_TOL, seed=range(50), + name='test_one_qubit_random_{basis}_basis_simplify_{simp_tol[0]}_{seed}') + def test_one_qubit_random_all_basis(self, basis, simp_tol, seed): """Verify for {basis} basis and random_unitary (seed={seed}).""" unitary = random_unitary(2, seed=seed) - self.check_one_qubit_euler_angles(unitary, basis) + self.check_one_qubit_euler_angles( + unitary, basis, simplify=simp_tol[0], tolerance=simp_tol[1]) def test_psx_zsx_special_cases(self): """Test decompositions of psx and zsx at special values of parameters""" oqed_psx = OneQubitEulerDecomposer(basis='PSX') oqed_zsx = OneQubitEulerDecomposer(basis='ZSX') + oqed_zsxx = OneQubitEulerDecomposer(basis='ZSXX') theta = np.pi / 3 phi = np.pi / 5 lam = np.pi / 7 @@ -227,8 +441,10 @@ def test_psx_zsx_special_cases(self): unitary = gate.to_matrix() qc_psx = oqed_psx(unitary) qc_zsx = oqed_zsx(unitary) + qc_zsxx = oqed_zsxx(unitary) self.assertTrue(np.allclose(unitary, Operator(qc_psx).data)) self.assertTrue(np.allclose(unitary, Operator(qc_zsx).data)) + self.assertTrue(np.allclose(unitary, Operator(qc_zsxx).data)) # FIXME: streamline the set of test cases @@ -236,6 +452,12 @@ class TestTwoQubitWeylDecomposition(CheckDecompositions): """Test TwoQubitWeylDecomposition() """ + def test_TwoQubitWeylDecomposition_repr(self, seed=42): + """Check that eval(__repr__) is exact round trip""" + target = random_unitary(4, seed=seed) + weyl1 = TwoQubitWeylDecomposition(target, fidelity=0.99) + self.assertRoundTrip(weyl1) + def test_two_qubit_weyl_decomposition_cnot(self): """Verify Weyl KAK decomposition for U~CNOT""" for k1l, k1r, k2l, k2r in K1K2S: @@ -361,7 +583,7 @@ def test_two_qubit_weyl_decomposition_aac(self, smallest=1e-18, factor=9.8, step self.check_two_qubit_weyl_decomposition(k1 @ a @ k2) def test_two_qubit_weyl_decomposition_abc(self, smallest=1e-18, factor=9.8, steps=11): - """Verify Weyl KAK decomposition for U~Ud(a,a,b)""" + """Verify Weyl KAK decomposition for U~Ud(a,b,c)""" for aaa in ([smallest * factor ** i for i in range(steps)] + [np.pi / 4 - smallest * factor ** i for i in range(steps)] + [np.pi / 8, 0.113 * np.pi, 0.1972 * np.pi]): @@ -370,15 +592,152 @@ def test_two_qubit_weyl_decomposition_abc(self, smallest=1e-18, factor=9.8, step for k1l, k1r, k2l, k2r in K1K2S: k1 = np.kron(k1l.data, k1r.data) k2 = np.kron(k2l.data, k2r.data) - a = Ud(aaa, aaa, ccc) + a = Ud(aaa, bbb, ccc) self.check_two_qubit_weyl_decomposition(k1 @ a @ k2) +K1K2SB = [[Operator(U3Gate(*xyz)) for xyz in xyzs] for xyzs in + [[(0.2, 0.3, 0.1), (0.7, 0.15, 0.22), (0.1, 0.97, 2.2), (3.14, 2.1, 0.9)], + [(0.21, 0.13, 0.45), (2.1, 0.77, 0.88), (1.5, 2.3, 2.3), (2.1, 0.4, 1.7)]]] +DELTAS = [(-0.019, 0.018, 0.021), (0.01, 0.015, 0.02), (-0.01, -0.009, 0.011), + (-0.002, -0.003, -0.004)] + + +class TestTwoQubitWeylDecompositionSpecialization(CheckDecompositions): + """Check TwoQubitWeylDecomposition specialized subclasses""" + + def test_weyl_specialize_id(self): + """Weyl specialization for Id gate""" + a, b, c = 0., 0., 0. + for da, db, dc in DELTAS: + for k1l, k1r, k2l, k2r in K1K2SB: + k1 = np.kron(k1l.data, k1r.data) + k2 = np.kron(k2l.data, k2r.data) + self.check_two_qubit_weyl_specialization(k1 @ Ud(a+da, b+db, c+dc) @ k2, + 0.999, TwoQubitWeylIdEquiv, + {'rz': 4, 'ry': 2}) + + def test_weyl_specialize_swap(self): + """Weyl specialization for swap gate""" + a, b, c = np.pi/4, np.pi/4, np.pi/4 + for da, db, dc in DELTAS: + for k1l, k1r, k2l, k2r in K1K2SB: + k1 = np.kron(k1l.data, k1r.data) + k2 = np.kron(k2l.data, k2r.data) + self.check_two_qubit_weyl_specialization(k1 @ Ud(a+da, b+db, c+dc) @ k2, + 0.999, TwoQubitWeylSWAPEquiv, + {'rz': 4, 'ry': 2, 'swap': 1}) + + def test_weyl_specialize_flip_swap(self): + """Weyl specialization for flip swap gate""" + a, b, c = np.pi/4, np.pi/4, -np.pi/4 + for da, db, dc in DELTAS: + for k1l, k1r, k2l, k2r in K1K2SB: + k1 = np.kron(k1l.data, k1r.data) + k2 = np.kron(k2l.data, k2r.data) + self.check_two_qubit_weyl_specialization(k1 @ Ud(a+da, b+db, c+dc) @ k2, + 0.999, TwoQubitWeylSWAPEquiv, + {'rz': 4, 'ry': 2, 'swap': 1}) + + def test_weyl_specialize_pswap(self, theta=0.123): + """Weyl specialization for partial swap gate""" + a, b, c = theta, theta, theta + for da, db, dc in DELTAS: + for k1l, k1r, k2l, k2r in K1K2SB: + k1 = np.kron(k1l.data, k1r.data) + k2 = np.kron(k2l.data, k2r.data) + self.check_two_qubit_weyl_specialization(k1 @ Ud(a+da, b+db, c+dc) @ k2, + 0.999, TwoQubitWeylPartialSWAPEquiv, + {'rz': 6, 'ry': 3, + 'rxx': 1, 'ryy': 1, 'rzz': 1}) + + def test_weyl_specialize_flip_pswap(self, theta=0.123): + """Weyl specialization for flipped partial swap gate""" + a, b, c = theta, theta, -theta + for da, db, dc in DELTAS: + for k1l, k1r, k2l, k2r in K1K2SB: + k1 = np.kron(k1l.data, k1r.data) + k2 = np.kron(k2l.data, k2r.data) + self.check_two_qubit_weyl_specialization(k1 @ Ud(a+da, b+db, c+dc) @ k2, + 0.999, TwoQubitWeylPartialSWAPFlipEquiv, + {'rz': 6, 'ry': 3, + 'rxx': 1, 'ryy': 1, 'rzz': 1}) + + def test_weyl_specialize_fsim_aab(self, aaa=0.456, bbb=0.132): + """Weyl specialization for partial swap gate""" + a, b, c = aaa, aaa, bbb + for da, db, dc in DELTAS: + for k1l, k1r, k2l, k2r in K1K2SB: + k1 = np.kron(k1l.data, k1r.data) + k2 = np.kron(k2l.data, k2r.data) + self.check_two_qubit_weyl_specialization(k1 @ Ud(a+da, b+db, c+dc) @ k2, + 0.999, TwoQubitWeylfSimaabEquiv, + {'rz': 7, 'ry': 4, + 'rxx': 1, 'ryy': 1, 'rzz': 1}) + + def test_weyl_specialize_fsim_abb(self, aaa=0.456, bbb=0.132): + """Weyl specialization for partial swap gate""" + a, b, c = aaa, bbb, bbb + for da, db, dc in DELTAS: + for k1l, k1r, k2l, k2r in K1K2SB: + k1 = np.kron(k1l.data, k1r.data) + k2 = np.kron(k2l.data, k2r.data) + self.check_two_qubit_weyl_specialization(k1 @ Ud(a+da, b+db, c+dc) @ k2, + 0.999, TwoQubitWeylfSimabbEquiv, + {'rx': 7, 'ry': 4, + 'rxx': 1, 'ryy': 1, 'rzz': 1}) + + def test_weyl_specialize_fsim_abmb(self, aaa=0.456, bbb=0.132): + """Weyl specialization for partial swap gate""" + a, b, c = aaa, bbb, -bbb + for da, db, dc in DELTAS: + for k1l, k1r, k2l, k2r in K1K2SB: + k1 = np.kron(k1l.data, k1r.data) + k2 = np.kron(k2l.data, k2r.data) + self.check_two_qubit_weyl_specialization(k1 @ Ud(a+da, b+db, c+dc) @ k2, + 0.999, TwoQubitWeylfSimabmbEquiv, + {'rx': 7, 'ry': 4, + 'rxx': 1, 'ryy': 1, 'rzz': 1}) + + def test_weyl_specialize_ctrl(self, aaa=0.456): + """Weyl specialization for partial swap gate""" + a, b, c = aaa, 0., 0. + for da, db, dc in DELTAS: + for k1l, k1r, k2l, k2r in K1K2SB: + k1 = np.kron(k1l.data, k1r.data) + k2 = np.kron(k2l.data, k2r.data) + self.check_two_qubit_weyl_specialization(k1 @ Ud(a+da, b+db, c+dc) @ k2, + 0.999, TwoQubitWeylControlledEquiv, + {'rx': 6, 'ry': 4, 'rxx': 1}) + + def test_weyl_specialize_mirror_ctrl(self, aaa=-0.456): + """Weyl specialization for partial swap gate""" + a, b, c = np.pi/4, np.pi/4, aaa + for da, db, dc in DELTAS: + for k1l, k1r, k2l, k2r in K1K2SB: + k1 = np.kron(k1l.data, k1r.data) + k2 = np.kron(k2l.data, k2r.data) + self.check_two_qubit_weyl_specialization(k1 @ Ud(a+da, b+db, c+dc) @ k2, + 0.999, TwoQubitWeylMirrorControlledEquiv, + {'rz': 6, 'ry': 4, 'rzz': 1, 'swap': 1}) + + def test_weyl_specialize_general(self, aaa=0.456, bbb=0.345, ccc=0.123): + """Weyl specialization for partial swap gate""" + a, b, c = aaa, bbb, ccc + for da, db, dc in DELTAS: + for k1l, k1r, k2l, k2r in K1K2SB: + k1 = np.kron(k1l.data, k1r.data) + k2 = np.kron(k2l.data, k2r.data) + self.check_two_qubit_weyl_specialization(k1 @ Ud(a+da, b+db, c+dc) @ k2, + 0.999, TwoQubitWeylGeneral, + {'rz': 8, 'ry': 4, + 'rxx': 1, 'ryy': 1, 'rzz': 1}) + + @ddt -class TestTwoQubitDecomposeExact(CheckDecompositions): - """Test TwoQubitBasisDecomposer() for exact decompositions +class TestTwoQubitDecompose(CheckDecompositions): + """Test TwoQubitBasisDecomposer() for exact/approx decompositions """ - def test_cnot_rxx_decompose(self): """Verify CNOT decomposition into RXX gate is correct""" cnot = Operator(CXGate()) @@ -403,20 +762,109 @@ def test_exact_two_qubit_cnot_decompose_paulis(self): unitary = Operator.from_label('XZ') self.check_exact_decomposition(unitary.data, two_qubit_cnot_decompose) + def make_random_supercontrolled_decomposer(self, seed): + """Return a random supercontrolled unitary given a seed""" + state = np.random.default_rng(seed) + basis_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_phase = state.random() * 2 * np.pi + basis_b = state.random() * np.pi / 4 + basis_unitary = np.exp(1j * basis_phase) * basis_k1 @ Ud(np.pi / 4, basis_b, 0) @ basis_k2 + decomposer = TwoQubitBasisDecomposer(UnitaryGate(basis_unitary)) + return decomposer + @combine(seed=range(10), name='test_exact_supercontrolled_decompose_random_{seed}') def test_exact_supercontrolled_decompose_random(self, seed): """Exact decomposition for random supercontrolled basis and random target (seed={seed})""" - k1 = np.kron(random_unitary(2, seed=seed).data, random_unitary(2, seed=seed + 1).data) - k2 = np.kron(random_unitary(2, seed=seed + 2).data, random_unitary(2, seed=seed + 3).data) - basis_unitary = k1 @ Ud(np.pi / 4, 0, 0) @ k2 + # pylint: disable=invalid-name + state = np.random.default_rng(seed) + decomposer = self.make_random_supercontrolled_decomposer(state) + self.check_exact_decomposition(random_unitary(4, seed=state).data, decomposer) + + @combine(seed=range(10), name='seed_{seed}') + def test_exact_supercontrolled_decompose_phase_0_use_random(self, seed): + """Exact decomposition supercontrolled basis, random target (0 basis uses) (seed={seed})""" + state = np.random.default_rng(seed) + decomposer = self.make_random_supercontrolled_decomposer(state) + + tgt_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_phase = state.random() * 2 * np.pi + tgt_unitary = np.exp(1j * tgt_phase) * tgt_k1 @ Ud(0, 0, 0) @ tgt_k2 + self.check_exact_decomposition(tgt_unitary, decomposer, num_basis_uses=0) + + @combine(seed=range(10), name='seed_{seed}') + def test_exact_supercontrolled_decompose_phase_1_use_random(self, seed): + """Exact decomposition supercontrolled basis, random tgt (1 basis uses) (seed={seed})""" + state = np.random.default_rng(seed) + basis_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_phase = state.random() * 2 * np.pi + basis_b = state.random() * np.pi / 4 + basis_unitary = np.exp(1j * basis_phase) * basis_k1 @ Ud(np.pi / 4, basis_b, 0) @ basis_k2 decomposer = TwoQubitBasisDecomposer(UnitaryGate(basis_unitary)) - self.check_exact_decomposition(random_unitary(4, seed=seed + 4).data, decomposer) + + tgt_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_phase = state.random() * 2 * np.pi + + tgt_unitary = np.exp(1j * tgt_phase) * tgt_k1 @ Ud(np.pi / 4, basis_b, 0) @ tgt_k2 + self.check_exact_decomposition(tgt_unitary, decomposer, num_basis_uses=1) + + @combine(seed=range(10), name='seed_{seed}') + def test_exact_supercontrolled_decompose_phase_2_use_random(self, seed): + """Exact decomposition supercontrolled basis, random tgt (2 basis uses) (seed={seed})""" + state = np.random.default_rng(seed) + decomposer = self.make_random_supercontrolled_decomposer(state) + + tgt_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_phase = state.random() * 2 * np.pi + tgt_a, tgt_b = state.random(size=2) * np.pi / 4 + tgt_unitary = np.exp(1j * tgt_phase) * tgt_k1 @ Ud(tgt_a, tgt_b, 0) @ tgt_k2 + self.check_exact_decomposition(tgt_unitary, decomposer, num_basis_uses=2) + + @combine(seed=range(10), name='seed_{seed}') + def test_exact_supercontrolled_decompose_phase_3_use_random(self, seed): + """Exact decomposition supercontrolled basis, random tgt (3 basis uses) (seed={seed})""" + state = np.random.default_rng(seed) + decomposer = self.make_random_supercontrolled_decomposer(state) + + tgt_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_phase = state.random() * 2 * np.pi + + tgt_a, tgt_b = state.random(size=2) * np.pi / 4 + tgt_c = state.random() * np.pi/2 - np.pi/4 + tgt_unitary = np.exp(1j * tgt_phase) * tgt_k1 @ Ud(tgt_a, tgt_b, tgt_c) @ tgt_k2 + self.check_exact_decomposition(tgt_unitary, decomposer, num_basis_uses=3) def test_exact_nonsupercontrolled_decompose(self): """Check that the nonsupercontrolled basis throws a warning""" with self.assertWarns(UserWarning, msg="Supposed to warn when basis non-supercontrolled"): TwoQubitBasisDecomposer(UnitaryGate(Ud(np.pi / 4, 0.2, 0.1))) + @combine(seed=range(10), name='seed_{seed}') + def test_approx_supercontrolled_decompose_random(self, seed): + """Check that n-uses of supercontrolled basis give the expected trace distance""" + state = np.random.default_rng(seed) + decomposer = self.make_random_supercontrolled_decomposer(state) + + tgt_phase = state.random() * 2 * np.pi + tgt = random_unitary(4, seed=state).data + tgt *= np.exp(1j * tgt_phase) + + with self.assertDebugOnly(): + traces_pred = decomposer.traces(TwoQubitWeylDecomposition(tgt)) + + for i in range(4): + decomp_circuit = decomposer(tgt, _num_basis_uses=i) + result = execute(decomp_circuit, UnitarySimulatorPy(), optimization_level=0).result() + decomp_unitary = result.get_unitary() + tr_actual = np.trace(decomp_unitary.conj().T @ tgt) + self.assertAlmostEqual(traces_pred[i], tr_actual, places=13, + msg=f"Trace doesn't match for {i}-basis decomposition") + def test_cx_equivalence_0cx(self, seed=0): """Check circuits with 0 cx gates locally equivalent to identity """ @@ -525,7 +973,8 @@ def test_seed_289(self): euler_bases=[('U321', ['u3', 'u2', 'u1']), ('U3', ['u3']), ('U', ['u']), ('U1X', ['u1', 'rx']), ('RR', ['r']), ('PSX', ['p', 'sx']), ('ZYZ', ['rz', 'ry']), ('ZXZ', ['rz', 'rx']), - ('XYX', ['rx', 'ry']), ('ZSX', ['rz', 'sx'])], + ('XYX', ['rx', 'ry']), ('ZSX', ['rz', 'sx']), + ('ZSXX', ['rz', 'sx', 'x'])], kak_gates=[(CXGate(), 'cx'), (CZGate(), 'cz'), (iSwapGate(), 'iswap'), (RXXGate(np.pi / 2), 'rxx')], name='test_euler_basis_selection_{seed}_{euler_bases[0]}_{kak_gates[1]}') @@ -545,7 +994,121 @@ def test_euler_basis_selection(self, euler_bases, kak_gates, seed): decomposition_basis.issubset(requested_basis)) -# FIXME: need to write tests for the approximate decompositions +@ddt +class TestTwoQubitDecomposeApprox(CheckDecompositions): + """Smoke tests for automatically-chosen approximate decompositions""" + + def check_approx_decomposition(self, target_unitary, decomposer, num_basis_uses): + """Check approx decomposition for a particular target""" + self.assertEqual(decomposer.num_basis_gates(target_unitary), num_basis_uses) + decomp_circuit = decomposer(target_unitary) + self.assertEqual(num_basis_uses, decomp_circuit.count_ops().get('unitary', 0)) + + @combine(seed=range(10), name='seed_{seed}') + def test_approx_supercontrolled_decompose_phase_0_use_random(self, seed, delta=0.01): + """Approx decomposition supercontrolled basis, random target (0 basis uses) (seed={seed}) + """ + state = np.random.default_rng(seed) + basis_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_phase = state.random() * 2 * np.pi + basis_b = 0.4 # how to safely randomize? + basis_unitary = np.exp(1j * basis_phase) * basis_k1 @ Ud(np.pi / 4, basis_b, 0) @ basis_k2 + decomposer = TwoQubitBasisDecomposer(UnitaryGate(basis_unitary), basis_fidelity=0.99) + + tgt_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_phase = state.random() * 2 * np.pi + d1, d2, d3 = state.random(size=3) * delta + tgt_unitary = np.exp(1j * tgt_phase) * tgt_k1 @ Ud(d1, d2, d3) @ tgt_k2 + self.check_approx_decomposition(tgt_unitary, decomposer, num_basis_uses=0) + + @combine(seed=range(10), name='seed_{seed}') + def test_approx_supercontrolled_decompose_phase_1_use_random(self, seed, delta=0.01): + """Approximate decomposition supercontrolled basis, random tgt (1 basis uses) (seed={seed}) + """ + state = np.random.default_rng(seed) + basis_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_phase = state.random() * 2 * np.pi + basis_b = 0.4 # how to safely randomize? + basis_unitary = np.exp(1j * basis_phase) * basis_k1 @ Ud(np.pi / 4, basis_b, 0) @ basis_k2 + decomposer = TwoQubitBasisDecomposer(UnitaryGate(basis_unitary), basis_fidelity=0.99) + + tgt_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_phase = state.random() * 2 * np.pi + d1, d2, d3 = state.random(size=3) * delta + tgt_unitary = np.exp(1j * tgt_phase) * tgt_k1 @ Ud(np.pi / 4-d1, basis_b+d2, d3) @ tgt_k2 + self.check_approx_decomposition(tgt_unitary, decomposer, num_basis_uses=1) + + @combine(seed=range(10), name='seed_{seed}') + def test_approx_supercontrolled_decompose_phase_2_use_random(self, seed, delta=0.01): + """Approximate decomposition supercontrolled basis, random tgt (2 basis uses) (seed={seed}) + """ + state = np.random.default_rng(seed) + basis_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_phase = state.random() * 2 * np.pi + basis_b = 0.4 # how to safely randomize? + basis_unitary = np.exp(1j * basis_phase) * basis_k1 @ Ud(np.pi / 4, basis_b, 0) @ basis_k2 + decomposer = TwoQubitBasisDecomposer(UnitaryGate(basis_unitary), basis_fidelity=0.99) + + tgt_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_phase = state.random() * 2 * np.pi + tgt_a, tgt_b = 0.3, 0.2 # how to safely randomize? + d1, d2, d3 = state.random(size=3) * delta + tgt_unitary = np.exp(1j * tgt_phase) * tgt_k1 @ Ud(tgt_a+d1, tgt_b+d2, d3) @ tgt_k2 + self.check_approx_decomposition(tgt_unitary, decomposer, num_basis_uses=2) + + @combine(seed=range(10), name='seed_{seed}') + def test_approx_supercontrolled_decompose_phase_3_use_random(self, seed, delta=0.01): + """Approximate decomposition supercontrolled basis, random tgt (3 basis uses) (seed={seed}) + """ + state = np.random.default_rng(seed) + basis_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + basis_phase = state.random() * 2 * np.pi + basis_b = state.random() * np.pi / 4 + basis_unitary = np.exp(1j * basis_phase) * basis_k1 @ Ud(np.pi / 4, basis_b, 0) @ basis_k2 + decomposer = TwoQubitBasisDecomposer(UnitaryGate(basis_unitary), basis_fidelity=0.99) + + tgt_k1 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_k2 = np.kron(random_unitary(2, seed=state).data, random_unitary(2, seed=state).data) + tgt_phase = state.random() * 2 * np.pi + tgt_a, tgt_b, tgt_c = 0.5, 0.4, 0.3 + d1, d2, d3 = state.random(size=3) * delta + tgt_unitary = np.exp(1j * tgt_phase) * tgt_k1 @ Ud(tgt_a+d1, tgt_b+d2, tgt_c+d3) @ tgt_k2 + self.check_approx_decomposition(tgt_unitary, decomposer, num_basis_uses=3) + + +class TestDecomposeProductRaises(QiskitTestCase): + """Check that exceptions are raised when 2q matrix is not a product of 1q unitaries""" + def test_decompose_two_qubit_product_gate_detr_too_small(self): + """Check that exception raised for too-small right component""" + kl = np.eye(2) + kr = 0.05*np.eye(2) + klkr = np.kron(kl, kr) + with self.assertRaises(QiskitError) as exc: + decompose_two_qubit_product_gate(klkr) + self.assertIn("detR <", exc.exception.message) + + def test_decompose_two_qubit_product_gate_detl_too_small(self): + """Check that exception raised for too-small left component""" + kl = np.array([[1, 0], [0, 0]]) + kr = np.eye(2) + klkr = np.kron(kl, kr) + with self.assertRaises(QiskitError) as exc: + decompose_two_qubit_product_gate(klkr) + self.assertIn("detL <", exc.exception.message) + + def test_decompose_two_qubit_product_gate_not_product(self): + """Check that exception raised for non-product unitary""" + klkr = Ud(1.E-6, 0, 0) + with self.assertRaises(QiskitError) as exc: + decompose_two_qubit_product_gate(klkr) + self.assertIn("decomposition failed", exc.exception.message) if __name__ == '__main__': diff --git a/test/python/transpiler/test_basis_translator.py b/test/python/transpiler/test_basis_translator.py index c1a21fccd2c9..ec6a9c732e9c 100644 --- a/test/python/transpiler/test_basis_translator.py +++ b/test/python/transpiler/test_basis_translator.py @@ -796,7 +796,8 @@ def test_global_phase(self): expected.global_phase = circ_angle - gate_angle / 2 expected_dag = circuit_to_dag(expected) self.assertEqual(out_dag, expected_dag) - self.assertEqual(float(out_dag.global_phase), float(expected_dag.global_phase)) + self.assertAlmostEqual(float(out_dag.global_phase), float(expected_dag.global_phase), + places=14) self.assertEqual(Operator(dag_to_circuit(out_dag)), Operator(expected)) def test_condition_set_substitute_node(self): diff --git a/test/python/transpiler/test_commutative_cancellation.py b/test/python/transpiler/test_commutative_cancellation.py index 6a7434dc0cb9..3af7d203cbb7 100644 --- a/test/python/transpiler/test_commutative_cancellation.py +++ b/test/python/transpiler/test_commutative_cancellation.py @@ -417,7 +417,8 @@ def test_commutative_circuit3(self): expected.append(RZGate(np.pi * 2 / 3), [qr[3]]) expected.cx(qr[2], qr[1]) - self.assertEqual(expected, new_circuit) + self.assertEqual(expected, new_circuit, + msg=f'expected:\n{expected}\nnew_circuit:\n{new_circuit}') def test_cnot_cascade(self): """ diff --git a/test/python/transpiler/test_optimize_1q_decomposition.py b/test/python/transpiler/test_optimize_1q_decomposition.py index b59c47bcb72b..6cce7a03dca1 100644 --- a/test/python/transpiler/test_optimize_1q_decomposition.py +++ b/test/python/transpiler/test_optimize_1q_decomposition.py @@ -346,7 +346,7 @@ def test_euler_decomposition_worse(self): result = passmanager.run(circuit) # decomposition of circuit will result in 3 gates instead of 2 # assert optimization pass doesn't use it. - self.assertEqual(result, circuit) + self.assertEqual(circuit, result, f"Circuit:\n{circuit}\nResult:\n{result}") def test_optimize_u_to_phase_gate(self): """U(0, 0, pi/4) -> p(pi/4). Basis [p, sx].""" @@ -363,7 +363,8 @@ def test_optimize_u_to_phase_gate(self): passmanager.append(Optimize1qGatesDecomposition(basis)) result = passmanager.run(circuit) - self.assertEqual(expected, result) + msg = f"expected:\n{expected}\nresult:\n{result}" + self.assertEqual(expected, result, msg=msg) def test_optimize_u_to_p_sx_p(self): """U(pi/2, 0, pi/4) -> p(-pi/4)-sx-p(p/2). Basis [p, sx].""" @@ -382,7 +383,8 @@ def test_optimize_u_to_p_sx_p(self): passmanager.append(Optimize1qGatesDecomposition(basis)) result = passmanager.run(circuit) - self.assertEqual(expected, result) + msg = f"expected:\n{expected}\nresult:\n{result}" + self.assertEqual(expected, result, msg=msg) def test_optimize_u3_to_u1(self): """U3(0, 0, pi/4) -> U1(pi/4). Basis [u1, u2, u3].""" @@ -399,7 +401,8 @@ def test_optimize_u3_to_u1(self): passmanager.append(Optimize1qGatesDecomposition(basis)) result = passmanager.run(circuit) - self.assertEqual(expected, result) + msg = f"expected:\n{expected}\nresult:\n{result}" + self.assertEqual(expected, result, msg=msg) def test_optimize_u3_to_u2(self): """U3(pi/2, 0, pi/4) -> U2(0, pi/4). Basis [u1, u2, u3].""" @@ -415,8 +418,25 @@ def test_optimize_u3_to_u2(self): passmanager.append(BasisTranslator(sel, basis)) passmanager.append(Optimize1qGatesDecomposition(basis)) result = passmanager.run(circuit) - self.assertEqual(expected, result) + msg = f"expected:\n{expected}\nresult:\n{result}" + self.assertEqual(expected, result, msg=msg) + + def test_y_simplification_rz_sx_x(self): + """Test that a y gate gets decomposed to x-zx with ibmq basis.""" + qc = QuantumCircuit(1) + qc.y(0) + basis = ["id", "rz", "sx", "x", "cx"] + passmanager = PassManager() + passmanager.append(BasisTranslator(sel, basis)) + passmanager.append(Optimize1qGatesDecomposition(basis)) + result = passmanager.run(qc) + expected = QuantumCircuit(1) + expected.x(0) + expected.rz(-np.pi, 0) + expected.global_phase += np.pi + msg = f"expected:\n{expected}\nresult:\n{result}" + self.assertEqual(expected, result, msg=msg) if __name__ == '__main__':