# This code is part of Qiskit.
#
# (C) Copyright IBM 2017, 2019.
#
# This code is licensed under the Apache License, Version 2.0. You may
# obtain a copy of this license in the LICENSE.txt file in the root directory
# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.
#
# Any modifications or derivative works of this code must retain this
# copyright notice, and modified files need to carry a notice indicating
# that they have been altered from the originals.
"""
Matrix Operator class.
"""
from __future__ import annotations
import copy
import re
from numbers import Number
from typing import TYPE_CHECKING
import numpy as np
from qiskit.circuit.instruction import Instruction
from qiskit.circuit.library.standard_gates import HGate, IGate, SGate, TGate, XGate, YGate, ZGate
from qiskit.circuit.operation import Operation
from qiskit.circuit.quantumcircuit import QuantumCircuit
from qiskit.exceptions import QiskitError
from qiskit.quantum_info.operators.base_operator import BaseOperator
from qiskit.quantum_info.operators.linear_op import LinearOp
from qiskit.quantum_info.operators.mixins import generate_apidocs
from qiskit.quantum_info.operators.predicates import is_unitary_matrix, matrix_equal
if TYPE_CHECKING:
from qiskit.transpiler.layout import Layout
[문서]class Operator(LinearOp):
r"""Matrix operator class
This represents a matrix operator :math:`M` that will
:meth:`~Statevector.evolve` a :class:`Statevector` :math:`|\psi\rangle`
by matrix-vector multiplication
.. math::
|\psi\rangle \mapsto M|\psi\rangle,
and will :meth:`~DensityMatrix.evolve` a :class:`DensityMatrix` :math:`\rho`
by left and right multiplication
.. math::
\rho \mapsto M \rho M^\dagger.
"""
def __init__(
self,
data: QuantumCircuit | Operation | BaseOperator | np.ndarray,
input_dims: tuple | None = None,
output_dims: tuple | None = None,
):
"""Initialize an operator object.
Args:
data (QuantumCircuit or Operation or BaseOperator or matrix):
data to initialize operator.
input_dims (tuple): the input subsystem dimensions.
[Default: None]
output_dims (tuple): the output subsystem dimensions.
[Default: None]
Raises:
QiskitError: if input data cannot be initialized as an operator.
Additional Information:
If the input or output dimensions are None, they will be
automatically determined from the input data. If the input data is
a Numpy array of shape (2**N, 2**N) qubit systems will be used. If
the input operator is not an N-qubit operator, it will assign a
single subsystem with dimension specified by the shape of the input.
"""
op_shape = None
if isinstance(data, (list, np.ndarray)):
# Default initialization from list or numpy array matrix
self._data = np.asarray(data, dtype=complex)
elif isinstance(data, (QuantumCircuit, Operation)):
# If the input is a Terra QuantumCircuit or Operation we
# perform a simulation to construct the unitary operator.
# This will only work if the circuit or instruction can be
# defined in terms of unitary gate instructions which have a
# 'to_matrix' method defined. Any other instructions such as
# conditional gates, measure, or reset will cause an
# exception to be raised.
self._data = self._init_instruction(data).data
elif hasattr(data, "to_operator"):
# If the data object has a 'to_operator' attribute this is given
# higher preference than the 'to_matrix' method for initializing
# an Operator object.
data = data.to_operator()
self._data = data.data
op_shape = data._op_shape
elif hasattr(data, "to_matrix"):
# If no 'to_operator' attribute exists we next look for a
# 'to_matrix' attribute to a matrix that will be cast into
# a complex numpy matrix.
self._data = np.asarray(data.to_matrix(), dtype=complex)
else:
raise QiskitError("Invalid input data format for Operator")
super().__init__(
op_shape=op_shape,
input_dims=input_dims,
output_dims=output_dims,
shape=self._data.shape,
)
def __array__(self, dtype=None):
if dtype:
return np.asarray(self.data, dtype=dtype)
return self.data
def __repr__(self):
prefix = "Operator("
pad = len(prefix) * " "
return "{}{},\n{}input_dims={}, output_dims={})".format(
prefix,
np.array2string(self.data, separator=", ", prefix=prefix),
pad,
self.input_dims(),
self.output_dims(),
)
def __eq__(self, other):
"""Test if two Operators are equal."""
if not super().__eq__(other):
return False
return np.allclose(self.data, other.data, rtol=self.rtol, atol=self.atol)
@property
def data(self):
"""The underlying Numpy array."""
return self._data
@property
def settings(self):
"""Return operator settings."""
return {
"data": self._data,
"input_dims": self.input_dims(),
"output_dims": self.output_dims(),
}
[문서] @classmethod
def from_label(cls, label: str) -> Operator:
"""Return a tensor product of single-qubit operators.
Args:
label (string): single-qubit operator string.
Returns:
Operator: The N-qubit operator.
Raises:
QiskitError: if the label contains invalid characters, or the
length of the label is larger than an explicitly
specified num_qubits.
Additional Information:
The labels correspond to the single-qubit matrices:
'I': [[1, 0], [0, 1]]
'X': [[0, 1], [1, 0]]
'Y': [[0, -1j], [1j, 0]]
'Z': [[1, 0], [0, -1]]
'H': [[1, 1], [1, -1]] / sqrt(2)
'S': [[1, 0], [0 , 1j]]
'T': [[1, 0], [0, (1+1j) / sqrt(2)]]
'0': [[1, 0], [0, 0]]
'1': [[0, 0], [0, 1]]
'+': [[0.5, 0.5], [0.5 , 0.5]]
'-': [[0.5, -0.5], [-0.5 , 0.5]]
'r': [[0.5, -0.5j], [0.5j , 0.5]]
'l': [[0.5, 0.5j], [-0.5j , 0.5]]
"""
# Check label is valid
label_mats = {
"I": IGate().to_matrix(),
"X": XGate().to_matrix(),
"Y": YGate().to_matrix(),
"Z": ZGate().to_matrix(),
"H": HGate().to_matrix(),
"S": SGate().to_matrix(),
"T": TGate().to_matrix(),
"0": np.array([[1, 0], [0, 0]], dtype=complex),
"1": np.array([[0, 0], [0, 1]], dtype=complex),
"+": np.array([[0.5, 0.5], [0.5, 0.5]], dtype=complex),
"-": np.array([[0.5, -0.5], [-0.5, 0.5]], dtype=complex),
"r": np.array([[0.5, -0.5j], [0.5j, 0.5]], dtype=complex),
"l": np.array([[0.5, 0.5j], [-0.5j, 0.5]], dtype=complex),
}
if re.match(r"^[IXYZHST01rl\-+]+$", label) is None:
raise QiskitError("Label contains invalid characters.")
# Initialize an identity matrix and apply each gate
num_qubits = len(label)
op = Operator(np.eye(2**num_qubits, dtype=complex))
for qubit, char in enumerate(reversed(label)):
if char != "I":
op = op.compose(label_mats[char], qargs=[qubit])
return op
[문서] def apply_permutation(self, perm: list, front: bool = False) -> Operator:
"""Modifies operator's data by composing it with a permutation.
Args:
perm (list): permutation pattern, describing which qubits
occupy the positions 0, 1, 2, etc. after applying the permutation.
front (bool): When set to ``True`` the permutation is applied before the
operator, when set to ``False`` the permutation is applied after the
operator.
Returns:
Operator: The modified operator.
Raises:
QiskitError: if the size of the permutation pattern does not match the
dimensions of the operator.
"""
# See https://github.com/Qiskit/qiskit-terra/pull/9403 for the math
# behind the following code.
inv_perm = np.argsort(perm)
raw_shape_l = self._op_shape.dims_l()
n_dims_l = len(raw_shape_l)
raw_shape_r = self._op_shape.dims_r()
n_dims_r = len(raw_shape_r)
if front:
# The permutation is applied first, the operator is applied after;
# however, in terms of matrices, we compute [O][P].
if len(perm) != n_dims_r:
raise QiskitError(
"The size of the permutation pattern does not match dimensions of the operator."
)
# shape: original on left, permuted on right
shape_l = self._op_shape.dims_l()
shape_r = tuple(raw_shape_r[n_dims_r - n - 1] for n in reversed(perm))
# axes order: id on left, inv-permuted on right
axes_l = tuple(x for x in range(self._op_shape._num_qargs_l))
axes_r = tuple(self._op_shape._num_qargs_l + x for x in (np.argsort(perm[::-1]))[::-1])
# updated shape: original on left, permuted on right
new_shape_l = self._op_shape.dims_l()
new_shape_r = tuple(raw_shape_r[n_dims_r - n - 1] for n in reversed(inv_perm))
else:
# The operator is applied first, the permutation is applied after;
# however, in terms of matrices, we compute [P][O].
if len(perm) != n_dims_l:
raise QiskitError(
"The size of the permutation pattern does not match dimensions of the operator."
)
# shape: inv-permuted on left, original on right
shape_l = tuple(raw_shape_l[n_dims_l - n - 1] for n in reversed(inv_perm))
shape_r = self._op_shape.dims_r()
# axes order: permuted on left, id on right
axes_l = tuple((np.argsort(inv_perm[::-1]))[::-1])
axes_r = tuple(
self._op_shape._num_qargs_l + x for x in range(self._op_shape._num_qargs_r)
)
# updated shape: permuted on left, original on right
new_shape_l = tuple(raw_shape_l[n_dims_l - n - 1] for n in reversed(perm))
new_shape_r = self._op_shape.dims_r()
# Computing the new operator
split_shape = shape_l + shape_r
axes_order = axes_l + axes_r
new_mat = (
self._data.reshape(split_shape).transpose(axes_order).reshape(self._op_shape.shape)
)
new_op = Operator(new_mat, input_dims=new_shape_r, output_dims=new_shape_l)
return new_op
[문서] @classmethod
def from_circuit(
cls,
circuit: QuantumCircuit,
ignore_set_layout: bool = False,
layout: Layout | None = None,
final_layout: Layout | None = None,
) -> Operator:
"""Create a new Operator object from a :class:`.QuantumCircuit`
While a :class:`~.QuantumCircuit` object can passed directly as ``data``
to the class constructor this provides no options on how the circuit
is used to create an :class:`.Operator`. This constructor method lets
you control how the :class:`.Operator` is created so it can be adjusted
for a particular use case.
By default this constructor method will permute the qubits based on a
configured initial layout (i.e. after it was transpiled). It also
provides an option to manually provide a :class:`.Layout` object
directly.
Args:
circuit (QuantumCircuit): The :class:`.QuantumCircuit` to create an Operator
object from.
ignore_set_layout (bool): When set to ``True`` if the input ``circuit``
has a layout set it will be ignored
layout (Layout): If specified this kwarg can be used to specify a
particular layout to use to permute the qubits in the created
:class:`.Operator`. If this is specified it will be used instead
of a layout contained in the ``circuit`` input. If specified
the virtual bits in the :class:`~.Layout` must be present in the
``circuit`` input.
final_layout (Layout): If specified this kwarg can be used to represent the
output permutation caused by swap insertions during the routing stage
of the transpiler.
Returns:
Operator: An operator representing the input circuit
"""
dimension = 2**circuit.num_qubits
op = cls(np.eye(dimension))
if layout is None:
if not ignore_set_layout:
layout = getattr(circuit, "_layout", None)
else:
from qiskit.transpiler.layout import TranspileLayout # pylint: disable=cyclic-import
layout = TranspileLayout(
initial_layout=layout,
input_qubit_mapping={qubit: index for index, qubit in enumerate(circuit.qubits)},
)
if final_layout is None:
if not ignore_set_layout and layout is not None:
final_layout = getattr(layout, "final_layout", None)
qargs = None
# If there was a layout specified (either from the circuit
# or via user input) use that to set qargs to permute qubits
# based on that layout
if layout is not None:
physical_to_virtual = layout.initial_layout.get_physical_bits()
qargs = [
layout.input_qubit_mapping[physical_to_virtual[physical_bit]]
for physical_bit in range(len(physical_to_virtual))
]
# Convert circuit to an instruction
instruction = circuit.to_instruction()
op._append_instruction(instruction, qargs=qargs)
# If final layout is set permute output indices based on layout
if final_layout is not None:
perm_pattern = [final_layout._v2p[v] for v in circuit.qubits]
op = op.apply_permutation(perm_pattern, front=False)
return op
[문서] def is_unitary(self, atol=None, rtol=None):
"""Return True if operator is a unitary matrix."""
if atol is None:
atol = self.atol
if rtol is None:
rtol = self.rtol
return is_unitary_matrix(self._data, rtol=rtol, atol=atol)
[문서] def to_operator(self) -> Operator:
"""Convert operator to matrix operator class"""
return self
[문서] def to_instruction(self):
"""Convert to a UnitaryGate instruction."""
# pylint: disable=cyclic-import
from qiskit.extensions.unitary import UnitaryGate
return UnitaryGate(self.data)
[문서] def conjugate(self):
# Make a shallow copy and update array
ret = copy.copy(self)
ret._data = np.conj(self._data)
return ret
[문서] def transpose(self):
# Make a shallow copy and update array
ret = copy.copy(self)
ret._data = np.transpose(self._data)
ret._op_shape = self._op_shape.transpose()
return ret
[문서] def compose(self, other: Operator, qargs: list | None = None, front: bool = False) -> Operator:
if qargs is None:
qargs = getattr(other, "qargs", None)
if not isinstance(other, Operator):
other = Operator(other)
# Validate dimensions are compatible and return the composed
# operator dimensions
new_shape = self._op_shape.compose(other._op_shape, qargs, front)
input_dims = new_shape.dims_r()
output_dims = new_shape.dims_l()
# Full composition of operators
if qargs is None:
if front:
# Composition self * other
data = np.dot(self._data, other.data)
else:
# Composition other * self
data = np.dot(other.data, self._data)
ret = Operator(data, input_dims, output_dims)
ret._op_shape = new_shape
return ret
# Compose with other on subsystem
num_qargs_l, num_qargs_r = self._op_shape.num_qargs
if front:
num_indices = num_qargs_r
shift = num_qargs_l
right_mul = True
else:
num_indices = num_qargs_l
shift = 0
right_mul = False
# Reshape current matrix
# Note that we must reverse the subsystem dimension order as
# qubit 0 corresponds to the right-most position in the tensor
# product, which is the last tensor wire index.
tensor = np.reshape(self.data, self._op_shape.tensor_shape)
mat = np.reshape(other.data, other._op_shape.tensor_shape)
indices = [num_indices - 1 - qubit for qubit in qargs]
final_shape = [int(np.prod(output_dims)), int(np.prod(input_dims))]
data = np.reshape(
Operator._einsum_matmul(tensor, mat, indices, shift, right_mul), final_shape
)
ret = Operator(data, input_dims, output_dims)
ret._op_shape = new_shape
return ret
[문서] def power(self, n: float) -> Operator:
"""Return the matrix power of the operator.
Args:
n (float): the power to raise the matrix to.
Returns:
Operator: the resulting operator ``O ** n``.
Raises:
QiskitError: if the input and output dimensions of the operator
are not equal.
"""
if self.input_dims() != self.output_dims():
raise QiskitError("Can only power with input_dims = output_dims.")
ret = copy.copy(self)
ret._data = np.linalg.matrix_power(self.data, n)
return ret
[문서] def tensor(self, other: Operator) -> Operator:
if not isinstance(other, Operator):
other = Operator(other)
return self._tensor(self, other)
[문서] def expand(self, other: Operator) -> Operator:
if not isinstance(other, Operator):
other = Operator(other)
return self._tensor(other, self)
@classmethod
def _tensor(cls, a, b):
ret = copy.copy(a)
ret._op_shape = a._op_shape.tensor(b._op_shape)
ret._data = np.kron(a.data, b.data)
return ret
def _add(self, other, qargs=None):
"""Return the operator self + other.
If ``qargs`` are specified the other operator will be added
assuming it is identity on all other subsystems.
Args:
other (Operator): an operator object.
qargs (None or list): optional subsystems to add on
(Default: None)
Returns:
Operator: the operator self + other.
Raises:
QiskitError: if other is not an operator, or has incompatible
dimensions.
"""
# pylint: disable=cyclic-import
from qiskit.quantum_info.operators.scalar_op import ScalarOp
if qargs is None:
qargs = getattr(other, "qargs", None)
if not isinstance(other, Operator):
other = Operator(other)
self._op_shape._validate_add(other._op_shape, qargs)
other = ScalarOp._pad_with_identity(self, other, qargs)
ret = copy.copy(self)
ret._data = self.data + other.data
return ret
def _multiply(self, other):
"""Return the operator self * other.
Args:
other (complex): a complex number.
Returns:
Operator: the operator other * self.
Raises:
QiskitError: if other is not a valid complex number.
"""
if not isinstance(other, Number):
raise QiskitError("other is not a number")
ret = copy.copy(self)
ret._data = other * self._data
return ret
[문서] def equiv(self, other: Operator, rtol: float | None = None, atol: float | None = None) -> bool:
"""Return True if operators are equivalent up to global phase.
Args:
other (Operator): an operator object.
rtol (float): relative tolerance value for comparison.
atol (float): absolute tolerance value for comparison.
Returns:
bool: True if operators are equivalent up to global phase.
"""
if not isinstance(other, Operator):
try:
other = Operator(other)
except QiskitError:
return False
if self.dim != other.dim:
return False
if atol is None:
atol = self.atol
if rtol is None:
rtol = self.rtol
return matrix_equal(self.data, other.data, ignore_phase=True, rtol=rtol, atol=atol)
[문서] def reverse_qargs(self) -> Operator:
r"""Return an Operator with reversed subsystem ordering.
For a tensor product operator this is equivalent to reversing
the order of tensor product subsystems. For an operator
:math:`A = A_{n-1} \otimes ... \otimes A_0`
the returned operator will be
:math:`A_0 \otimes ... \otimes A_{n-1}`.
Returns:
Operator: the operator with reversed subsystem order.
"""
ret = copy.copy(self)
axes = tuple(range(self._op_shape._num_qargs_l - 1, -1, -1))
axes = axes + tuple(len(axes) + i for i in axes)
ret._data = np.reshape(
np.transpose(np.reshape(self.data, self._op_shape.tensor_shape), axes),
self._op_shape.shape,
)
ret._op_shape = self._op_shape.reverse()
return ret
[문서] def to_matrix(self):
"""Convert operator to NumPy matrix."""
return self.data
@classmethod
def _einsum_matmul(cls, tensor, mat, indices, shift=0, right_mul=False):
"""Perform a contraction using Numpy.einsum
Args:
tensor (np.array): a vector or matrix reshaped to a rank-N tensor.
mat (np.array): a matrix reshaped to a rank-2M tensor.
indices (list): tensor indices to contract with mat.
shift (int): shift for indices of tensor to contract [Default: 0].
right_mul (bool): if True right multiply tensor by mat
(else left multiply) [Default: False].
Returns:
Numpy.ndarray: the matrix multiplied rank-N tensor.
Raises:
QiskitError: if mat is not an even rank tensor.
"""
rank = tensor.ndim
rank_mat = mat.ndim
if rank_mat % 2 != 0:
raise QiskitError("Contracted matrix must have an even number of indices.")
# Get einsum indices for tensor
indices_tensor = list(range(rank))
for j, index in enumerate(indices):
indices_tensor[index + shift] = rank + j
# Get einsum indices for mat
mat_contract = list(reversed(range(rank, rank + len(indices))))
mat_free = [index + shift for index in reversed(indices)]
if right_mul:
indices_mat = mat_contract + mat_free
else:
indices_mat = mat_free + mat_contract
return np.einsum(tensor, indices_tensor, mat, indices_mat)
@classmethod
def _init_instruction(cls, instruction):
"""Convert a QuantumCircuit or Operation to an Operator."""
# Initialize an identity operator of the correct size of the circuit
if hasattr(instruction, "__array__"):
return Operator(np.array(instruction, dtype=complex))
dimension = 2**instruction.num_qubits
op = Operator(np.eye(dimension))
# Convert circuit to an instruction
if isinstance(instruction, QuantumCircuit):
instruction = instruction.to_instruction()
op._append_instruction(instruction)
return op
@classmethod
def _instruction_to_matrix(cls, obj):
"""Return Operator for instruction if defined or None otherwise."""
# Note: to_matrix() is not a required method for Operations, so for now
# we do not allow constructing matrices for general Operations.
# However, for backward compatibility we need to support constructing matrices
# for Cliffords, which happen to have a to_matrix() method.
# pylint: disable=cyclic-import
from qiskit.quantum_info import Clifford
if not isinstance(obj, (Instruction, Clifford)):
raise QiskitError("Input is neither an Instruction nor Clifford.")
mat = None
if hasattr(obj, "to_matrix"):
# If instruction is a gate first we see if it has a
# `to_matrix` definition and if so use that.
try:
mat = obj.to_matrix()
except QiskitError:
pass
return mat
def _append_instruction(self, obj, qargs=None):
"""Update the current Operator by apply an instruction."""
from qiskit.circuit.barrier import Barrier
from .scalar_op import ScalarOp
mat = self._instruction_to_matrix(obj)
if mat is not None:
# Perform the composition and inplace update the current state
# of the operator
op = self.compose(mat, qargs=qargs)
self._data = op.data
elif isinstance(obj, Barrier):
return
else:
# If the instruction doesn't have a matrix defined we use its
# circuit decomposition definition if it exists, otherwise we
# cannot compose this gate and raise an error.
if obj.definition is None:
raise QiskitError(f"Cannot apply Operation: {obj.name}")
if not isinstance(obj.definition, QuantumCircuit):
raise QiskitError(
'Operation "{}" '
"definition is {} but expected QuantumCircuit.".format(
obj.name, type(obj.definition)
)
)
if obj.definition.global_phase:
dimension = 2**obj.num_qubits
op = self.compose(
ScalarOp(dimension, np.exp(1j * float(obj.definition.global_phase))),
qargs=qargs,
)
self._data = op.data
flat_instr = obj.definition
bit_indices = {
bit: index
for bits in [flat_instr.qubits, flat_instr.clbits]
for index, bit in enumerate(bits)
}
for instruction in flat_instr:
if instruction.clbits:
raise QiskitError(
"Cannot apply operation with classical bits:"
f" {instruction.operation.name}"
)
# Get the integer position of the flat register
if qargs is None:
new_qargs = [bit_indices[tup] for tup in instruction.qubits]
else:
new_qargs = [qargs[bit_indices[tup]] for tup in instruction.qubits]
self._append_instruction(instruction.operation, qargs=new_qargs)
# Update docstrings for API docs
generate_apidocs(Operator)