Repository: quantumlib/OpenFermion

Issue:
Implement Bravyi-Kitaev Transform for MajoranaOperator

Discussion:
kevinsung: Another coveralls false positive.

Patch:
diff --git a/src/openfermion/transforms/_bravyi_kitaev.py b/src/openfermion/transforms/_bravyi_kitaev.py
index 6007dcf..f999465 100644
--- a/src/openfermion/transforms/_bravyi_kitaev.py
+++ b/src/openfermion/transforms/_bravyi_kitaev.py
@@ -12,7 +12,7 @@
 
 """Bravyi-Kitaev transform on fermionic operators."""
 
-from openfermion.ops import QubitOperator
+from openfermion.ops import FermionOperator, MajoranaOperator, QubitOperator
 from openfermion.utils import count_qubits, inline_sum
 
 
@@ -40,20 +40,12 @@ def bravyi_kitaev(operator, n_qubits=None):
     Raises:
         ValueError: Invalid number of qubits specified.
     """
-    # Compute the number of qubits.
-    if n_qubits is None:
-        n_qubits = count_qubits(operator)
-    if n_qubits < count_qubits(operator):
-        raise ValueError('Invalid number of qubits specified.')
-
-    # Compute transformed operator.
-    transformed_terms = (
-        _transform_operator_term(term=term,
-                                 coefficient=operator.terms[term],
-                                 n_qubits=n_qubits)
-        for term in operator.terms
-    )
-    return inline_sum(summands=transformed_terms, seed=QubitOperator())
+    if isinstance(operator, FermionOperator):
+        return _bravyi_kitaev_fermion_operator(operator, n_qubits)
+    if isinstance(operator, MajoranaOperator):
+        return _bravyi_kitaev_majorana_operator(operator, n_qubits)
+    raise TypeError("Couldn't apply the Bravyi-Kitaev Transform to object "
+                    "of type {}.".format(type(operator)))
 
 
 def _update_set(index, n_qubits):
@@ -105,6 +97,52 @@ def _parity_set(index):
     return indices
 
 
+def _bravyi_kitaev_majorana_operator(operator, n_qubits):
+    # Compute the number of qubits.
+    N = count_qubits(operator)
+    if n_qubits is None:
+        n_qubits = N
+    if n_qubits < N:
+        raise ValueError('Invalid number of qubits specified.')
+
+    # Compute transformed operator.
+    transformed_terms = (
+        _transform_majorana_term(term=term,
+                                 coefficient=coeff,
+                                 n_qubits=n_qubits)
+        for term, coeff in operator.terms.items()
+    )
+    return inline_sum(summands=transformed_terms, seed=QubitOperator())
+
+
+def _transform_majorana_term(term, coefficient, n_qubits):
+    # Build the Bravyi-Kitaev transformed operators.
+    transformed_ops = (
+        _transform_majorana_operator(majorana_index, n_qubits)
+        for majorana_index in term
+    )
+    return inline_product(factors=transformed_ops,
+                          seed=QubitOperator((), coefficient))
+
+
+def _transform_majorana_operator(majorana_index, n_qubits):
+    q, b = divmod(majorana_index, 2)
+
+    update_set = _update_set(q, n_qubits)
+    occupation_set = _occupation_set(q)
+    parity_set = _parity_set(q - 1)
+
+    if b:
+        return QubitOperator(
+                [(q, 'Y')] +
+                [(i, 'X') for i in update_set - {q}] +
+                [(i, 'Z') for i in (parity_set ^ occupation_set) - {q}])
+    else:
+        return QubitOperator(
+                [(i, 'X') for i in update_set] +
+                [(i, 'Z') for i in parity_set])
+
+
 def _transform_operator_term(term, coefficient, n_qubits):
     """
     Args:
@@ -123,7 +161,24 @@ def _transform_operator_term(term, coefficient, n_qubits):
     )
     return inline_product(factors=transformed_ladder_ops,
                           seed=QubitOperator((), coefficient))
-                          
+
+
+def _bravyi_kitaev_fermion_operator(operator, n_qubits):
+    # Compute the number of qubits.
+    N = count_qubits(operator)
+    if n_qubits is None:
+        n_qubits = N
+    if n_qubits < N:
+        raise ValueError('Invalid number of qubits specified.')
+
+    # Compute transformed operator.
+    transformed_terms = (
+        _transform_operator_term(term=term,
+                                 coefficient=operator.terms[term],
+                                 n_qubits=n_qubits)
+        for term in operator.terms
+    )
+    return inline_sum(summands=transformed_terms, seed=QubitOperator())
 
 
 def _transform_ladder_operator(ladder_operator, n_qubits):
diff --git a/src/openfermion/utils/_operator_utils.py b/src/openfermion/utils/_operator_utils.py
index dc7af40..f36bdf5 100644
--- a/src/openfermion/utils/_operator_utils.py
+++ b/src/openfermion/utils/_operator_utils.py
@@ -11,7 +11,7 @@
 #   limitations under the License.
 
 """This module provides generic tools for classes in ops/"""
-from __future__ import absolute_import
+from __future__ import absolute_import, division
 from builtins import map, zip
 import copy
 import itertools
@@ -28,6 +28,7 @@ from openfermion.ops import (BosonOperator,
                              FermionOperator,
                              InteractionOperator,
                              InteractionRDM,
+                             MajoranaOperator,
                              QuadOperator,
                              QubitOperator,
                              PolynomialTensor)
@@ -326,6 +327,15 @@ def count_qubits(operator):
                     num_qubits = term[-1][0] + 1
         return num_qubits
 
+    # Handle MajoranaOperator.
+    if isinstance(operator, MajoranaOperator):
+        num_qubits = 0
+        for term in operator.terms:
+            for majorana_index in term:
+                if numpy.ceil((majorana_index+1) / 2) > num_qubits:
+                    num_qubits = int(numpy.ceil((majorana_index+1) / 2))
+        return num_qubits
+
     # Handle DiagonalCoulombHamiltonian
     elif isinstance(operator, DiagonalCoulombHamiltonian):
         return operator.one_body.shape[0]


Test patch:
diff --git a/src/openfermion/transforms/_bravyi_kitaev_test.py b/src/openfermion/transforms/_bravyi_kitaev_test.py
index 33db57d..834c5b4 100644
--- a/src/openfermion/transforms/_bravyi_kitaev_test.py
+++ b/src/openfermion/transforms/_bravyi_kitaev_test.py
@@ -18,8 +18,10 @@ import unittest
 import numpy
 
 from openfermion.ops import (FermionOperator,
+                             MajoranaOperator,
                              QubitOperator)
-from openfermion.transforms import bravyi_kitaev, jordan_wigner
+from openfermion.transforms import (
+        bravyi_kitaev, get_fermion_operator, jordan_wigner)
 from openfermion.utils import eigenspectrum, number_operator
 
 
@@ -68,6 +70,8 @@ class BravyiKitaevTransformTest(unittest.TestCase):
     def test_bk_n_qubits_too_small(self):
         with self.assertRaises(ValueError):
             bravyi_kitaev(FermionOperator('2^ 3^ 5 0'), n_qubits=4)
+        with self.assertRaises(ValueError):
+            bravyi_kitaev(MajoranaOperator((2, 3, 9, 0)), n_qubits=4)
 
     def test_bk_jw_number_operator(self):
         # Check if number operator has the same spectrum in both
@@ -186,3 +190,13 @@ class BravyiKitaevTransformTest(unittest.TestCase):
         bk_spectrum = eigenspectrum(bk_qubit_operator)
         self.assertAlmostEqual(0., numpy.amax(numpy.absolute(jw_spectrum -
                                               bk_spectrum)), places=5)
+    def test_bk_bad_type(self):
+        with self.assertRaises(TypeError):
+            bravyi_kitaev(QubitOperator())
+
+
+def test_bravyi_kitaev_majorana_op_consistent():
+    op = (MajoranaOperator((1, 3, 4), 0.5)
+          + MajoranaOperator((3, 7, 8, 9, 10, 12), 1.8)
+          + MajoranaOperator((0, 4)))
+    assert bravyi_kitaev(op) == bravyi_kitaev(get_fermion_operator(op))
diff --git a/src/openfermion/utils/_operator_utils_test.py b/src/openfermion/utils/_operator_utils_test.py
index 612467b..ae4a2d4 100644
--- a/src/openfermion/utils/_operator_utils_test.py
+++ b/src/openfermion/utils/_operator_utils_test.py
@@ -38,6 +38,7 @@ class OperatorUtilsTest(unittest.TestCase):
 
     def setUp(self):
         self.n_qubits = 5
+        self.majorana_operator = MajoranaOperator((1, 4, 9))
         self.fermion_term = FermionOperator('1^ 2^ 3 4', -3.17)
         self.fermion_operator = self.fermion_term + hermitian_conjugated(
             self.fermion_term)
@@ -45,6 +46,10 @@ class OperatorUtilsTest(unittest.TestCase):
         self.interaction_operator = get_interaction_operator(
             self.fermion_operator)
 
+    def test_n_qubits_majorana_operator(self):
+        self.assertEqual(self.n_qubits,
+                         count_qubits(self.majorana_operator))
+
     def test_n_qubits_single_fermion_term(self):
         self.assertEqual(self.n_qubits,
                          count_qubits(self.fermion_term))