Implement the Tensor class (#1033)

* Implement the Tensor class

* Always wrap scalar numbers into numpy array

This matches what we did in the PolynomialTensor and avoids unnecessary
complicated code to deal with Number objects where otherwise one always
encounters an array.

* Add qiskit_nature.settings.tensor_wrapping

* Add release note

* fix: remove jupyter-execute blocks

The Jupyter Sphinx integration is currently unable to properly suppress
warnings resulting in prints to stderr [1]. This causes the docs to fail
building properly.
Qiskit Terra already removed the usage of `jupyter-execute` statements a
while ago [2], so we are following suit in this regard, too.

[1]: jupyter/jupyter-sphinx#178
[2]: Qiskit/qiskit#9346

* fix: rework settings.tensor_wrapping as settings.tensor_unwrapping

* docs: update documentaiton of PolynomialTensor now that Tensor exists

---------

Co-authored-by: mergify[bot] <37929162+mergify[bot]@users.noreply.github.com>

.pylintdict

@@ -242,6 +242,7 @@ intel
 intelvem
 interatomic
 internuclear
+interoperability
 interpretable
 ints
 io
@@ -354,6 +355,7 @@ natom
 nbasis
 nd
 ndarray
+ndim
 nelec
 neq
 networkx
@@ -583,6 +585,7 @@ uccsd
 ucrx
 ucry
 ucrz
+ufuncs
 uhf
 ulimit
 un

qiskit_nature/second_q/hamiltonians/lattices/hyper_cubic_lattice.py

@@ -1,6 +1,6 @@
 # This code is part of Qiskit.
 #
-# (C) Copyright IBM 2021, 2022.
+# (C) Copyright IBM 2021, 2023.
 #
 # 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
@@ -30,7 +30,7 @@ class HyperCubicLattice(Lattice):
     tuples of `size`, `edge_parameters`, and `boundary_conditions`.
     For example,
 
-    .. jupyter-execute::
+    .. code-block:: python
 
         from qiskit_nature.second_q.hamiltonians.lattices import (
             BoundaryCondition,

qiskit_nature/second_q/operators/__init__.py

@@ -1,6 +1,6 @@
 # This code is part of Qiskit.
 #
-# (C) Copyright IBM 2022.
+# (C) Copyright IBM 2022, 2023.
 #
 # 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
@@ -28,6 +28,7 @@
    VibrationalOp
    VibrationalIntegrals
    PolynomialTensor
+   Tensor
 
 Modules
 -------
@@ -46,6 +47,7 @@
 from .vibrational_integrals import VibrationalIntegrals
 from .polynomial_tensor import PolynomialTensor
 from .sparse_label_op import SparseLabelOp
+from .tensor import Tensor
 
 __all__ = [
     "ElectronicIntegrals",
@@ -55,4 +57,5 @@
     "VibrationalIntegrals",
     "PolynomialTensor",
     "SparseLabelOp",
+    "Tensor",
 ]

qiskit_nature/second_q/operators/electronic_integrals.py

@@ -1,6 +1,6 @@
 # This code is part of Qiskit.
 #
-# (C) Copyright IBM 2022.
+# (C) Copyright IBM 2022, 2023.
 #
 # 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
@@ -25,7 +25,8 @@
 from qiskit_nature.exceptions import QiskitNatureError
 import qiskit_nature.optionals as _optionals
 
-from .polynomial_tensor import ARRAY_TYPE, PolynomialTensor
+from .polynomial_tensor import PolynomialTensor
+from .tensor import Tensor
 from .tensor_ordering import (
     IndexType,
     find_index_order,
@@ -256,8 +257,8 @@ def alpha_beta(self) -> PolynomialTensor:
         if self.beta_alpha.is_empty():
             return self.beta_alpha
 
-        beta_alpha = cast(ARRAY_TYPE, self.beta_alpha["++--"])
-        alpha_beta = np.moveaxis(beta_alpha, (0, 1), (2, 3))
+        # TODO: remove extra-wrapping of Tensor once settings.tensor_unwrapping is removed
+        alpha_beta = Tensor(np.moveaxis(Tensor(self.beta_alpha["++--"]), (0, 1), (2, 3)))
         return PolynomialTensor({"++--": alpha_beta}, validate=False)
 
     @property
@@ -378,7 +379,7 @@ def _add(self, other: ElectronicIntegrals, qargs=None) -> ElectronicIntegrals:
     @classmethod
     def apply(
         cls,
-        function: Callable[..., np.ndarray | SparseArray | Number],
+        function: Callable[..., np.ndarray | SparseArray | complex],
         *operands: ElectronicIntegrals,
         validate: bool = True,
     ) -> ElectronicIntegrals:
@@ -568,9 +569,7 @@ def second_q_coeffs(self) -> PolynomialTensor:
 
         kron_one_body = np.zeros((2, 2))
         kron_two_body = np.zeros((2, 2, 2, 2))
-        kron_tensor = PolynomialTensor(
-            {"": cast(Number, 1.0), "+-": kron_one_body, "++--": kron_two_body}
-        )
+        kron_tensor = PolynomialTensor({"": 1.0, "+-": kron_one_body, "++--": kron_two_body})
 
         if beta_empty and beta_alpha_empty:
             kron_one_body[(0, 0)] = 1

qiskit_nature/second_q/operators/fermionic_op.py

@@ -1,6 +1,6 @@
 # This code is part of Qiskit.
 #
-# (C) Copyright IBM 2021, 2022.
+# (C) Copyright IBM 2021, 2023.
 #
 # 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
@@ -46,7 +46,7 @@ class FermionicOp(SparseLabelOp):
     A ``FermionicOp`` is initialized with a dictionary, mapping terms to their respective
     coefficients:
 
-    .. jupyter-execute::
+    .. code-block:: python
 
         from qiskit_nature.second_q.operators import FermionicOp
 
@@ -62,7 +62,7 @@ class FermionicOp(SparseLabelOp):
     If you have very restricted memory resources available, or would like to avoid the additional
     copy, the dictionary will be stored by reference if you disable ``copy`` like so:
 
-    .. jupyter-execute::
+    .. code-block:: python
 
         some_big_data = {
             "+_0 -_0": 1.0,
@@ -91,40 +91,40 @@ class FermionicOp(SparseLabelOp):
 
     Addition
 
-    .. jupyter-execute::
+    .. code-block:: python
 
       FermionicOp({"+_1": 1}, num_spin_orbitals=2) + FermionicOp({"+_0": 1}, num_spin_orbitals=2)
 
     Sum
 
-    .. jupyter-execute::
+    .. code-block:: python
 
       sum(FermionicOp({label: 1}, num_spin_orbitals=3) for label in ["+_0", "-_1", "+_2 -_2"])
 
     Scalar multiplication
 
-    .. jupyter-execute::
+    .. code-block:: python
 
       0.5 * FermionicOp({"+_1": 1}, num_spin_orbitals=2)
 
     Operator multiplication
 
-    .. jupyter-execute::
+    .. code-block:: python
 
       op1 = FermionicOp({"+_0 -_1": 1}, num_spin_orbitals=2)
       op2 = FermionicOp({"-_0 +_0 +_1": 1}, num_spin_orbitals=2)
       print(op1 @ op2)
 
     Tensor multiplication
 
-    .. jupyter-execute::
+    .. code-block:: python
 
       op = FermionicOp({"+_0 -_1": 1}, num_spin_orbitals=2)
       print(op ^ op)
 
     Adjoint
 
-    .. jupyter-execute::
+    .. code-block:: python
 
       FermionicOp({"+_0 -_1": 1j}, num_spin_orbitals=2).adjoint()
 
@@ -258,6 +258,7 @@ def from_polynomial_tensor(cls, tensor: PolynomialTensor) -> FermionicOp:
                 data[""] = cast(float, tensor[key])
                 continue
 
+            # TODO: extract label_template into Tensor class
             label_template = " ".join(f"{op}_{{}}" for op in key)
 
             # PERF: the following matrix unpacking is a performance bottleneck!