Merge 7a23aee into e6de198

src/fermilib/ops/_fermion_operator.py

@@ -263,6 +263,29 @@ def __init__(self, term=None, coefficient=1.):
                         'Invalid action in FermionOperator: '
                         'Must be 0 (lowering) or 1 (raising).')
 
+    @staticmethod
+    def annihilate_create(annihilate_mode,
+                          create_mode,
+                          coefficient=1.0,
+                          anti_hermitian=False):
+        """
+        Args:
+            annihilate_mode (int)
+            create_mode (int)
+            coefficient (float|complex)
+            anti_hermitian (bool)
+        Returns:
+            FermionOperator: A paired annihilation-creation operator.
+        """
+        result = FermionOperator(term=((annihilate_mode, 1),
+                                       (create_mode, 0)),
+                                 coefficient=coefficient)
+        if anti_hermitian:
+            result += FermionOperator(term=((create_mode, 1),
+                                            (annihilate_mode, 0)),
+                                      coefficient=-coefficient)
+        return result
+
     def compress(self, abs_tol=EQ_TOLERANCE):
         """
         Eliminates all terms with coefficients close to zero and removes

src/fermilib/ops/_fermion_operator_test.py

@@ -427,7 +427,7 @@ def test_add(self):
     def test_add_bad_addend(self):
         op = FermionOperator((), 1.0)
         with self.assertRaises(TypeError):
-            op = op + "0.5"
+            _ = op + "0.5"
 
     def test_sub(self):
         term_a = ((1, 1), (3, 1), (8, 1))
@@ -446,7 +446,29 @@ def test_sub(self):
     def test_sub_bad_subtrahend(self):
         op = FermionOperator((), 1.0)
         with self.assertRaises(TypeError):
-            op = op - "0.5"
+            _ = op - "0.5"
+
+    def test_annihilate_create(self):
+        at = FermionOperator(term=((100, 1),))
+        a = FermionOperator(term=((100, 0),))
+
+        bt = FermionOperator(term=((101, 1),))
+        b = FermionOperator(term=((101, 0),))
+
+        ct = FermionOperator(term=((102, 1),))
+        c = FermionOperator(term=((102, 0),))
+
+        self.assertTrue(FermionOperator.annihilate_create(
+            100, 101).isclose(at * b))
+        self.assertTrue(FermionOperator.annihilate_create(
+            100, 101, 0.5j).isclose(at * b * 0.5j))
+        self.assertTrue(FermionOperator.annihilate_create(
+            100, 102, 0.2j).isclose(at * c * 0.2j))
+        self.assertTrue(FermionOperator.annihilate_create(
+            100, 102, anti_hermitian=True).isclose(at * c - ct * a))
+        self.assertTrue(FermionOperator.annihilate_create(
+            100, 101, 0.2j, anti_hermitian=True).isclose(
+                at * b * 0.2j - bt * a * 0.2j))
 
     def test_isub_different_term(self):
         term_a = ((1, 1), (3, 1), (8, 0))
@@ -468,7 +490,7 @@ def test_isub_bad_addend(self):
 
     def test_neg(self):
         op = FermionOperator(((1, 1), (3, 1), (8, 1)), 0.5)
-        -op
+        _ = -op
         # out of place
         self.assertTrue(op.isclose(FermionOperator(((1, 1), (3, 1), (8, 1)),
                                                    0.5)))

src/fermilib/transforms/_bravyi_kitaev_test.py

@@ -179,7 +179,7 @@ def test_bk_jw_integration_original(self):
         n_qubits = 5
         fermion_operator = FermionOperator(((3, 1), (2, 1), (1, 0), (0, 0)),
                                            -4.3)
-        fermion_operator += FermionOperator(((3, 1), (1, 0)), 8.17)
+        fermion_operator += FermionOperator.annihilate_create(3, 1, 8.17)
         fermion_operator += 3.2 * FermionOperator()
 
         # Map to qubits and compare matrix versions.

src/fermilib/transforms/_conversion.py

@@ -12,24 +12,22 @@
 
 """Transformations acting on operators and RDMs."""
 from __future__ import absolute_import
-from future.utils import iteritems
 
-import copy
 import itertools
+
 import numpy
+from future.utils import iteritems
+from projectq.ops import QubitOperator
 
 from fermilib.ops import (FermionOperator,
                           normal_ordered,
-                          number_operator,
                           InteractionOperator,
                           InteractionRDM)
 from fermilib.ops._interaction_operator import InteractionOperatorError
 from fermilib.utils import (count_qubits,
                             jordan_wigner_sparse,
                             qubit_operator_sparse)
 
-from projectq.ops import QubitOperator
-
 
 def get_sparse_operator(operator, n_qubits=None):
     """Map a Fermion, Qubit, or InteractionOperator to a SparseOperator."""
@@ -161,17 +159,15 @@ def get_fermion_operator(interaction_operator):
     n_qubits = count_qubits(interaction_operator)
 
     # Add one-body terms.
-    for p in range(n_qubits):
-        for q in range(n_qubits):
-            coefficient = interaction_operator[p, q]
-            fermion_operator += FermionOperator(((p, 1), (q, 0)), coefficient)
-
-            # Add two-body terms.
-            for r in range(n_qubits):
-                for s in range(n_qubits):
-                    coefficient = interaction_operator[p, q, r, s]
-                    fermion_operator += FermionOperator(((p, 1), (q, 1),
-                                                         (r, 0), (s, 0)),
-                                                        coefficient)
+    for p, q in itertools.product(range(n_qubits), repeat=2):
+        fermion_operator += FermionOperator.annihilate_create(
+            p, q, coefficient=interaction_operator[p, q])
+
+    # Add two-body terms.
+    for p, q, r, s in itertools.product(range(n_qubits), repeat=4):
+        coefficient = interaction_operator[p, q, r, s]
+        fermion_operator += FermionOperator(((p, 1), (q, 1),
+                                             (r, 0), (s, 0)),
+                                            coefficient)
 
     return fermion_operator

src/fermilib/utils/_sparse_tools_test.py

@@ -36,7 +36,7 @@ def test_qubit_jw_fermion_integration(self):
         # Initialize a random fermionic operator.
         fermion_operator = FermionOperator(((3, 1), (2, 1), (1, 0), (0, 0)),
                                            -4.3)
-        fermion_operator += FermionOperator(((3, 1), (1, 0)), 8.17)
+        fermion_operator += FermionOperator.annihilate_create(3, 1, 8.17)
         fermion_operator += 3.2 * FermionOperator()
 
         # Map to qubits and compare matrix versions.

src/fermilib/utils/_unitary_cc.py

@@ -60,11 +60,8 @@ def uccsd_operator(single_amplitudes, double_amplitudes, anti_hermitian=True):
     for i, j in itertools.product(range(n_orbitals), repeat=2):
         if single_amplitudes[i, j] == 0.:
             continue
-        uccsd_generator += FermionOperator(
-            ((i, 1), (j, 0)), single_amplitudes[i, j])
-        if (anti_hermitian):
-            uccsd_generator += FermionOperator(
-                ((j, 1), (i, 0)), -single_amplitudes[i, j])
+        uccsd_generator += FermionOperator.annihilate_create(
+            i, j, single_amplitudes[i, j], anti_hermitian)
 
         # Add double excitations
     for i, j, k, l in itertools.product(range(n_orbitals), repeat=4):
@@ -77,7 +74,6 @@ def uccsd_operator(single_amplitudes, double_amplitudes, anti_hermitian=True):
             uccsd_generator += FermionOperator(
                 ((l, 1), (k, 0), (j, 1), (i, 0)),
                 -double_amplitudes[i, j, k, l])
-
     return uccsd_generator
 
 
@@ -126,7 +122,6 @@ def uccsd_singlet_operator(packed_amplitudes,
     n_occupied = int(numpy.ceil(n_electrons / 2.))
     n_virtual = int(n_qubits / 2 - n_occupied)  # Virtual Spatial Orbitals
     n_t1 = int(n_occupied * n_virtual)
-    n_t2 = int(n_t1 ** 2)
 
     t1 = packed_amplitudes[:n_t1]
     t2 = packed_amplitudes[n_t1:]
@@ -142,38 +137,30 @@ def t2_ind(i, j, k, l):
 
     uccsd_generator = FermionOperator()
 
-    for i in range(n_virtual):
-        for j in range(n_occupied):
-            for s1 in range(2):
-                uccsd_generator += FermionOperator((
-                    (2 * (i + n_occupied) + s1, 1),
-                    (2 * j + s1, 0)),
-                    t1[t1_ind(i, j)])
-
-                uccsd_generator += FermionOperator((
-                    (2 * j + s1, 1),
-                    (2 * (i + n_occupied) + s1, 0)),
-                    -t1[t1_ind(i, j)])
-
-    for i in range(n_virtual):
-        for j in range(n_occupied):
-            for s1 in range(2):
-                for k in range(n_virtual):
-                    for l in range(n_occupied):
-                        for s2 in range(2):
-                            uccsd_generator += FermionOperator((
-                                (2 * (i + n_occupied) + s1, 1),
-                                (2 * j + s1, 0),
-                                (2 * (k + n_occupied) + s2, 1),
-                                (2 * l + s2, 0)),
-                                t2[t2_ind(i, j, k, l)])
-
-                            uccsd_generator += FermionOperator((
-                                (2 * l + s2, 1),
-                                (2 * (k + n_occupied) + s2, 0),
-                                (2 * j + s1, 1),
-                                (2 * (i + n_occupied) + s1, 0)),
-                                -t2[t2_ind(i, j, k, l)])
+    spaces = range(n_virtual), range(n_occupied), range(2)
+
+    for i, j, s in itertools.product(*spaces):
+        uccsd_generator += FermionOperator.annihilate_create(
+            2 * (i + n_occupied) + s,
+            2 * j + s,
+            coefficient=t1[t1_ind(i, j)],
+            anti_hermitian=True)
+
+    for i, j, s, i2, j2, s2 in itertools.product(*spaces, repeat=2):
+        uccsd_generator += FermionOperator((
+            (2 * (i + n_occupied) + s, 1),
+            (2 * j + s, 0),
+            (2 * (i2 + n_occupied) + s2, 1),
+            (2 * j2 + s2, 0)),
+            t2[t2_ind(i, j, i2, j2)])
+
+        uccsd_generator += FermionOperator((
+            (2 * j2 + s2, 1),
+            (2 * (i2 + n_occupied) + s2, 0),
+            (2 * j + s, 1),
+            (2 * (i + n_occupied) + s, 0)),
+            -t2[t2_ind(i, j, i2, j2)])
+
     return uccsd_generator