Operator.compose should allow indices

qiskit/aqua/operators/evolutions/evolved_op.py

@@ -86,13 +86,23 @@ def tensor(self, other: OperatorBase) -> OperatorBase:
 
         return TensoredOp([self, other])
 
-    def compose(self, other: OperatorBase) -> OperatorBase:
-        other = self._check_zero_for_composition_and_expand(other)
+    def _expand_dim(self, num_qubits: int) -> OperatorBase:
+        from qiskit.aqua.operators import I
+        return self.tensor(I ^ num_qubits)
 
+    def permute(self, permutation: List[int]) -> OperatorBase:
+        return EvolvedOp(self.primitive.permute(permutation), coeff=self.coeff)  # type: ignore
+
+    def compose(self, other: OperatorBase,
+                permutation: Optional[List[int]] = None, front: bool = False) -> OperatorBase:
+
+        new_self, other = self._expand_shorter_operator_and_permute(other, permutation)
+        if front:
+            return other.compose(new_self)
         if isinstance(other, ComposedOp):
-            return ComposedOp([self] + other.oplist)  # type: ignore
+            return ComposedOp([new_self] + other.oplist)
 
-        return ComposedOp([self, other])
+        return ComposedOp([new_self, other])
 
     def __str__(self) -> str:
         prim_str = str(self.primitive)

qiskit/aqua/operators/list_ops/composed_op.py

@@ -12,7 +12,7 @@
 
 """ ComposedOp Class """
 
-from typing import List, Union, cast
+from typing import List, Union, cast, Optional
 from functools import reduce, partial
 import numpy as np
 
@@ -82,24 +82,31 @@ def to_circuit(self) -> QuantumCircuit:
     def adjoint(self) -> OperatorBase:
         return ComposedOp([op.adjoint() for op in reversed(self.oplist)], coeff=self.coeff)
 
-    def compose(self, other: OperatorBase) -> OperatorBase:
+    def compose(self, other: OperatorBase,
+                permutation: Optional[List[int]] = None, front: bool = False) -> OperatorBase:
+
+        new_self, other = self._expand_shorter_operator_and_permute(other, permutation)
+        new_self = cast(ComposedOp, new_self)
+
+        if front:
+            return other.compose(new_self)
         # Try composing with last element in list
         if isinstance(other, ComposedOp):
-            return ComposedOp(self.oplist + other.oplist, coeff=self.coeff * other.coeff)
+            return ComposedOp(new_self.oplist + other.oplist, coeff=new_self.coeff * other.coeff)
 
         # Try composing with last element of oplist. We only try
         # this if that last element isn't itself an
         # ComposedOp, so we can tell whether composing the
         # two elements directly worked. If it doesn't,
         # continue to the final return statement below, appending other to the oplist.
-        if not isinstance(self.oplist[-1], ComposedOp):
-            comp_with_last = self.oplist[-1].compose(other)
+        if not isinstance(new_self.oplist[-1], ComposedOp):
+            comp_with_last = new_self.oplist[-1].compose(other)
             # Attempt successful
             if not isinstance(comp_with_last, ComposedOp):
-                new_oplist = self.oplist[0:-1] + [comp_with_last]
-                return ComposedOp(new_oplist, coeff=self.coeff)
+                new_oplist = new_self.oplist[0:-1] + [comp_with_last]
+                return ComposedOp(new_oplist, coeff=new_self.coeff)
 
-        return ComposedOp(self.oplist + [other], coeff=self.coeff)
+        return ComposedOp(new_self.oplist + [other], coeff=new_self.coeff)
 
     def eval(self,
              front: Union[str, dict, np.ndarray,

qiskit/aqua/operators/list_ops/list_op.py

@@ -13,15 +13,18 @@
 """ ListOp Operator Class """
 
 from functools import reduce
-from typing import List, Union, Optional, Callable, Iterator, Set, Dict
+from typing import List, Union, Optional, Callable, Iterator, Set, Dict, cast
 from numbers import Number
 
 import numpy as np
 from scipy.sparse import spmatrix
 
-from qiskit.circuit import ParameterExpression
+from qiskit.circuit import QuantumCircuit, ParameterExpression
+
 from ..legacy.base_operator import LegacyBaseOperator
 from ..operator_base import OperatorBase
+from ... import AquaError
+from ...utils import arithmetic
 
 
 class ListOp(OperatorBase):
@@ -212,11 +215,61 @@ def tensorpower(self, other: int) -> Union[OperatorBase, int]:
         from .tensored_op import TensoredOp
         return TensoredOp([self] * other)
 
-    def compose(self, other: OperatorBase) -> OperatorBase:
+    def _expand_dim(self, num_qubits: int) -> 'ListOp':
+        return ListOp([op._expand_dim(num_qubits + self.num_qubits - op.num_qubits)
+                       for op in self.oplist], combo_fn=self.combo_fn, coeff=self.coeff)
+
+    def permute(self, permutation: List[int]) -> 'ListOp':
+        """Permute the qubits of the operator.
+
+        Args:
+            permutation: A list defining where each qubit should be permuted. The qubit at index
+                j should be permuted to position permutation[j].
+
+        Returns:
+            A new ListOp representing the permuted operator.
+
+        Raises:
+            AquaError: if indices do not define a new index for each qubit.
+        """
+        new_self = self
+        circuit_size = max(permutation) + 1
+
+        try:
+            if self.num_qubits != len(permutation):
+                raise AquaError("New index must be defined for each qubit of the operator.")
+        except ValueError:
+            raise AquaError("Permute is only possible if all operators in the ListOp have the "
+                            "same number of qubits.") from ValueError
+        if self.num_qubits < circuit_size:
+            # pad the operator with identities
+            new_self = self._expand_dim(circuit_size - self.num_qubits)
+        qc = QuantumCircuit(circuit_size)
+        # extend the indices to match the size of the circuit
+        permutation \
+            = list(filter(lambda x: x not in permutation, range(circuit_size))) + permutation
+
+        # decompose permutation into sequence of transpositions
+        transpositions = arithmetic.transpositions(permutation)
+        for trans in transpositions:
+            qc.swap(trans[0], trans[1])
+
+        from qiskit.aqua.operators import CircuitOp
+
+        return CircuitOp(qc.reverse_ops()) @ new_self @ CircuitOp(qc)
+
+    def compose(self, other: OperatorBase,
+                permutation: Optional[List[int]] = None, front: bool = False) -> OperatorBase:
+
+        new_self, other = self._expand_shorter_operator_and_permute(other, permutation)
+        new_self = cast(ListOp, new_self)
+
+        if front:
+            return other.compose(new_self)
         # Avoid circular dependency
         # pylint: disable=cyclic-import,import-outside-toplevel
         from .composed_op import ComposedOp
-        return ComposedOp([self, other])
+        return ComposedOp([new_self, other])
 
     def power(self, exponent: int) -> OperatorBase:
         if not isinstance(exponent, int) or exponent <= 0:
@@ -258,8 +311,8 @@ def to_spmatrix(self) -> Union[spmatrix, List[spmatrix]]:
             [op.to_spmatrix() for op in self.oplist]) * self.coeff  # type: ignore
 
     def eval(self,
-             front: Optional[Union[str, Dict[str, complex], 'OperatorBase']] = None
-             ) -> Union['OperatorBase', float, complex, list]:
+             front: Optional[Union[str, Dict[str, complex], OperatorBase]] = None
+             ) -> Union[OperatorBase, float, complex, list]:
         """
         Evaluate the Operator's underlying function, either on a binary string or another Operator.
         A square binary Operator can be defined as a function taking a binary function to another

qiskit/aqua/operators/list_ops/tensored_op.py

@@ -54,6 +54,16 @@ def num_qubits(self) -> int:
     def distributive(self) -> bool:
         return False
 
+    def _expand_dim(self, num_qubits: int) -> 'TensoredOp':
+        """Appends I ^ num_qubits to ``oplist``. Choice of PauliOp as
+        identity is arbitrary and can be substituted for other PrimitiveOp identity.
+
+        Returns:
+            TensoredOp expanded with identity operator.
+        """
+        from qiskit.aqua.operators import I
+        return TensoredOp(self.oplist + [I ^ num_qubits], coeff=self.coeff)
+
     def tensor(self, other: OperatorBase) -> OperatorBase:
         if isinstance(other, TensoredOp):
             return TensoredOp(self.oplist + other.oplist, coeff=self.coeff * other.coeff)

qiskit/aqua/operators/operator_base.py

@@ -12,7 +12,7 @@
 
 """ OperatorBase Class """
 
-from typing import Set, Union, Dict, Optional, List, cast
+from typing import Set, Union, Dict, Optional, List, cast, Tuple
 from numbers import Number
 from abc import ABC, abstractmethod
 import numpy as np
@@ -452,6 +452,32 @@ def assign_parameters(self,
         """
         raise NotImplementedError
 
+    @abstractmethod
+    def _expand_dim(self, num_qubits: int) -> 'OperatorBase':
+        """Expands the operator with identity operator of dimension 2**num_qubits.
+
+        Returns:
+            Operator corresponding to self.tensor(identity_operator), where dimension of identity
+            operator is 2 ** num_qubits.
+        """
+        raise NotImplementedError
+
+    @abstractmethod
+    def permute(self, permutation: List[int]) -> 'OperatorBase':
+        """Permutes the qubits of the operator.
+
+        Args:
+            permutation: A list defining where each qubit should be permuted. The qubit at index
+                j should be permuted to position permutation[j].
+
+        Returns:
+            A new OperatorBase containing the permuted operator.
+
+        Raises:
+            AquaError: if indices do not define a new index for each qubit.
+        """
+        raise NotImplementedError
+
     def bind_parameters(self,
                         param_dict: Dict[ParameterExpression,
                                          Union[Number,
@@ -501,6 +527,24 @@ def _get_param_dict_for_index(unrolled_dict: Dict[ParameterExpression, List[Numb
         """ Gets a single non-list-nested param_dict for a given list index from a nested one. """
         return {k: v[i] for (k, v) in unrolled_dict.items()}
 
+    def _expand_shorter_operator_and_permute(self, other: 'OperatorBase',
+                                             permutation: Optional[List[int]] = None) \
+            -> Tuple['OperatorBase', 'OperatorBase']:
+        if permutation is not None:
+            other = other.permute(permutation)
+        new_self = self
+        if not self.num_qubits == other.num_qubits:
+            # pylint: disable=cyclic-import,import-outside-toplevel
+            from .operator_globals import Zero
+            if other == Zero:
+                # Zero is special - we'll expand it to the correct qubit number.
+                other = Zero.__class__('0' * self.num_qubits)
+            elif other.num_qubits < self.num_qubits:
+                other = other._expand_dim(self.num_qubits - other.num_qubits)
+            elif other.num_qubits > self.num_qubits:
+                new_self = self._expand_dim(other.num_qubits - self.num_qubits)
+        return new_self, other
+
     # Composition
 
     def __matmul__(self, other: 'OperatorBase') -> 'OperatorBase':
@@ -515,7 +559,8 @@ def __matmul__(self, other: 'OperatorBase') -> 'OperatorBase':
         return self.compose(other)
 
     @abstractmethod
-    def compose(self, other: 'OperatorBase') -> 'OperatorBase':
+    def compose(self, other: 'OperatorBase',
+                permutation: Optional[List[int]] = None, front: bool = False) -> 'OperatorBase':
         r""" Return Operator Composition between self and other (linear algebra-style:
         A@B(x) = A(B(x))), overloaded by ``@``.
 
@@ -529,6 +574,8 @@ def compose(self, other: 'OperatorBase') -> 'OperatorBase':
 
         Args:
             other: The ``OperatorBase`` with which to compose self.
+            permutation: ``List[int]`` which defines permutation on other operator.
+            front: If front==True, return ``other.compose(self)``.
 
         Returns:
             An ``OperatorBase`` equivalent to the function composition of self and other.

qiskit/aqua/operators/primitive_ops/circuit_op.py

@@ -108,30 +108,37 @@ def tensor(self, other: OperatorBase) -> OperatorBase:
 
         return TensoredOp([self, other])
 
-    def compose(self, other: OperatorBase) -> OperatorBase:
-        other = self._check_zero_for_composition_and_expand(other)
+    def compose(self, other: OperatorBase,
+                permutation: Optional[List[int]] = None, front: bool = False) -> OperatorBase:
+
+        new_self, other = self._expand_shorter_operator_and_permute(other, permutation)
+        new_self = cast(CircuitOp, new_self)
+
+        if front:
+            return other.compose(new_self)
+        # ignore
         # pylint: disable=cyclic-import,import-outside-toplevel
         from ..operator_globals import Zero
         from ..state_fns import CircuitStateFn
         from .pauli_op import PauliOp
         from .matrix_op import MatrixOp
 
-        if other == Zero ^ self.num_qubits:
-            return CircuitStateFn(self.primitive, coeff=self.coeff)
+        if other == Zero ^ new_self.num_qubits:
+            return CircuitStateFn(new_self.primitive, coeff=new_self.coeff)
 
         if isinstance(other, (PauliOp, CircuitOp, MatrixOp)):
             other = other.to_circuit_op()
 
         if isinstance(other, (CircuitOp, CircuitStateFn)):
-            new_qc = other.primitive.compose(self.primitive)  # type: ignore
+            new_qc = other.primitive.compose(new_self.primitive)  # type: ignore
             if isinstance(other, CircuitStateFn):
                 return CircuitStateFn(new_qc,
                                       is_measurement=other.is_measurement,
-                                      coeff=self.coeff * other.coeff)
+                                      coeff=new_self.coeff * other.coeff)
             else:
-                return CircuitOp(new_qc, coeff=self.coeff * other.coeff)
+                return CircuitOp(new_qc, coeff=new_self.coeff * other.coeff)
 
-        return super().compose(other)
+        return super(CircuitOp, new_self).compose(other)
 
     def to_matrix(self, massive: bool = False) -> np.ndarray:
         if self.num_qubits > 16 and not massive:
@@ -219,6 +226,9 @@ def reduce(self) -> OperatorBase:
                     del self.primitive.data[i]  # type: ignore
         return self
 
+    def _expand_dim(self, num_qubits: int) -> 'CircuitOp':
+        return self.permute(list(range(num_qubits, num_qubits + self.num_qubits)))
+
     def permute(self, permutation: List[int]) -> 'CircuitOp':
         r"""
         Permute the qubits of the circuit.
@@ -230,5 +240,5 @@ def permute(self, permutation: List[int]) -> 'CircuitOp':
         Returns:
             A new CircuitOp containing the permuted circuit.
         """
-        new_qc = QuantumCircuit(self.num_qubits).compose(self.primitive, qubits=permutation)
+        new_qc = QuantumCircuit(max(permutation) + 1).compose(self.primitive, qubits=permutation)
         return CircuitOp(new_qc, coeff=self.coeff)