Repository: quantumlib/OpenFermion
Pull request #90 (2017-10-28 20:04:55)

Incorporate constant into n_body_tensors dictionary of PolynomialTensor
Now `n_body_tensors[()]` is a complex number, the constant term.


Patch:
diff --git a/src/openfermion/ops/_interaction_operator.py b/src/openfermion/ops/_interaction_operator.py
index 145fe3a..9d08f18 100644
--- a/src/openfermion/ops/_interaction_operator.py
+++ b/src/openfermion/ops/_interaction_operator.py
@@ -61,8 +61,10 @@ class InteractionOperator(PolynomialTensor):
         """
         # Make sure nonzero elements are only for normal ordered terms.
         super(InteractionOperator, self).__init__(
-                constant,
-                {(1, 0): one_body_tensor, (1, 1, 0, 0): two_body_tensor})
+                {(): constant,
+                 (1, 0): one_body_tensor,
+                 (1, 1, 0, 0): two_body_tensor})
+        self.constant = self.n_body_tensors[()]
         self.one_body_tensor = self.n_body_tensors[1, 0]
         self.two_body_tensor = self.n_body_tensors[1, 1, 0, 0]
 
diff --git a/src/openfermion/ops/_interaction_rdm.py b/src/openfermion/ops/_interaction_rdm.py
index d602aa4..bf0adf3 100644
--- a/src/openfermion/ops/_interaction_rdm.py
+++ b/src/openfermion/ops/_interaction_rdm.py
@@ -45,7 +45,6 @@ class InteractionRDM(PolynomialTensor):
                 <a^\dagger_p a^\dagger_q a_r a_s>.
         """
         super(InteractionRDM, self).__init__(
-                None,
                 {(1, 0): one_body_tensor, (1, 1, 0, 0): two_body_tensor})
         self.one_body_tensor = self.n_body_tensors[1, 0]
         self.two_body_tensor = self.n_body_tensors[1, 1, 0, 0]
diff --git a/src/openfermion/ops/_polynomial_tensor.py b/src/openfermion/ops/_polynomial_tensor.py
index 489095d..ffdac71 100644
--- a/src/openfermion/ops/_polynomial_tensor.py
+++ b/src/openfermion/ops/_polynomial_tensor.py
@@ -100,8 +100,6 @@ class PolynomialTensor(object):
 
     Attributes:
         n_qubits(int): The number of sites on which the tensor acts.
-        constant(complex or float): A constant term in the operator given
-            as a complex number. For instance, the nuclear repulsion energy.
         n_body_tensors(dict): A dictionary storing the tensors describing
             n-body interactions. The keys are tuples that indicate the
             type of tensor. For instance, n_body_tensors[(1, 0)] would
@@ -112,21 +110,21 @@ class PolynomialTensor(object):
             of the form a_i a^\dagger_j.
     """
 
-    def __init__(self, constant, n_body_tensors):
+    def __init__(self, n_body_tensors):
         """Initialize the PolynomialTensor class.
 
         Args:
-            constant(complex or float): A constant term in the operator given
-                as a complex number. For instance, the nuclear repulsion
-                energy.
             n_body_tensors(dict): A dictionary storing the tensors describing
                 n-body interactions.
         """
-        if constant is None:
-            constant = 0.
-        self.constant = constant
         self.n_body_tensors = n_body_tensors
-        self.n_qubits = list(n_body_tensors.values())[0].shape[0]
+
+        # Set n_qubits
+        key_iterator = iter(n_body_tensors.keys())
+        key = next(key_iterator)
+        if key == ():
+            key = next(key_iterator)
+        self.n_qubits = n_body_tensors[key].shape[0]
 
     def __getitem__(self, args):
         """Look up matrix element.
@@ -138,7 +136,7 @@ class PolynomialTensor(object):
                 `my_tensor.n_body_tensors[1, 1, 0][6, 8, 2]`
         """
         if len(args) == 0:
-            return self.constant
+            return self.n_body_tensors[()]
         else:
             index = tuple([operator[0] for operator in args])
             key = tuple([operator[1] for operator in args])
@@ -151,7 +149,7 @@ class PolynomialTensor(object):
             args: Tuples indicating which coefficient to set.
         """
         if len(args) == 0:
-            self.constant = value
+            self.n_body_tensors[()] = value
         else:
             key = tuple([operator[1] for operator in args])
             index = tuple([operator[0] for operator in args])
@@ -162,7 +160,7 @@ class PolynomialTensor(object):
             return False
         if self.n_body_tensors.keys() != other_operator.n_body_tensors.keys():
             return False
-        diff = abs(other_operator.constant - self.constant)
+        diff = 0.
         for key in self.n_body_tensors:
             self_tensor = self.n_body_tensors[key]
             other_tensor = other_operator.n_body_tensors[key]
@@ -184,7 +182,6 @@ class PolynomialTensor(object):
         if self.n_body_tensors.keys() != addend.n_body_tensors.keys():
             raise TypeError('Invalid tensor type.')
 
-        self.constant += addend.constant
         for key in self.n_body_tensors:
             self.n_body_tensors[key] = numpy.add(self.n_body_tensors[key],
                                                  addend.n_body_tensors[key])
@@ -199,8 +196,7 @@ class PolynomialTensor(object):
         neg_n_body_tensors = dict()
         for key in self.n_body_tensors:
             neg_n_body_tensors[key] = numpy.negative(self.n_body_tensors[key])
-        return PolynomialTensor(-self.constant,
-                                n_body_tensors=neg_n_body_tensors)
+        return PolynomialTensor(neg_n_body_tensors)
 
     def __isub__(self, subtrahend):
         if not issubclass(type(subtrahend), PolynomialTensor):
@@ -212,7 +208,6 @@ class PolynomialTensor(object):
         if self.n_body_tensors.keys() != subtrahend.n_body_tensors.keys():
             raise TypeError('Invalid tensor type.')
 
-        self.constant -= subtrahend.constant
         for key in self.n_body_tensors:
             self.n_body_tensors[key] = numpy.subtract(
                     self.n_body_tensors[key], subtrahend.n_body_tensors[key])
@@ -233,7 +228,6 @@ class PolynomialTensor(object):
         if self.n_body_tensors.keys() != multiplier.n_body_tensors.keys():
             raise TypeError('Invalid tensor type.')
 
-        self.constant *= multiplier.constant
         for key in self.n_body_tensors:
             self.n_body_tensors[key] = numpy.multiply(
                     self.n_body_tensors[key], multiplier.n_body_tensors[key])
@@ -246,33 +240,32 @@ class PolynomialTensor(object):
 
     def __iter__(self):
         """Iterate over non-zero elements of PolynomialTensor."""
-        # Constant.
-        if self.constant:
-            yield ()
-
-        # n-body elements
         def sort_key(key):
             """This determines how the keys to n_body_tensors
             should be sorted."""
             # Interpret key as an integer written in binary
-            key_int = int(''.join(map(str, key)))
-            return len(key), key_int
+            if key == ():
+                return 0, 0
+            else:
+                key_int = int(''.join(map(str, key)))
+                return len(key), key_int
+
         for key in sorted(self.n_body_tensors.keys(), key=sort_key):
-            n_body_tensor = self.n_body_tensors[key]
-            for index in itertools.product(
-                    range(self.n_qubits), repeat=len(key)):
-                if n_body_tensor[index]:
-                    yield tuple(zip(index, key))
+            if key == ():
+                yield ()
+            else:
+                n_body_tensor = self.n_body_tensors[key]
+                for index in itertools.product(
+                        range(self.n_qubits), repeat=len(key)):
+                    if n_body_tensor[index]:
+                        yield tuple(zip(index, key))
 
     def __str__(self):
         """Print out the non-zero elements of PolynomialTensor."""
-        string = ''
+        strings = []
         for key in self:
-            if len(key) == 0:
-                string += '() {}\n'.format(self[key])
-            else:
-                string += '{} {}\n'.format(key, self[key])
-        return string if string else '0'
+            strings.append('{} {}\n'.format(key, self[key]))
+        return ''.join(strings) if strings else '0'
 
     def rotate_basis(self, rotation_matrix):
         """
@@ -283,10 +276,12 @@ class PolynomialTensor(object):
                 dimensions of n_qubits by n_qubits. Assumed to be real and
                 invertible.
         """
-        self.n_body_tensors[1, 0] = one_body_basis_change(
-            self.n_body_tensors[1, 0], rotation_matrix)
-        self.n_body_tensors[1, 1, 0, 0] = two_body_basis_change(
-            self.n_body_tensors[1, 1, 0, 0], rotation_matrix)
+        if (1, 0) in self.n_body_tensors:
+            self.n_body_tensors[1, 0] = one_body_basis_change(
+                self.n_body_tensors[1, 0], rotation_matrix)
+        if (1, 1, 0, 0) in self.n_body_tensors:
+            self.n_body_tensors[1, 1, 0, 0] = two_body_basis_change(
+                self.n_body_tensors[1, 1, 0, 0], rotation_matrix)
 
     def __repr__(self):
         return str(self)
diff --git a/src/openfermion/ops/_quadratic_hamiltonian.py b/src/openfermion/ops/_quadratic_hamiltonian.py
index 285e9ab..0254b30 100644
--- a/src/openfermion/ops/_quadratic_hamiltonian.py
+++ b/src/openfermion/ops/_quadratic_hamiltonian.py
@@ -76,8 +76,8 @@ class QuadraticHamiltonian(PolynomialTensor):
             antisymmetric_part = numpy.zeros((n_qubits, n_qubits), complex)
 
         super(QuadraticHamiltonian, self).__init__(
-                constant,
-                {(1, 0): combined_hermitian_part,
+                {(): constant,
+                 (1, 0): combined_hermitian_part,
                  (1, 1): .5 * antisymmetric_part,
                  (0, 0): -.5 * antisymmetric_part.conj()})
         self.chemical_potential = chemical_potential


Test patch:
diff --git a/src/openfermion/ops/_polynomial_tensor_test.py b/src/openfermion/ops/_polynomial_tensor_test.py
index fbc6a3d..38fac0c 100644
--- a/src/openfermion/ops/_polynomial_tensor_test.py
+++ b/src/openfermion/ops/_polynomial_tensor_test.py
@@ -36,8 +36,9 @@ class PolynomialTensorTest(unittest.TestCase):
         two_body_a[1, 1, 0, 0] = 5
 
         self.polynomial_tensor_a = PolynomialTensor(
-                self.constant,
-                {(1, 0): one_body_a, (1, 1, 0, 0): two_body_a})
+                {(): self.constant,
+                 (1, 0): one_body_a,
+                 (1, 1, 0, 0): two_body_a})
 
         one_body_na = numpy.zeros((self.n_qubits, self.n_qubits))
         two_body_na = numpy.zeros((self.n_qubits, self.n_qubits,
@@ -47,8 +48,9 @@ class PolynomialTensorTest(unittest.TestCase):
         two_body_na[0, 1, 0, 1] = -4
         two_body_na[1, 1, 0, 0] = -5
         self.polynomial_tensor_na = PolynomialTensor(
-                -self.constant,
-                {(1, 0): one_body_na, (1, 1, 0, 0): two_body_na})
+                {(): -self.constant,
+                 (1, 0): one_body_na,
+                 (1, 1, 0, 0): two_body_na})
 
         one_body_b = numpy.zeros((self.n_qubits, self.n_qubits))
         two_body_b = numpy.zeros((self.n_qubits, self.n_qubits,
@@ -58,8 +60,9 @@ class PolynomialTensorTest(unittest.TestCase):
         two_body_b[0, 1, 0, 1] = 3
         two_body_b[1, 0, 0, 1] = 4
         self.polynomial_tensor_b = PolynomialTensor(
-                self.constant,
-                {(1, 0): one_body_b, (1, 1, 0, 0): two_body_b})
+                {(): self.constant,
+                 (1, 0): one_body_b,
+                 (1, 1, 0, 0): two_body_b})
 
         one_body_ab = numpy.zeros((self.n_qubits, self.n_qubits))
         two_body_ab = numpy.zeros((self.n_qubits, self.n_qubits,
@@ -70,8 +73,9 @@ class PolynomialTensorTest(unittest.TestCase):
         two_body_ab[1, 0, 0, 1] = 4
         two_body_ab[1, 1, 0, 0] = 5
         self.polynomial_tensor_ab = PolynomialTensor(
-                2.0 * self.constant,
-                {(1, 0): one_body_ab, (1, 1, 0, 0): two_body_ab})
+                {(): 2.0 * self.constant,
+                 (1, 0): one_body_ab,
+                 (1, 1, 0, 0): two_body_ab})
 
         constant_axb = self.constant * self.constant
         one_body_axb = numpy.zeros((self.n_qubits, self.n_qubits))
@@ -81,8 +85,9 @@ class PolynomialTensorTest(unittest.TestCase):
         one_body_axb[1, 0] = 6
         two_body_axb[0, 1, 0, 1] = 12
         self.polynomial_tensor_axb = PolynomialTensor(
-                constant_axb,
-                {(1, 0): one_body_axb, (1, 1, 0, 0): two_body_axb})
+                {(): constant_axb,
+                 (1, 0): one_body_axb,
+                 (1, 1, 0, 0): two_body_axb})
 
         self.n_qubits_plus_one = self.n_qubits + 1
         one_body_c = numpy.zeros((self.n_qubits_plus_one,
@@ -96,8 +101,9 @@ class PolynomialTensorTest(unittest.TestCase):
         two_body_c[0, 1, 0, 1] = 3
         two_body_c[1, 0, 0, 1] = 4
         self.polynomial_tensor_c = PolynomialTensor(
-                self.constant,
-                {(1, 0): one_body_c, (1, 1, 0, 0): two_body_c})
+                {(): self.constant,
+                 (1, 0): one_body_c,
+                 (1, 1, 0, 0): two_body_c})
 
         one_body_hole = numpy.zeros((self.n_qubits, self.n_qubits))
         two_body_hole = numpy.zeros((self.n_qubits, self.n_qubits,
@@ -108,8 +114,9 @@ class PolynomialTensorTest(unittest.TestCase):
         two_body_hole[1, 1, 0, 0] = 5
 
         self.polynomial_tensor_hole = PolynomialTensor(
-                self.constant,
-                {(0, 1): one_body_hole, (0, 0, 1, 1): two_body_hole})
+                {(): self.constant,
+                 (0, 1): one_body_hole,
+                 (0, 0, 1, 1): two_body_hole})
 
         one_body_spinful = numpy.zeros((2 * self.n_qubits, 2 * self.n_qubits))
         two_body_spinful = numpy.zeros((2 * self.n_qubits, 2 * self.n_qubits,
@@ -124,8 +131,9 @@ class PolynomialTensorTest(unittest.TestCase):
         two_body_spinful[3, 3, 2, 2] = 9
 
         self.polynomial_tensor_spinful = PolynomialTensor(
-                self.constant,
-                {(1, 0): one_body_spinful, (1, 1, 0, 0): two_body_spinful})
+                {(): self.constant,
+                 (1, 0): one_body_spinful,
+                 (1, 1, 0, 0): two_body_spinful})
 
     def test_setitem_1body(self):
         expected_one_body_tensor = numpy.array([[0, 3], [2, 0]])
@@ -254,33 +262,14 @@ class PolynomialTensorTest(unittest.TestCase):
         one_body[0, 1] = 11.0
         two_body[0, 1, 1, 0] = 22.0
         polynomial_tensor = PolynomialTensor(
-                self.constant,
-                {(1, 0): one_body, (1, 1, 0, 0): two_body})
+                {(): self.constant,
+                 (1, 0): one_body,
+                 (1, 1, 0, 0): two_body})
         want_str = ('() 23.0\n((0, 1), (1, 0)) 11.0\n'
                     '((0, 1), (1, 1), (1, 0), (0, 0)) 22.0\n')
         self.assertEqual(polynomial_tensor.__str__(), want_str)
         self.assertEqual(polynomial_tensor.__repr__(), want_str)
 
-    def test_init_none(self):
-        n_qubits = 2
-        constant = 23.0
-
-        one_body = numpy.zeros((self.n_qubits, self.n_qubits))
-        two_body = numpy.zeros((self.n_qubits, self.n_qubits,
-                                self.n_qubits, self.n_qubits))
-        one_body[0, 1] = 2
-        one_body[1, 0] = 3
-        two_body[0, 1, 0, 1] = 4
-        two_body[1, 1, 0, 0] = 5
-        polynomial_tensor = PolynomialTensor(
-                None,
-                {(1, 0): one_body, (1, 1, 0, 0): two_body})
-        self.assertAlmostEqual(polynomial_tensor.constant, 0.)
-
-        self.polynomial_tensor_a[()] = 0.
-        self.assertAlmostEqual(self.polynomial_tensor_a, polynomial_tensor)
-        self.polynomial_tensor_a[()] = self.constant
-
     def test_rotate_basis_identical(self):
         rotation_matrix_identical = numpy.zeros((self.n_qubits, self.n_qubits))
         rotation_matrix_identical[0, 0] = 1
@@ -310,17 +299,21 @@ class PolynomialTensorTest(unittest.TestCase):
                                          s + self.n_qubits] = j
                         j = j + 1
         polynomial_tensor = PolynomialTensor(
-                self.constant,
-                {(1, 0): one_body, (1, 1, 0, 0): two_body})
+                {(): self.constant,
+                 (1, 0): one_body,
+                 (1, 1, 0, 0): two_body})
         want_polynomial_tensor = PolynomialTensor(
-                self.constant,
-                {(1, 0): one_body, (1, 1, 0, 0): two_body})
+                {(): self.constant,
+                 (1, 0): one_body,
+                 (1, 1, 0, 0): two_body})
         polynomial_tensor_spinful = PolynomialTensor(
-                self.constant,
-                {(1, 0): one_body_spinful, (1, 1, 0, 0): two_body_spinful})
+                {(): self.constant,
+                 (1, 0): one_body_spinful,
+                 (1, 1, 0, 0): two_body_spinful})
         want_polynomial_tensor_spinful = PolynomialTensor(
-                self.constant,
-                {(1, 0): one_body_spinful, (1, 1, 0, 0): two_body_spinful})
+                {(): self.constant,
+                 (1, 0): one_body_spinful,
+                 (1, 1, 0, 0): two_body_spinful})
 
         polynomial_tensor.rotate_basis(rotation_matrix_identical)
         polynomial_tensor_spinful.rotate_basis(rotation_matrix_identical)
@@ -356,10 +349,12 @@ class PolynomialTensorTest(unittest.TestCase):
                         two_body_reverse[p, q, r, s] = j_reverse
                         j_reverse = j_reverse - 1
         polynomial_tensor = PolynomialTensor(
-                self.constant,
-                {(1, 0): one_body, (1, 1, 0, 0): two_body})
+                {(): self.constant,
+                 (1, 0): one_body,
+                 (1, 1, 0, 0): two_body})
         want_polynomial_tensor = PolynomialTensor(
-                self.constant,
-                {(1, 0): one_body_reverse, (1, 1, 0, 0): two_body_reverse})
+                {(): self.constant,
+                 (1, 0): one_body_reverse,
+                 (1, 1, 0, 0): two_body_reverse})
         polynomial_tensor.rotate_basis(rotation_matrix_reverse)
         self.assertEqual(polynomial_tensor, want_polynomial_tensor)