From 0a703da8e0db7fbceb588e389043834a1e8c4816 Mon Sep 17 00:00:00 2001
From: chmwzc <70616433+chmwzc@users.noreply.github.com>
Date: Tue, 27 Feb 2024 11:38:21 +0800
Subject: [PATCH 01/11] Test out pytest parametrize

---
 tests/test_expectation_samples.py | 157 ++++++++++++++----------------
 1 file changed, 73 insertions(+), 84 deletions(-)

diff --git a/tests/test_expectation_samples.py b/tests/test_expectation_samples.py
index 57d2aeb..5eed036 100644
--- a/tests/test_expectation_samples.py
+++ b/tests/test_expectation_samples.py
@@ -16,39 +16,43 @@
 )
 
 
-def test_expectation_z0():
-    """Test from_samples functionality of expectation function"""
-    hamiltonian = SymbolicHamiltonian(Z(0))
-    circuit = Circuit(2)
-    circuit.add(gates.X(0))
-    result = expectation(circuit, hamiltonian, from_samples=True)
-    assert pytest.approx(result) == -1.0
-
-
-def test_expectation_z0z1():
-    """Tests expectation_from_samples for diagonal Hamiltonians (only Z's)"""
-    hamiltonian = SymbolicHamiltonian(Z(0) * Z(1))
-    circuit = Circuit(2)
-    circuit.add(gates.X(0))
-    result = expectation(circuit, hamiltonian, from_samples=True)
-    assert pytest.approx(result) == -1.0
-
-
-def test_expectation_x0():
-    """Tests expectation_from_samples for Hamiltonians with X"""
-    hamiltonian = SymbolicHamiltonian(X(0))
+@pytest.mark.parametrize(
+    "terms,gates_to_add,expected,threshold",
+    [
+        (
+            Z(0),
+            [
+                gates.X(0),
+            ],
+            -1.0,
+            None,
+        ),
+        (
+            Z(0) * Z(1),
+            [
+                gates.X(0),
+            ],
+            -1.0,
+            None,
+        ),
+        (
+            X(0),
+            [
+                gates.H(0),
+            ],
+            1.0,
+            None,
+        ),
+        (X(0), [gates.X(0), gates.X(0)], 0.0, 0.05),
+    ],
+)
+def test_expectation_samples(terms, gates_to_add, expected, threshold):
+    """Test from_samples functionality of expectation function with various Hamiltonians"""
+    hamiltonian = SymbolicHamiltonian(terms)
     circuit = Circuit(2)
-    circuit.add(gates.H(0))
+    circuit.add(gates_to_add)
     result = expectation(circuit, hamiltonian, from_samples=True)
-    assert pytest.approx(result) == 1.0
-
-
-def test_expectation_x0_2():
-    """Test 2 of expectation_from_samples for Hamiltonians with X"""
-    hamiltonian = SymbolicHamiltonian(X(0))
-    circuit = Circuit(2)
-    result = expectation(circuit, hamiltonian, from_samples=True, n_shots=10000)
-    assert pytest.approx(result, abs=0.05) == 0.00
+    assert pytest.approx(result, abs=threshold) == expected
 
 
 def test_measurement_basis_rotations_error():
@@ -58,38 +62,26 @@ def test_measurement_basis_rotations_error():
         _ = measurement_basis_rotations(hamiltonian, 2, grouping="test")
 
 
-def test_allocate_shots_uniform():
+@pytest.mark.parametrize(
+    "method,max_shots_per_term,expected",
+    [
+        ("u", None, [100, 100]),  # Control test; i.e. working normally
+        (None, None, [190, 10]),  # Default arguments test
+        (None, 100, [100, 100]),  # max_shots_per_term error
+        (None, 25, [100, 100]),  # If max_shots_per_term is too small
+        (None, 1000, [190, 10]),  # If max_shots_per_term is too large
+    ],
+)
+def test_allocate_shots(method, max_shots_per_term, expected):
     hamiltonian = SymbolicHamiltonian(94 * Z(0) + Z(1) + 5 * X(0))
     grouped_terms = measurement_basis_rotations(hamiltonian, 1)
     n_shots = 200
-    # Control test; i.e. working normally
-    assert allocate_shots(grouped_terms, method="u", n_shots=n_shots) == [100, 100]
+    assert (
+        allocate_shots(grouped_terms, method=method, n_shots=n_shots, max_shots_per_term=max_shots_per_term) == expected
+    )
 
 
-def test_allocate_shots_coefficient():
-    hamiltonian = SymbolicHamiltonian(94 * Z(0) + Z(1) + 5 * X(0))
-    grouped_terms = measurement_basis_rotations(hamiltonian, 1)
-    n_shots = 200
-    # Default arguments
-    assert allocate_shots(grouped_terms, n_shots=n_shots) == [190, 10], "Default arguments error!"
-    # Reasonable max_shots_per_term test
-    assert allocate_shots(grouped_terms, n_shots=n_shots, max_shots_per_term=100) == [
-        100,
-        100,
-    ], "max_shots_per_term error!"
-    # Too small max_shots_per_term test
-    assert allocate_shots(grouped_terms, n_shots=n_shots, max_shots_per_term=25) == [
-        100,
-        100,
-    ], "max_shots_per_term error: Too small test"
-    # Too big max_shots_per_term test
-    assert allocate_shots(grouped_terms, n_shots=n_shots, max_shots_per_term=1000) == [
-        190,
-        10,
-    ], "max_shots_per_term error: Too big test"
-
-
-def test_allocate_shots_coefficient_edges():
+def test_allocate_shots_coefficient_edge_case():
     """Edge cases of allocate_shots"""
     hamiltonian = SymbolicHamiltonian(Z(0) + X(0))
     grouped_terms = measurement_basis_rotations(hamiltonian, 1)
@@ -104,29 +96,21 @@ def test_allocate_shots_input_validity():
         _ = allocate_shots(grouped_terms, n_shots=1, method="wrong")
 
 
-def test_expectation_manual_shot_allocation():
-    # State vector: -|1>
-    circuit = Circuit(1)
-    circuit.add(gates.X(0))
-    circuit.add(gates.Z(0))
-    hamiltonian = SymbolicHamiltonian(Z(0) + X(0))
-    shot_allocation = (10, 0)
-    result = expectation(
-        circuit, hamiltonian, from_samples=True, n_shots_per_pauli_term=False, shot_allocation=shot_allocation
-    )
-    assert pytest.approx(result) == -1.0
-
-
-def test_expectation_manual_shot_allocation2():
-    # State vector: |1>
+@pytest.mark.parametrize(
+    "gates_to_add,shot_allocation,threshold,expected",
+    [
+        ([gates.X(0), gates.Z(0)], [10, 0], None, -1.0),  # State vector: -|1>, Measuring X
+        ([gates.X(0)], [0, 1000], 0.1, 0.0),  # State vector: |1>, Measuring X
+    ],
+)
+def test_expectation_manual_shot_allocation(gates_to_add, shot_allocation, threshold, expected):
     circuit = Circuit(1)
-    circuit.add(gates.X(0))
+    circuit.add(gates_to_add)
     hamiltonian = SymbolicHamiltonian(Z(0) + X(0))
-    shot_allocation = (0, 1000)
     result = expectation(
         circuit, hamiltonian, from_samples=True, n_shots_per_pauli_term=False, shot_allocation=shot_allocation
     )
-    assert pytest.approx(result, abs=0.1) == 0.0
+    assert pytest.approx(result, abs=threshold) == expected
 
 
 def test_expectation_invalid_shot_allocation():
@@ -139,7 +123,14 @@ def test_expectation_invalid_shot_allocation():
         )
 
 
-def test_h2_hf_energy():
+@pytest.mark.parametrize(
+    "n_shots_per_pauli_term,threshold",
+    [
+        (True, 0.005),  # 10000 shots used for each term in Hamiltonian
+        (False, 0.015),  # 10000 shots divided between each Pauli string in Hamiltonian
+    ],
+)
+def test_h2_hf_energy(n_shots_per_pauli_term, threshold):
     """Test HF energy of H2 molecule"""
     # Hardcoded benchmark results
     h2_ref_energy = -1.117349035
@@ -153,10 +144,8 @@ def test_h2_hf_energy():
     # Molecular Hamiltonian and the HF expectation value
     hamiltonian = h2.hamiltonian()
 
-    # n_shots (=10000) allocated to each term
-    hf_energy = expectation(circuit, hamiltonian, from_samples=True, n_shots=10000)
-    assert pytest.approx(hf_energy, abs=0.005) == h2_ref_energy
-
-    # n_shots divided amongst every term
-    hf_energy = expectation(circuit, hamiltonian, from_samples=True, n_shots_per_pauli_term=False, n_shots=10000)
-    assert pytest.approx(hf_energy, abs=0.01) == h2_ref_energy
+    n_shots = 10000
+    hf_energy = expectation(
+        circuit, hamiltonian, from_samples=True, n_shots_per_pauli_term=n_shots_per_pauli_term, n_shots=n_shots
+    )
+    assert pytest.approx(hf_energy, abs=threshold) == h2_ref_energy

From aa64bea8d0d6f7d7ecd7b4803fc0255ef178116e Mon Sep 17 00:00:00 2001
From: chmwzc <70616433+chmwzc@users.noreply.github.com>
Date: Tue, 27 Feb 2024 15:02:53 +0800
Subject: [PATCH 02/11] Refactor HF tests

---
 tests/test_hf_circuit.py | 107 ++++++++++++++-------------------------
 1 file changed, 37 insertions(+), 70 deletions(-)

diff --git a/tests/test_hf_circuit.py b/tests/test_hf_circuit.py
index 839f6ad..bd5bbe3 100644
--- a/tests/test_hf_circuit.py
+++ b/tests/test_hf_circuit.py
@@ -1,70 +1,37 @@
-"""
-Test HF reference circuit ansatz
-"""
-
-# import numpy as np
-import pytest
-
-from qibochem.ansatz.hf_reference import hf_circuit
-from qibochem.driver.molecule import Molecule
-from qibochem.measurement.expectation import expectation
-
-
-def test_jw_circuit():
-    """Tests the HF circuit with the Jordan-Wigner mapping"""
-    # Hardcoded benchmark results
-    h2_ref_energy = -1.117349035
-
-    h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
-    h2.run_pyscf()
-
-    # JW-HF circuit
-    circuit = hf_circuit(h2.nso, h2.nelec, ferm_qubit_map=None)
-
-    # Molecular Hamiltonian and the HF expectation value
-    hamiltonian = h2.hamiltonian()
-    hf_energy = expectation(circuit, hamiltonian)
-
-    # assert h2.e_hf == pytest.approx(hf_energy)
-    assert h2_ref_energy == pytest.approx(hf_energy)
-
-
-def test_bk_circuit_1():
-    """Tests the HF circuit with the Brayvi-Kitaev mapping for H2"""
-    # Hardcoded benchmark results
-    h2_ref_energy = -1.117349035
-
-    h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
-    h2.run_pyscf()
-
-    # JW-HF circuit
-    circuit = hf_circuit(h2.nso, h2.nelec, ferm_qubit_map="bk")
-
-    # Molecular Hamiltonian and the HF expectation value
-    hamiltonian = h2.hamiltonian(ferm_qubit_map="bk")
-    hf_energy = expectation(circuit, hamiltonian)
-
-    # assert h2.e_hf == pytest.approx(hf_energy)
-    assert h2_ref_energy == pytest.approx(hf_energy)
-
-
-def test_bk_circuit_2():
-    """Tests the HF circuit with the Brayvi-Kitaev mapping for LiH"""
-    # Hardcoded benchmark results
-    lih = Molecule([("Li", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 1.3))])
-    lih.run_pyscf()
-
-    # JW-HF circuit
-    circuit = hf_circuit(lih.nso, lih.nelec, ferm_qubit_map="bk")
-
-    # Molecular Hamiltonian and the HF expectation value
-    hamiltonian = lih.hamiltonian(ferm_qubit_map="bk")
-    hf_energy = expectation(circuit, hamiltonian)
-
-    assert lih.e_hf == pytest.approx(hf_energy)
-
-
-def test_mapping_error():
-    """Tests the HF circuit with an incorrect mapping"""
-    with pytest.raises(KeyError):
-        hf_circuit(4, 2, ferm_qubit_map="incorrect")
+"""
+Test HF reference circuit ansatz
+"""
+
+import pytest
+
+from qibochem.ansatz.hf_reference import hf_circuit
+from qibochem.driver.molecule import Molecule
+from qibochem.measurement.expectation import expectation
+
+
+@pytest.mark.parametrize(
+    "mapping,",
+    [
+        None,  # JW mapping
+        "bk",  # BK mapping
+    ],
+)
+def test_h2(mapping):
+    """Tests the HF circuit for H2"""
+    # Hardcoded benchmark results
+    h2_ref_energy = -1.117349035
+
+    h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
+    h2.run_pyscf()
+    hamiltonian = h2.hamiltonian(ferm_qubit_map=mapping)
+    circuit = hf_circuit(h2.nso, h2.nelec, ferm_qubit_map=mapping)
+    hf_energy = expectation(circuit, hamiltonian)
+
+    # assert h2.e_hf == pytest.approx(hf_energy)
+    assert pytest.approx(hf_energy) == h2_ref_energy
+
+
+def test_mapping_error():
+    """Tests the HF circuit with an incorrect mapping"""
+    with pytest.raises(KeyError):
+        hf_circuit(4, 2, ferm_qubit_map="incorrect")

From a52ad25c16744b7c5e9c5843f8838c4b703ac34f Mon Sep 17 00:00:00 2001
From: chmwzc <70616433+chmwzc@users.noreply.github.com>
Date: Tue, 27 Feb 2024 15:05:12 +0800
Subject: [PATCH 03/11] Clean up spacing

---
 tests/test_expectation_samples.py | 27 +++------------------------
 1 file changed, 3 insertions(+), 24 deletions(-)

diff --git a/tests/test_expectation_samples.py b/tests/test_expectation_samples.py
index 5eed036..087c6df 100644
--- a/tests/test_expectation_samples.py
+++ b/tests/test_expectation_samples.py
@@ -19,30 +19,9 @@
 @pytest.mark.parametrize(
     "terms,gates_to_add,expected,threshold",
     [
-        (
-            Z(0),
-            [
-                gates.X(0),
-            ],
-            -1.0,
-            None,
-        ),
-        (
-            Z(0) * Z(1),
-            [
-                gates.X(0),
-            ],
-            -1.0,
-            None,
-        ),
-        (
-            X(0),
-            [
-                gates.H(0),
-            ],
-            1.0,
-            None,
-        ),
+        (Z(0), [gates.X(0)], -1.0, None),
+        (Z(0) * Z(1), [gates.X(0)], -1.0, None),
+        (X(0), [gates.H(0)], 1.0, None),
         (X(0), [gates.X(0), gates.X(0)], 0.0, 0.05),
     ],
 )

From 027c92a2bd03b29aa174ff3944d6637202775195 Mon Sep 17 00:00:00 2001
From: Wong Zi Cheng <70616433+chmwzc@users.noreply.github.com>
Date: Wed, 28 Feb 2024 02:31:17 +0000
Subject: [PATCH 04/11] Refactor UCC tests (in-progress)

More for testing git push/pull after the problems in the week so far
---
 tests/test_ucc.py | 81 +++++++++++++++--------------------------------
 1 file changed, 25 insertions(+), 56 deletions(-)

diff --git a/tests/test_ucc.py b/tests/test_ucc.py
index 9a5d688..a28c935 100644
--- a/tests/test_ucc.py
+++ b/tests/test_ucc.py
@@ -17,64 +17,33 @@
 from qibochem.measurement.expectation import expectation
 
 
-def test_generate_excitations_1():
-    ex1 = generate_excitations(1, [0, 1], [2, 3])
-    ref_ex1 = np.array([[0, 2], [1, 3]])
-
-    assert np.allclose(ex1, ref_ex1)
-
-
-def test_generate_excitations_2():
-    ex2 = generate_excitations(2, [0, 1], [2, 3])
-    ref_ex2 = np.array([[0, 1, 2, 3]])
-
-    assert np.allclose(ex2, ref_ex2)
-
-
-def test_generate_excitations_3():
-    ex3 = generate_excitations(3, [0, 1], [2, 3])
-
-    assert np.allclose(ex3, [[]])
-
-
-def test_sort_excitations_1():
-    ex1 = generate_excitations(1, [0, 1], [2, 3, 4, 5])
-    sorted_excitations = sort_excitations(ex1)
-    ref_sorted_ex1 = np.array([[0, 2], [1, 3], [0, 4], [1, 5]])
-
-    assert np.allclose(sorted_excitations, ref_sorted_ex1)
-
-
-def test_sort_excitations_2():
-    ex2 = generate_excitations(2, [0, 1, 2, 3], [4, 5, 6, 7])
-    sorted_excitations = sort_excitations(ex2)
-    ref_sorted_ex2 = np.array(
-        [
-            [0, 1, 4, 5],
-            [0, 1, 6, 7],
-            [2, 3, 4, 5],
-            [2, 3, 6, 7],
-            [0, 1, 4, 7],
-            [0, 1, 5, 6],
-            [2, 3, 4, 7],
-            [2, 3, 5, 6],
-            [0, 3, 4, 5],
-            [1, 2, 4, 5],
-            [0, 3, 6, 7],
-            [1, 2, 6, 7],
-            [0, 2, 4, 6],
-            [1, 3, 5, 7],
-            [0, 3, 4, 7],
-            [0, 3, 5, 6],
-            [1, 2, 4, 7],
-            [1, 2, 5, 6],
-        ]
-    )
-
-    assert np.allclose(sorted_excitations, ref_sorted_ex2)
+@pytest.mark.parametrize(
+    "order,excite_from,excite_to,expected",
+    [
+        (1, [0, 1], [2, 3], [[0, 2], [1, 3]]),
+        (2, [2, 3], [4, 5], [[2, 3, 4, 5]]),
+        (3, [0, 1], [2, 3], [[]]),
+    ],
+)
+def test_generate_excitations(order, excite_from, excite_to, expected):
+    """Test generation of all possible excitations between two lists of orbitals"""
+    test = generate_excitations(order, excite_from, excite_to)
+    assert test == expected
+
+
+@pytest.mark.parametrize(
+    "order,excite_from,excite_to,expected",
+    [
+        (1, [0, 1], [2, 3, 4, 5], [[0, 2], [1, 3], [0, 4], [1, 5]]),
+        (2, [0, 1], [2, 3, 4, 5], [[0, 1, 2, 3], [0, 1, 4, 5], [0, 1, 2, 5], [0, 1, 3, 4]]),
+    ],
+)
+def test_sort_excitations(order, excite_from, excite_to, expected):
+    test = sort_excitations(generate_excitations(order, excite_from, excite_to))
+    assert test == expected
 
 
-def test_sort_excitations_3():
+def test_sort_excitations_triples():
     with pytest.raises(NotImplementedError):
         sort_excitations([[1, 2, 3, 4, 5, 6]])
 

From 23ce69f1fff57554c0cd37aa41862d55b9ea7948 Mon Sep 17 00:00:00 2001
From: Wong Zi Cheng <70616433+chmwzc@users.noreply.github.com>
Date: Wed, 28 Feb 2024 02:37:51 +0000
Subject: [PATCH 05/11] Add docstring to test_ucc.py

Still testing whether can push without issues...
---
 tests/test_ucc.py | 4 ++++
 1 file changed, 4 insertions(+)

diff --git a/tests/test_ucc.py b/tests/test_ucc.py
index a28c935..d8c74e1 100644
--- a/tests/test_ucc.py
+++ b/tests/test_ucc.py
@@ -1,3 +1,7 @@
+"""
+Tests for the UCC ansatz and related functions
+"""
+
 from functools import partial
 
 import numpy as np

From 4a876ed72c4127ef2f2b1b30ccd6ca037ec78237 Mon Sep 17 00:00:00 2001
From: Wong Zi Cheng <70616433+chmwzc@users.noreply.github.com>
Date: Wed, 28 Feb 2024 03:45:24 +0000
Subject: [PATCH 06/11] Refactor UCC tests

---
 tests/test_ucc.py | 106 ++++++++++++++++------------------------------
 1 file changed, 36 insertions(+), 70 deletions(-)

diff --git a/tests/test_ucc.py b/tests/test_ucc.py
index d8c74e1..b3d5581 100644
--- a/tests/test_ucc.py
+++ b/tests/test_ucc.py
@@ -7,7 +7,6 @@
 import numpy as np
 import pytest
 from qibo import gates
-from scipy.optimize import minimize
 
 from qibochem.ansatz.hf_reference import hf_circuit
 from qibochem.ansatz.ucc import (
@@ -17,8 +16,7 @@
     ucc_ansatz,
     ucc_circuit,
 )
-from qibochem.driver.molecule import Molecule
-from qibochem.measurement.expectation import expectation
+from qibochem.driver import Molecule
 
 
 @pytest.mark.parametrize(
@@ -65,84 +63,52 @@ def test_mp2_amplitude_doubles():
     assert np.isclose(l, ref_l)
 
 
-def test_ucc_jw_singles():
-    """Build a UCC singles (JW) circuit"""
-    rx_gate = partial(gates.RX, theta=-0.5 * np.pi, trainable=False)
-    cnot_cascade = [gates.CNOT(2, 1), gates.CNOT(1, 0), gates.RZ(0, 0.0), gates.CNOT(1, 0), gates.CNOT(2, 1)]
-    basis_rotation_gates = ([rx_gate(0), gates.H(2)], [gates.H(0), rx_gate(2)])
-
-    # Build control list of gates
-    nested_gate_list = [gate_list + cnot_cascade + gate_list for gate_list in basis_rotation_gates]
-    gate_list = [_gate for gate_list in nested_gate_list for _gate in gate_list]
-
-    # Test ucc_function
-    circuit = ucc_circuit(4, [0, 2], ferm_qubit_map="jw")
-    # Check gates are correct
-    assert all(
-        control.name == test.name and control.target_qubits == test.target_qubits
-        for control, test in zip(gate_list, list(circuit.queue))
-    )
-    # Check that only two parametrised gates
-    assert len(circuit.get_parameters()) == 2
-
-
-def test_ucc_jw_doubles():
-    """Build a UCC doubles (JW) circuit"""
-    rx_gate = partial(gates.RX, theta=-0.5 * np.pi, trainable=False)
-    cnot_cascade = [
-        gates.CNOT(3, 2),
-        gates.CNOT(2, 1),
-        gates.CNOT(1, 0),
-        gates.RZ(0, 0.0),
-        gates.CNOT(1, 0),
-        gates.CNOT(2, 1),
-        gates.CNOT(3, 2),
+@pytest.mark.parametrize(
+    "excitation,mapping,basis_rotations",
+    [
+        ([0, 2], None, ([("Y", 0), ("X", 2)], [("X", 0), ("Y", 2)])),  # JW singles
+        (
+            [0, 1, 2, 3],
+            None,
+            (
+                [("X", 0), ("X", 1), ("Y", 2), ("X", 3)],
+                [("Y", 0), ("Y", 1), ("Y", 2), ("X", 3)],
+                [("Y", 0), ("X", 1), ("X", 2), ("X", 3)],
+                [("X", 0), ("Y", 1), ("X", 2), ("X", 3)],
+                [("Y", 0), ("X", 1), ("Y", 2), ("Y", 3)],
+                [("X", 0), ("Y", 1), ("Y", 2), ("Y", 3)],
+                [("X", 0), ("X", 1), ("X", 2), ("Y", 3)],
+                [("Y", 0), ("Y", 1), ("X", 2), ("Y", 3)],
+            ),
+        ),  # JW doubles
+        ([0, 2], "bk", ([("X", 0), ("Y", 1), ("X", 2)], [("Y", 0), ("Y", 1), ("Y", 2)])),  # BK singles
+    ],
+)
+def test_ucc_circuit(excitation, mapping, basis_rotations):
+    """Build a UCC circuit with only one excitation"""
+    gate_dict = {"X": gates.H, "Y": partial(gates.RX, theta=-0.5 * np.pi, trainable=False)}
+    # Build the list of basis rotation gates
+    basis_rotation_gates = [
+        [gate_dict[_gate[0]](_gate[1]) for _gate in basis_rotation] for basis_rotation in basis_rotations
     ]
-    basis_rotation_gates = (
-        [gates.H(0), gates.H(1), rx_gate(2), gates.H(3)],
-        [rx_gate(0), rx_gate(1), rx_gate(2), gates.H(3)],
-        [rx_gate(0), gates.H(1), gates.H(2), gates.H(3)],
-        [gates.H(0), rx_gate(1), gates.H(2), gates.H(3)],
-        [rx_gate(0), gates.H(1), rx_gate(2), rx_gate(3)],
-        [gates.H(0), rx_gate(1), rx_gate(2), rx_gate(3)],
-        [gates.H(0), gates.H(1), gates.H(2), rx_gate(3)],
-        [rx_gate(0), rx_gate(1), gates.H(2), rx_gate(3)],
-    )
-
-    # Build control list of gates
-    nested_gate_list = [gate_list + cnot_cascade + gate_list for gate_list in basis_rotation_gates]
-    gate_list = [_gate for gate_list in nested_gate_list for _gate in gate_list]
-
-    # Test ucc_function
-    circuit = ucc_circuit(4, [0, 1, 2, 3], ferm_qubit_map="jw")
-    # Check gates are correct
-    assert all(
-        control.name == test.name and control.target_qubits == test.target_qubits
-        for control, test in zip(gate_list, list(circuit.queue))
-    )
-    # Check that only two parametrised gates
-    assert len(circuit.get_parameters()) == 8
-
-
-def test_ucc_bk_singles():
-    """Build a UCC doubles (BK) circuit"""
-    rx_gate = partial(gates.RX, theta=-0.5 * np.pi, trainable=False)
-    cnot_cascade = [gates.CNOT(2, 1), gates.CNOT(1, 0), gates.RZ(0, 0.0), gates.CNOT(1, 0), gates.CNOT(2, 1)]
-    basis_rotation_gates = ([gates.H(0), rx_gate(1), gates.H(2)], [rx_gate(0), rx_gate(1), rx_gate(2)])
+    # Build the CNOT cascade manually
+    cnot_cascade = [gates.CNOT(_i, _i - 1) for _i in range(excitation[-1], excitation[0], -1)]
+    cnot_cascade = cnot_cascade + [gates.RZ(excitation[0], 0.0)]
+    cnot_cascade = cnot_cascade + [gates.CNOT(_i + 1, _i) for _i in range(excitation[0], excitation[-1])]
 
     # Build control list of gates
     nested_gate_list = [gate_list + cnot_cascade + gate_list for gate_list in basis_rotation_gates]
     gate_list = [_gate for gate_list in nested_gate_list for _gate in gate_list]
 
     # Test ucc_function
-    circuit = ucc_circuit(4, [0, 2], ferm_qubit_map="bk")
+    circuit = ucc_circuit(4, excitation, ferm_qubit_map=mapping)
     # Check gates are correct
     assert all(
         control.name == test.name and control.target_qubits == test.target_qubits
         for control, test in zip(gate_list, list(circuit.queue))
     )
     # Check that only two parametrised gates
-    assert len(circuit.get_parameters()) == 2
+    assert len(circuit.get_parameters()) == len(basis_rotations)
 
 
 def test_ucc_ferm_qubit_map_error():
@@ -156,14 +122,14 @@ def test_ucc_parameter_coefficients():
     # UCC-JW singles
     control_values = (-1.0, 1.0)
     coeffs = []
-    circuit = ucc_circuit(2, [0, 1], coeffs=coeffs)
+    _circuit = ucc_circuit(2, [0, 1], coeffs=coeffs)
     # Check that the signs of the coefficients have been saved
     assert all(control == test for control, test in zip(control_values, coeffs))
 
     # UCC-JW doubles
     control_values = (-0.25, 0.25, 0.25, 0.25, -0.25, -0.25, -0.25, 0.25)
     coeffs = []
-    circuit = ucc_circuit(4, [0, 1, 2, 3], coeffs=coeffs)
+    _circuit = ucc_circuit(4, [0, 1, 2, 3], coeffs=coeffs)
     # Check that the signs of the coefficients have been saved
     assert all(control == test for control, test in zip(control_values, coeffs))
 

From cdeb1360d4b2463ca309d333c01a2cf1a128bd76 Mon Sep 17 00:00:00 2001
From: Wong Zi Cheng <70616433+chmwzc@users.noreply.github.com>
Date: Wed, 28 Feb 2024 07:10:22 +0000
Subject: [PATCH 07/11] Refactor test_molecule.py

---
 tests/test_expectation_samples.py |   3 +-
 tests/test_molecule.py            | 365 ++++++++++++++----------------
 2 files changed, 176 insertions(+), 192 deletions(-)

diff --git a/tests/test_expectation_samples.py b/tests/test_expectation_samples.py
index 087c6df..bd2ba4d 100644
--- a/tests/test_expectation_samples.py
+++ b/tests/test_expectation_samples.py
@@ -2,11 +2,10 @@
 Test expectation functionality
 """
 
-# import numpy as np
 import pytest
 from qibo import Circuit, gates
 from qibo.hamiltonians import SymbolicHamiltonian
-from qibo.symbols import X, Y, Z
+from qibo.symbols import X, Z
 
 from qibochem.driver import Molecule
 from qibochem.measurement import expectation
diff --git a/tests/test_molecule.py b/tests/test_molecule.py
index 198b45b..bbf765e 100644
--- a/tests/test_molecule.py
+++ b/tests/test_molecule.py
@@ -9,49 +9,55 @@
 import pytest
 from qibo import gates, models
 from qibo.hamiltonians import SymbolicHamiltonian
-from qibo.symbols import X, Y, Z
+from qibo.symbols import Z
 
-from qibochem.driver import hamiltonian
-from qibochem.driver.molecule import Molecule
+from qibochem.driver import Molecule
+from qibochem.driver.hamiltonian import (  # fermionic_hamiltonian,; qubit_hamiltonian,; symbolic_hamiltonian,
+    parse_pauli_string,
+)
 from qibochem.measurement.expectation import expectation
 
 
-def test_run_pyscf():
-    """PySCF driver"""
-    # Hardcoded benchmark results
-    # Change to run PySCF directly during a test?
-    h2_ref_energy = -1.117349035
-    h2_ref_hcore = np.array([[-1.14765024, -1.00692423], [-1.00692423, -1.14765024]])
-
-    h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
-    h2.run_pyscf()
-    assert h2.e_hf == pytest.approx(h2_ref_energy)
-    assert np.allclose(h2.hcore, h2_ref_hcore)
-
-
-def test_run_pyscf_molecule_xyz():
-    """Pyscf driver with xyz file"""
-    file_path = Path("tests", "data", "lih.xyz")
-    if not file_path.is_file():
-        with open(file_path, "w") as file_handler:
-            file_handler.write("2\n0 1\nLi 0.0 0.0 0.0\nH 0.0 0.0 1.2\n")
-    lih_ref_energy = -7.83561582555692
-    lih = Molecule(xyz_file=file_path)
-    lih.run_pyscf()
-
-    assert lih.e_hf == pytest.approx(lih_ref_energy)
-
+@pytest.mark.parametrize(
+    "xyz_file,expected",
+    [
+        (None, -1.117349035),
+        ("lih.xyz", -7.83561582555692),
+        ("h2.xyz", -1.117349035),
+    ],
+)
+def test_pyscf_driver(xyz_file, expected):
+    if xyz_file is None:
+        mol = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
+    else:
+        file_path = Path("tests", "data") / Path(xyz_file)
+        # In case .xyz files somehow not found
+        if not file_path.is_file():
+            with open(file_path, "w") as file_handler:
+                if xyz_file == "lih.xyz":
+                    file_handler.write("2\n0 1\nLi 0.0 0.0 0.0\nH 0.0 0.0 1.2\n")
+                elif xyz_file == "h2.xyz":
+                    file_handler.write("2\n \nH 0.0 0.0 0.0\nH 0.0 0.0 0.7\n")
+                else:
+                    file_handler.write("2\n \nH 0.0 0.0 0.0\nH 0.0 0.0 0.7\n")
+        # Define Molecule using .xyz file
+        mol = Molecule(xyz_file=file_path)
+    # Run PySCF and check that the HF energy matches
+    mol.run_pyscf()
+    assert mol.e_hf == pytest.approx(expected)
 
-def test_run_pyscf_molecule_xyz_charged():
-    file_path = Path("tests", "data", "h2.xyz")
-    if not file_path.is_file():
-        with open(file_path, "w") as file_handler:
-            file_handler.write("2\n \nH 0.0 0.0 0.0\nH 0.0 0.0 0.7\n")
-    h2_ref_energy = -1.117349035
-    h2 = Molecule(xyz_file=file_path)
-    h2.run_pyscf()
 
-    assert h2.e_hf == pytest.approx(h2_ref_energy)
+# Commenting out since not actively supporting PSI4 at the moment
+# @pytest.mark.skip(reason="Psi4 doesn't offer pip install, so needs to be installed through conda or manually.")
+# def test_run_psi4():
+#     """PSI4 driver"""
+#     # Hardcoded benchmark results
+#     h2_ref_energy = -1.117349035
+#
+#     h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
+#     h2.run_psi4()
+#
+#     assert h2.e_hf == pytest.approx(h2_ref_energy)
 
 
 def test_molecule_custom_basis():
@@ -60,20 +66,21 @@ def test_molecule_custom_basis():
     assert np.isclose(mol.e_hf, -7.94129296352493)
 
 
-def test_define_active_space():
+@pytest.mark.parametrize(
+    "active,frozen,expected",
+    [
+        (None, None, (list(range(6)), [])),  # Default arguments: Nothing given
+        ([1, 2, 5], None, ([1, 2, 5], [0])),  # Default frozen argument if active given
+        (None, [0], (list(range(1, 6)), [0])),  # Default active argument if frozen given
+        ([0, 1, 2, 3], [], (list(range(4)), [])),  # active, frozen arguments both given
+        ([1, 2, 3], [0], (list(range(1, 4)), [0])),  # active, frozen arguments both given
+    ],
+)
+def test_define_active_space(active, frozen, expected):
     mol = Molecule([("Li", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 1.2))])
     mol.nalpha = 2
     mol.norb = 6
-    # Default arguments: Nothing given
-    assert mol._active_space(None, None) == (list(range(mol.norb)), [])
-    # Default frozen argument if active given
-    assert mol._active_space([1, 2, 5], None) == ([1, 2, 5], [0])
-    # Default active argument if frozen given
-    assert mol._active_space(None, [0]) == (list(range(1, 6)), [0])
-    # active, frozen arguments both given
-    assert mol._active_space([0, 1, 2, 3], []) == (list(range(4)), [])
-    # active, frozen arguments both given
-    assert mol._active_space([1, 2, 3], [0]) == (list(range(1, 4)), [0])
+    assert mol._active_space(active, frozen) == expected
 
 
 def test_define_active_space_assertions():
@@ -107,160 +114,138 @@ def test_hf_embedding():
     assert np.allclose(mol.embed_tei, mol.tei[:dim, :dim, :dim, :dim])
 
 
-def test_fermionic_hamiltonian():
-    h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
-    h2.run_pyscf()
-    h2_ferm_ham = hamiltonian.fermionic_hamiltonian(h2.oei, h2.tei, h2.e_nuc)
-
-    # a^\dagger_0 a_0
-    ref_one_body_tensor = np.array(
-        [
-            [-1.27785301, 0.0, 0.0, 0.0],
-            [0.0, -1.27785301, 0.0, 0.0],
-            [0.0, 0.0, -0.4482997, 0.0],
-            [0.0, 0.0, 0.0, -0.4482997],
-        ]
-    )
-
-    assert np.isclose(h2_ferm_ham[()], 0.7559674441714287)
-    assert np.allclose(h2_ferm_ham.one_body_tensor, ref_one_body_tensor)
-
-
-def test_fermionic_hamiltonian_2():
-    h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7414))])
-    h2.run_pyscf()
-
-    h2_ferm_ham_1 = h2.hamiltonian("f", h2.oei, h2.tei)
-    h2_qub_ham_jw_1 = h2.hamiltonian("q", h2.oei, h2.tei, ferm_qubit_map="jw")
-    h2_qub_ham_bk_1 = h2.hamiltonian("q", h2.oei, h2.tei, ferm_qubit_map="bk")
-    # h2_mol_ham has format of InteractionOperator
-    h2_mol_ham = hamiltonian.fermionic_hamiltonian(h2.oei, h2.tei, h2.e_nuc)
-    h2_ferm_ham_2 = openfermion.transforms.get_fermion_operator(h2_mol_ham)
-    h2_qub_ham_jw_2 = openfermion.jordan_wigner(h2_mol_ham)
-    h2_qub_ham_bk_2 = openfermion.bravyi_kitaev(h2_mol_ham)
-
-    assert h2_ferm_ham_2.isclose(h2_ferm_ham_1)
-    assert h2_qub_ham_jw_2.isclose(h2_qub_ham_jw_1)
-    assert h2_qub_ham_bk_2.isclose(h2_qub_ham_bk_1)
-
+@pytest.mark.parametrize(
+    "pauli_string,coeff,expected",
+    [
+        (((0, "X"), (1, "Y")), 0.5, "0.5*X0*Y1"),
+        (None, 0.1, "0.1"),
+    ],
+)
+def test_parse_pauli_string(pauli_string, coeff, expected):
+    test = parse_pauli_string(pauli_string, coeff)
+    assert str(test) == expected
 
-def test_parse_pauli_string_1():
-    pauli_string = ((0, "X"), (1, "Y"))
-    qibo_pauli_string = hamiltonian.parse_pauli_string(pauli_string, 0.5)
-    ref_pauli_string = "0.5*X0*Y1"
-    assert str(qibo_pauli_string) == ref_pauli_string
-
-
-def test_parse_pauli_string_2():
-    qibo_pauli_string = hamiltonian.parse_pauli_string(None, 0.1)
-    assert str(qibo_pauli_string) == "0.1"
-
-
-def test_qubit_hamiltonian():
-    h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
-    h2.run_pyscf()
-    h2_ferm_ham = hamiltonian.fermionic_hamiltonian(h2.oei, h2.tei, h2.e_nuc)
-    h2_qubit_ham_jw = hamiltonian.qubit_hamiltonian(h2_ferm_ham, "jw")
-    h2_qubit_ham_bk = hamiltonian.qubit_hamiltonian(h2_ferm_ham, "bk")
-    ref_h2_qubit_ham_jw = {
-        (): -0.04207897647782238,
-        ((0, "Z"),): 0.17771287465139918,
-        ((1, "Z"),): 0.1777128746513992,
-        ((2, "Z"),): -0.24274280513140478,
-        ((3, "Z"),): -0.24274280513140478,
-        ((0, "Z"), (1, "Z")): 0.17059738328801052,
-        ((0, "Z"), (2, "Z")): 0.12293305056183809,
-        ((0, "Z"), (3, "Z")): 0.16768319457718972,
-        ((1, "Z"), (2, "Z")): 0.16768319457718972,
-        ((1, "Z"), (3, "Z")): 0.12293305056183809,
-        ((2, "Z"), (3, "Z")): 0.17627640804319608,
-        ((0, "X"), (1, "X"), (2, "Y"), (3, "Y")): -0.04475014401535165,
-        ((0, "X"), (1, "Y"), (2, "Y"), (3, "X")): 0.04475014401535165,
-        ((0, "Y"), (1, "X"), (2, "X"), (3, "Y")): 0.04475014401535165,
-        ((0, "Y"), (1, "Y"), (2, "X"), (3, "X")): -0.04475014401535165,
-    }
-    ref_h2_qubit_ham_bk = {
-        ((0, "Z"),): 0.17771287465139923,
-        (): -0.04207897647782244,
-        ((0, "Z"), (1, "Z")): 0.17771287465139918,
-        ((2, "Z"),): -0.24274280513140484,
-        ((1, "Z"), (2, "Z"), (3, "Z")): -0.24274280513140484,
-        ((0, "Y"), (1, "Z"), (2, "Y")): 0.04475014401535165,
-        ((0, "X"), (1, "Z"), (2, "X")): 0.04475014401535165,
-        ((0, "X"), (1, "Z"), (2, "X"), (3, "Z")): 0.04475014401535165,
-        ((0, "Y"), (1, "Z"), (2, "Y"), (3, "Z")): 0.04475014401535165,
-        ((1, "Z"),): 0.17059738328801052,
-        ((0, "Z"), (2, "Z")): 0.12293305056183809,
-        ((0, "Z"), (1, "Z"), (2, "Z")): 0.16768319457718972,
-        ((0, "Z"), (1, "Z"), (2, "Z"), (3, "Z")): 0.16768319457718972,
-        ((0, "Z"), (2, "Z"), (3, "Z")): 0.12293305056183809,
-        ((1, "Z"), (3, "Z")): 0.17627640804319608,
-    }
-
-    jw_array = np.array([terms for terms in h2_qubit_ham_jw.terms.values()])
-    bk_array = np.array([terms for terms in h2_qubit_ham_bk.terms.values()])
-
-    ref_jw_array = np.array([terms for terms in ref_h2_qubit_ham_jw.values()])
-    ref_bk_array = np.array([terms for terms in ref_h2_qubit_ham_bk.values()])
-
-    assert np.allclose(jw_array, ref_jw_array)
-    assert np.allclose(bk_array, ref_bk_array)
-
-    # incorrect mapping circuit
-    with pytest.raises(KeyError):
-        hamiltonian.qubit_hamiltonian(h2_ferm_ham, "incorrect")
 
+# def test_fermionic_hamiltonian():
+#     h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
+#     h2.run_pyscf()
+#     h2_ferm_ham = h2.hamiltonian("f")
+#
+#     # a^\dagger_0 a_0
+#     ref_one_body_tensor = np.diag([-1.27785301, -1.27785301, -0.4482997, -0.4482997])
+#
+#     assert np.isclose(h2_ferm_ham.constant, 0.7559674441714287) # N-N repulsion
+#     assert np.allclose(h2_ferm_ham.one_body_tensor, ref_one_body_tensor)
 
-def test_symbolic_hamiltonian():
-    h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
-    h2.run_pyscf()
-    h2_ferm_ham = hamiltonian.fermionic_hamiltonian(h2.oei, h2.tei, h2.e_nuc)
-    h2_qubit_ham = hamiltonian.qubit_hamiltonian(h2_ferm_ham, "jw")
-    h2_sym_ham = hamiltonian.symbolic_hamiltonian(h2_qubit_ham)
-    ref_ham = (
-        -0.0420789764778224
-        - 0.0447501440153516 * X(0) * X(1) * Y(2) * Y(3)
-        + 0.0447501440153516 * X(0) * Y(1) * Y(2) * X(3)
-        + 0.0447501440153516 * Y(0) * X(1) * X(2) * Y(3)
-        - 0.0447501440153516 * Y(0) * Y(1) * X(2) * X(3)
-        + 0.177712874651399 * Z(0)
-        + 0.170597383288011 * Z(0) * Z(1)
-        + 0.122933050561838 * Z(0) * Z(2)
-        + 0.16768319457719 * Z(0) * Z(3)
-        + 0.177712874651399 * Z(1)
-        + 0.16768319457719 * Z(1) * Z(2)
-        + 0.122933050561838 * Z(1) * Z(3)
-        - 0.242742805131405 * Z(2)
-        + 0.176276408043196 * Z(2) * Z(3)
-        - 0.242742805131405 * Z(3)
-    )
-    ref_sym_ham = SymbolicHamiltonian(ref_ham)
 
-    assert np.allclose(h2_sym_ham.matrix, ref_sym_ham.matrix)
+# def test_fermionic_hamiltonian_2():
+#     h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7414))])
+#     h2.run_pyscf()
+#
+#     h2_ferm_ham_1 = h2.hamiltonian("f", h2.oei, h2.tei)
+#     h2_qub_ham_jw_1 = h2.hamiltonian("q", h2.oei, h2.tei, ferm_qubit_map="jw")
+#     h2_qub_ham_bk_1 = h2.hamiltonian("q", h2.oei, h2.tei, ferm_qubit_map="bk")
+#     # h2_mol_ham has format of InteractionOperator
+#     h2_mol_ham = fermionic_hamiltonian(h2.oei, h2.tei, h2.e_nuc)
+#     h2_ferm_ham_2 = openfermion.transforms.get_fermion_operator(h2_mol_ham)
+#     h2_qub_ham_jw_2 = openfermion.jordan_wigner(h2_mol_ham)
+#     h2_qub_ham_bk_2 = openfermion.bravyi_kitaev(h2_mol_ham)
+#
+#     assert h2_ferm_ham_2.isclose(h2_ferm_ham_1)
+#     assert h2_qub_ham_jw_2.isclose(h2_qub_ham_jw_1)
+#     assert h2_qub_ham_bk_2.isclose(h2_qub_ham_bk_1)
+
+
+# @pytest.mark.parametrize(
+#     "mapping,expected_operators",
+#     [
+#         (
+#             None,
+#             [
+#                 ((), -0.04207897647782238),
+#                 (((0, "Z")), 0.17771287465139918),
+#                 (((1, "Z")), 0.1777128746513992),
+#                 (((2, "Z")), -0.24274280513140478),
+#                 (((3, "Z")), -0.24274280513140478),
+#                 (((0, "Z"), (1, "Z")), 0.17059738328801052),
+#                 (((0, "Z"), (2, "Z")), 0.12293305056183809),
+#                 (((0, "Z"), (3, "Z")), 0.16768319457718972),
+#                 (((1, "Z"), (2, "Z")), 0.16768319457718972),
+#                 (((1, "Z"), (3, "Z")), 0.12293305056183809),
+#                 (((2, "Z"), (3, "Z")), 0.17627640804319608),
+#                 (((0, "X"), (1, "X"), (2, "Y"), (3, "Y")), -0.04475014401535165),
+#                 (((0, "X"), (1, "Y"), (2, "Y"), (3, "X")), 0.04475014401535165),
+#                 (((0, "Y"), (1, "X"), (2, "X"), (3, "Y")), 0.04475014401535165),
+#                 (((0, "Y"), (1, "Y"), (2, "X"), (3, "X")), -0.04475014401535165),
+#             ]
+#         ), # H2 JW mapping
+#         (
+#             "bk",
+#             [
+#                 ((), -0.04207897647782244),
+#                 (((0, "Z"),), 0.17771287465139923),
+#                 (((0, "Z"), (1, "Z")), 0.17771287465139918),
+#                 (((2, "Z"),), -0.24274280513140484),
+#                 (((1, "Z"), (2, "Z"), (3, "Z")), -0.24274280513140484),
+#                 (((0, "Y"), (1, "Z"), (2, "Y")), 0.04475014401535165),
+#                 (((0, "X"), (1, "Z"), (2, "X")), 0.04475014401535165),
+#                 (((0, "X"), (1, "Z"), (2, "X"), (3, "Z")), 0.04475014401535165),
+#                 (((0, "Y"), (1, "Z"), (2, "Y"), (3, "Z")), 0.04475014401535165),
+#                 (((1, "Z"),), 0.17059738328801052),
+#                 (((0, "Z"), (2, "Z")), 0.12293305056183809),
+#                 (((0, "Z"), (1, "Z"), (2, "Z")), 0.16768319457718972),
+#                 (((0, "Z"), (1, "Z"), (2, "Z"), (3, "Z")), 0.16768319457718972),
+#                 (((0, "Z"), (2, "Z"), (3, "Z")), 0.12293305056183809),
+#                 (((1, "Z"), (3, "Z")), 0.17627640804319608),
+#             ],
+#         ), # H2 BK mapping
+#     ],
+# )
+# def test_qubit_hamiltonian(mapping, expected_operators):
+#     control = sum(openfermion.QubitOperator(_op, coeff) for _op, coeff in expected_operators)
+#
+#     h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
+#     h2.run_pyscf()
+#
+#     h2_qubit_hamiltonian = h2.hamiltonian("q", ferm_qubit_map=mapping)
+#     assert h2_qubit_hamiltonian.isclose(control)
 
 
-def test_hamiltonian_input_error():
+# def test_symbolic_hamiltonian():
+#     h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
+#     h2.run_pyscf()
+#     hamiltonian = h2.hamiltonian()
+#     ref_sym_ham = SymbolicHamiltonian(
+#         -0.0420789764778224
+#         - 0.0447501440153516 * X(0) * X(1) * Y(2) * Y(3)
+#         + 0.0447501440153516 * X(0) * Y(1) * Y(2) * X(3)
+#         + 0.0447501440153516 * Y(0) * X(1) * X(2) * Y(3)
+#         - 0.0447501440153516 * Y(0) * Y(1) * X(2) * X(3)
+#         + 0.177712874651399 * Z(0)
+#         + 0.170597383288011 * Z(0) * Z(1)
+#         + 0.122933050561838 * Z(0) * Z(2)
+#         + 0.16768319457719 * Z(0) * Z(3)
+#         + 0.177712874651399 * Z(1)
+#         + 0.16768319457719 * Z(1) * Z(2)
+#         + 0.122933050561838 * Z(1) * Z(3)
+#         - 0.242742805131405 * Z(2)
+#         + 0.176276408043196 * Z(2) * Z(3)
+#         - 0.242742805131405 * Z(3)
+#     )
+#     assert np.allclose(hamiltonian.matrix, ref_sym_ham.matrix)
+
+
+def test_hamiltonian_input_errors():
     h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
     h2.e_nuc = 0.0
     h2.oei = np.random.rand(4, 4)
     h2.tei = np.random.rand(4, 4, 4, 4)
+    # Hamiltonian type error
     with pytest.raises(NameError):
         h2.hamiltonian("ihpc")
-
-
-# Commenting out since not actively supporting PSI4 at the moment
-# @pytest.mark.skip(reason="Psi4 doesn't offer pip install, so needs to be installed through conda or manually.")
-# def test_run_psi4():
-#     """PSI4 driver"""
-#     # Hardcoded benchmark results
-#     h2_ref_energy = -1.117349035
-#     h2_ref_hcore = np.array([[-1.14765024, -1.00692423], [-1.00692423, -1.14765024]])
-#
-#     h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
-#     h2.run_psi4()
-#
-#     assert h2.e_hf == pytest.approx(h2_ref_energy)
-#     assert np.allclose(h2.hcore, h2_ref_hcore)
+    # Fermion to qubit mapping error
+    with pytest.raises(KeyError):
+        h2.hamiltonian(ferm_qubit_map="incorrect")
 
 
 def test_expectation_value():

From fc69264485e321db0007ddd117078b9ce4d49844 Mon Sep 17 00:00:00 2001
From: Wong Zi Cheng <70616433+chmwzc@users.noreply.github.com>
Date: Wed, 28 Feb 2024 07:42:14 +0000
Subject: [PATCH 08/11] Fill in missing test coverage

---
 tests/test_molecule.py | 206 ++++++++++++++++++++---------------------
 1 file changed, 100 insertions(+), 106 deletions(-)

diff --git a/tests/test_molecule.py b/tests/test_molecule.py
index bbf765e..010ef6d 100644
--- a/tests/test_molecule.py
+++ b/tests/test_molecule.py
@@ -12,9 +12,7 @@
 from qibo.symbols import Z
 
 from qibochem.driver import Molecule
-from qibochem.driver.hamiltonian import (  # fermionic_hamiltonian,; qubit_hamiltonian,; symbolic_hamiltonian,
-    parse_pauli_string,
-)
+from qibochem.driver.hamiltonian import parse_pauli_string
 from qibochem.measurement.expectation import expectation
 
 
@@ -126,113 +124,109 @@ def test_parse_pauli_string(pauli_string, coeff, expected):
     assert str(test) == expected
 
 
-# def test_fermionic_hamiltonian():
-#     h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
-#     h2.run_pyscf()
-#     h2_ferm_ham = h2.hamiltonian("f")
-#
-#     # a^\dagger_0 a_0
-#     ref_one_body_tensor = np.diag([-1.27785301, -1.27785301, -0.4482997, -0.4482997])
-#
-#     assert np.isclose(h2_ferm_ham.constant, 0.7559674441714287) # N-N repulsion
-#     assert np.allclose(h2_ferm_ham.one_body_tensor, ref_one_body_tensor)
+def test_fermionic_hamiltonian():
+    # Reference result
+    fermion_operator_list = [
+        ("", 0.7559674441714287),
+        ("0^ 0", -1.277853006156875),
+        ("0^ 0^ 0 0", 0.34119476657602105),
+        ("0^ 0^ 2 2", 0.08950028803070331),
+        ("0^ 1^ 1 0", 0.34119476657602105),
+        ("0^ 1^ 3 2", 0.08950028803070331),
+        ("0^ 2^ 0 2", 0.08950028803070331),
+        ("0^ 2^ 2 0", 0.33536638915437944),
+        ("0^ 3^ 1 2", 0.08950028803070331),
+        ("0^ 3^ 3 0", 0.33536638915437944),
+        ("1^ 0^ 0 1", 0.34119476657602105),
+        ("1^ 0^ 2 3", 0.08950028803070331),
+        ("1^ 1", -1.277853006156875),
+        ("1^ 1^ 1 1", 0.34119476657602105),
+        ("1^ 1^ 3 3", 0.08950028803070331),
+        ("1^ 2^ 0 3", 0.08950028803070331),
+        ("1^ 2^ 2 1", 0.33536638915437944),
+        ("1^ 3^ 1 3", 0.08950028803070331),
+        ("1^ 3^ 3 1", 0.33536638915437944),
+        ("2^ 0^ 0 2", 0.3353663891543795),
+        ("2^ 0^ 2 0", 0.08950028803070331),
+        ("2^ 1^ 1 2", 0.3353663891543795),
+        ("2^ 1^ 3 0", 0.08950028803070331),
+        ("2^ 2", -0.448299696101638),
+        ("2^ 2^ 0 0", 0.08950028803070331),
+        ("2^ 2^ 2 2", 0.35255281608639233),
+        ("2^ 3^ 1 0", 0.08950028803070331),
+        ("2^ 3^ 3 2", 0.35255281608639233),
+        ("3^ 0^ 0 3", 0.3353663891543795),
+        ("3^ 0^ 2 1", 0.08950028803070331),
+        ("3^ 1^ 1 3", 0.3353663891543795),
+        ("3^ 1^ 3 1", 0.08950028803070331),
+        ("3^ 2^ 0 1", 0.08950028803070331),
+        ("3^ 2^ 2 3", 0.35255281608639233),
+        ("3^ 3", -0.448299696101638),
+        ("3^ 3^ 1 1", 0.08950028803070331),
+        ("3^ 3^ 3 3", 0.35255281608639233),
+    ]
+    ref_h2_ferm_ham = sum(openfermion.FermionOperator(_op[0], _op[1]) for _op in fermion_operator_list)
+    # Test case
+    h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
+    h2.run_pyscf()
+    h2_ferm_ham = h2.hamiltonian("f")
 
+    assert h2_ferm_ham.isclose(ref_h2_ferm_ham)
 
-# def test_fermionic_hamiltonian_2():
-#     h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7414))])
-#     h2.run_pyscf()
-#
-#     h2_ferm_ham_1 = h2.hamiltonian("f", h2.oei, h2.tei)
-#     h2_qub_ham_jw_1 = h2.hamiltonian("q", h2.oei, h2.tei, ferm_qubit_map="jw")
-#     h2_qub_ham_bk_1 = h2.hamiltonian("q", h2.oei, h2.tei, ferm_qubit_map="bk")
-#     # h2_mol_ham has format of InteractionOperator
-#     h2_mol_ham = fermionic_hamiltonian(h2.oei, h2.tei, h2.e_nuc)
-#     h2_ferm_ham_2 = openfermion.transforms.get_fermion_operator(h2_mol_ham)
-#     h2_qub_ham_jw_2 = openfermion.jordan_wigner(h2_mol_ham)
-#     h2_qub_ham_bk_2 = openfermion.bravyi_kitaev(h2_mol_ham)
-#
-#     assert h2_ferm_ham_2.isclose(h2_ferm_ham_1)
-#     assert h2_qub_ham_jw_2.isclose(h2_qub_ham_jw_1)
-#     assert h2_qub_ham_bk_2.isclose(h2_qub_ham_bk_1)
-
-
-# @pytest.mark.parametrize(
-#     "mapping,expected_operators",
-#     [
-#         (
-#             None,
-#             [
-#                 ((), -0.04207897647782238),
-#                 (((0, "Z")), 0.17771287465139918),
-#                 (((1, "Z")), 0.1777128746513992),
-#                 (((2, "Z")), -0.24274280513140478),
-#                 (((3, "Z")), -0.24274280513140478),
-#                 (((0, "Z"), (1, "Z")), 0.17059738328801052),
-#                 (((0, "Z"), (2, "Z")), 0.12293305056183809),
-#                 (((0, "Z"), (3, "Z")), 0.16768319457718972),
-#                 (((1, "Z"), (2, "Z")), 0.16768319457718972),
-#                 (((1, "Z"), (3, "Z")), 0.12293305056183809),
-#                 (((2, "Z"), (3, "Z")), 0.17627640804319608),
-#                 (((0, "X"), (1, "X"), (2, "Y"), (3, "Y")), -0.04475014401535165),
-#                 (((0, "X"), (1, "Y"), (2, "Y"), (3, "X")), 0.04475014401535165),
-#                 (((0, "Y"), (1, "X"), (2, "X"), (3, "Y")), 0.04475014401535165),
-#                 (((0, "Y"), (1, "Y"), (2, "X"), (3, "X")), -0.04475014401535165),
-#             ]
-#         ), # H2 JW mapping
-#         (
-#             "bk",
-#             [
-#                 ((), -0.04207897647782244),
-#                 (((0, "Z"),), 0.17771287465139923),
-#                 (((0, "Z"), (1, "Z")), 0.17771287465139918),
-#                 (((2, "Z"),), -0.24274280513140484),
-#                 (((1, "Z"), (2, "Z"), (3, "Z")), -0.24274280513140484),
-#                 (((0, "Y"), (1, "Z"), (2, "Y")), 0.04475014401535165),
-#                 (((0, "X"), (1, "Z"), (2, "X")), 0.04475014401535165),
-#                 (((0, "X"), (1, "Z"), (2, "X"), (3, "Z")), 0.04475014401535165),
-#                 (((0, "Y"), (1, "Z"), (2, "Y"), (3, "Z")), 0.04475014401535165),
-#                 (((1, "Z"),), 0.17059738328801052),
-#                 (((0, "Z"), (2, "Z")), 0.12293305056183809),
-#                 (((0, "Z"), (1, "Z"), (2, "Z")), 0.16768319457718972),
-#                 (((0, "Z"), (1, "Z"), (2, "Z"), (3, "Z")), 0.16768319457718972),
-#                 (((0, "Z"), (2, "Z"), (3, "Z")), 0.12293305056183809),
-#                 (((1, "Z"), (3, "Z")), 0.17627640804319608),
-#             ],
-#         ), # H2 BK mapping
-#     ],
-# )
-# def test_qubit_hamiltonian(mapping, expected_operators):
-#     control = sum(openfermion.QubitOperator(_op, coeff) for _op, coeff in expected_operators)
-#
-#     h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
-#     h2.run_pyscf()
-#
-#     h2_qubit_hamiltonian = h2.hamiltonian("q", ferm_qubit_map=mapping)
-#     assert h2_qubit_hamiltonian.isclose(control)
 
+@pytest.mark.parametrize(
+    "mapping,expected_operators",
+    [
+        (
+            None,
+            [
+                ((), -0.04207897647782238),
+                (((0, "Z")), 0.17771287465139918),
+                (((1, "Z")), 0.1777128746513992),
+                (((2, "Z")), -0.24274280513140478),
+                (((3, "Z")), -0.24274280513140478),
+                (((0, "Z"), (1, "Z")), 0.17059738328801052),
+                (((0, "Z"), (2, "Z")), 0.12293305056183809),
+                (((0, "Z"), (3, "Z")), 0.16768319457718972),
+                (((1, "Z"), (2, "Z")), 0.16768319457718972),
+                (((1, "Z"), (3, "Z")), 0.12293305056183809),
+                (((2, "Z"), (3, "Z")), 0.17627640804319608),
+                (((0, "X"), (1, "X"), (2, "Y"), (3, "Y")), -0.04475014401535165),
+                (((0, "X"), (1, "Y"), (2, "Y"), (3, "X")), 0.04475014401535165),
+                (((0, "Y"), (1, "X"), (2, "X"), (3, "Y")), 0.04475014401535165),
+                (((0, "Y"), (1, "Y"), (2, "X"), (3, "X")), -0.04475014401535165),
+            ],
+        ),  # H2 JW mapping
+        (
+            "bk",
+            [
+                ((), -0.04207897647782244),
+                (((0, "Z"),), 0.17771287465139923),
+                (((0, "Z"), (1, "Z")), 0.17771287465139918),
+                (((2, "Z"),), -0.24274280513140484),
+                (((1, "Z"), (2, "Z"), (3, "Z")), -0.24274280513140484),
+                (((0, "Y"), (1, "Z"), (2, "Y")), 0.04475014401535165),
+                (((0, "X"), (1, "Z"), (2, "X")), 0.04475014401535165),
+                (((0, "X"), (1, "Z"), (2, "X"), (3, "Z")), 0.04475014401535165),
+                (((0, "Y"), (1, "Z"), (2, "Y"), (3, "Z")), 0.04475014401535165),
+                (((1, "Z"),), 0.17059738328801052),
+                (((0, "Z"), (2, "Z")), 0.12293305056183809),
+                (((0, "Z"), (1, "Z"), (2, "Z")), 0.16768319457718972),
+                (((0, "Z"), (1, "Z"), (2, "Z"), (3, "Z")), 0.16768319457718972),
+                (((0, "Z"), (2, "Z"), (3, "Z")), 0.12293305056183809),
+                (((1, "Z"), (3, "Z")), 0.17627640804319608),
+            ],
+        ),  # H2 BK mapping
+    ],
+)
+def test_qubit_hamiltonian(mapping, expected_operators):
+    control = sum(openfermion.QubitOperator(_op, coeff) for _op, coeff in expected_operators)
 
-# def test_symbolic_hamiltonian():
-#     h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
-#     h2.run_pyscf()
-#     hamiltonian = h2.hamiltonian()
-#     ref_sym_ham = SymbolicHamiltonian(
-#         -0.0420789764778224
-#         - 0.0447501440153516 * X(0) * X(1) * Y(2) * Y(3)
-#         + 0.0447501440153516 * X(0) * Y(1) * Y(2) * X(3)
-#         + 0.0447501440153516 * Y(0) * X(1) * X(2) * Y(3)
-#         - 0.0447501440153516 * Y(0) * Y(1) * X(2) * X(3)
-#         + 0.177712874651399 * Z(0)
-#         + 0.170597383288011 * Z(0) * Z(1)
-#         + 0.122933050561838 * Z(0) * Z(2)
-#         + 0.16768319457719 * Z(0) * Z(3)
-#         + 0.177712874651399 * Z(1)
-#         + 0.16768319457719 * Z(1) * Z(2)
-#         + 0.122933050561838 * Z(1) * Z(3)
-#         - 0.242742805131405 * Z(2)
-#         + 0.176276408043196 * Z(2) * Z(3)
-#         - 0.242742805131405 * Z(3)
-#     )
-#     assert np.allclose(hamiltonian.matrix, ref_sym_ham.matrix)
+    h2 = Molecule([("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7))])
+    h2.run_pyscf()
+
+    h2_qubit_hamiltonian = h2.hamiltonian("q", ferm_qubit_map=mapping)
+    assert h2_qubit_hamiltonian.isclose(control)
 
 
 def test_hamiltonian_input_errors():

From 20a086dea72bbfadac655a42123203dbeb3f493d Mon Sep 17 00:00:00 2001
From: Wong Zi Cheng <70616433+chmwzc@users.noreply.github.com>
Date: Wed, 28 Feb 2024 08:01:31 +0000
Subject: [PATCH 09/11] Update workflows

---
 .github/workflows/deploy.yml  | 4 ++--
 .github/workflows/publish.yml | 4 ++--
 .github/workflows/tests.yml   | 2 +-
 3 files changed, 5 insertions(+), 5 deletions(-)

diff --git a/.github/workflows/deploy.yml b/.github/workflows/deploy.yml
index d9a6d6d..d9ae6cb 100644
--- a/.github/workflows/deploy.yml
+++ b/.github/workflows/deploy.yml
@@ -12,11 +12,11 @@ jobs:
     strategy:
       matrix:
         os: [ubuntu-latest, macos-latest]
-        python-version: [3.9]
+        python-version: [3.9, "3.10", "3.11"]
     uses: qiboteam/workflows/.github/workflows/deploy-pip-poetry.yml@main
     with:
       os: ${{ matrix.os }}
       python-version: ${{ matrix.python-version }}
-      publish: ${{ github.event_name == 'release' && github.event.action == 'published' && matrix.os == 'ubuntu-latest' && matrix.python-version == '3.9' }}
+      publish: ${{ github.event_name == 'release' && github.event.action == 'published' && matrix.os == 'ubuntu-latest' && matrix.python-version == "3.10" }}
       poetry-extras: --with tests,docs --all-extras
     secrets: inherit
diff --git a/.github/workflows/publish.yml b/.github/workflows/publish.yml
index 6b3db8f..4ba9548 100644
--- a/.github/workflows/publish.yml
+++ b/.github/workflows/publish.yml
@@ -1,4 +1,4 @@
-name: Deploy docs
+name: docs
 
 on:
   workflow_dispatch:
@@ -29,7 +29,7 @@ jobs:
     needs: [evaluate-label]
     uses: qiboteam/workflows/.github/workflows/deploy-ghpages-latest-stable.yml@main
     with:
-      python-version: 3.9
+      python-version: 3.10
       package-manager: "poetry"
       dependency-path: "**/poetry.lock"
       trigger-label: "${{needs.evaluate-label.outputs.label}}"
diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml
index 8bf989c..b892833 100644
--- a/.github/workflows/tests.yml
+++ b/.github/workflows/tests.yml
@@ -9,7 +9,7 @@ jobs:
     strategy:
       matrix:
         os: [ubuntu-latest, macos-latest]
-        python-version: [3.9]
+        python-version: ["3.9", "3.10", "3.11"]
     uses: qiboteam/workflows/.github/workflows/rules-poetry.yml@main
     with:
       os: ${{ matrix.os }}

From e2effbf4f460c0a81e720efa3ea389eea9cfe534 Mon Sep 17 00:00:00 2001
From: Wong Zi Cheng <70616433+chmwzc@users.noreply.github.com>
Date: Wed, 28 Feb 2024 08:14:06 +0000
Subject: [PATCH 10/11] Revert "Update workflows"

This reverts commit 20a086dea72bbfadac655a42123203dbeb3f493d.
---
 .github/workflows/deploy.yml  | 4 ++--
 .github/workflows/publish.yml | 4 ++--
 .github/workflows/tests.yml   | 2 +-
 3 files changed, 5 insertions(+), 5 deletions(-)

diff --git a/.github/workflows/deploy.yml b/.github/workflows/deploy.yml
index d9ae6cb..d9a6d6d 100644
--- a/.github/workflows/deploy.yml
+++ b/.github/workflows/deploy.yml
@@ -12,11 +12,11 @@ jobs:
     strategy:
       matrix:
         os: [ubuntu-latest, macos-latest]
-        python-version: [3.9, "3.10", "3.11"]
+        python-version: [3.9]
     uses: qiboteam/workflows/.github/workflows/deploy-pip-poetry.yml@main
     with:
       os: ${{ matrix.os }}
       python-version: ${{ matrix.python-version }}
-      publish: ${{ github.event_name == 'release' && github.event.action == 'published' && matrix.os == 'ubuntu-latest' && matrix.python-version == "3.10" }}
+      publish: ${{ github.event_name == 'release' && github.event.action == 'published' && matrix.os == 'ubuntu-latest' && matrix.python-version == '3.9' }}
       poetry-extras: --with tests,docs --all-extras
     secrets: inherit
diff --git a/.github/workflows/publish.yml b/.github/workflows/publish.yml
index 4ba9548..6b3db8f 100644
--- a/.github/workflows/publish.yml
+++ b/.github/workflows/publish.yml
@@ -1,4 +1,4 @@
-name: docs
+name: Deploy docs
 
 on:
   workflow_dispatch:
@@ -29,7 +29,7 @@ jobs:
     needs: [evaluate-label]
     uses: qiboteam/workflows/.github/workflows/deploy-ghpages-latest-stable.yml@main
     with:
-      python-version: 3.10
+      python-version: 3.9
       package-manager: "poetry"
       dependency-path: "**/poetry.lock"
       trigger-label: "${{needs.evaluate-label.outputs.label}}"
diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml
index b892833..8bf989c 100644
--- a/.github/workflows/tests.yml
+++ b/.github/workflows/tests.yml
@@ -9,7 +9,7 @@ jobs:
     strategy:
       matrix:
         os: [ubuntu-latest, macos-latest]
-        python-version: ["3.9", "3.10", "3.11"]
+        python-version: [3.9]
     uses: qiboteam/workflows/.github/workflows/rules-poetry.yml@main
     with:
       os: ${{ matrix.os }}

From 3ebf614a698ffa1ecc900976a03efb60e8409bd3 Mon Sep 17 00:00:00 2001
From: "pre-commit-ci[bot]"
 <66853113+pre-commit-ci[bot]@users.noreply.github.com>
Date: Mon, 18 Mar 2024 22:29:52 +0000
Subject: [PATCH 11/11] [pre-commit.ci] pre-commit autoupdate
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

updates:
- [github.com/psf/black: 24.2.0 → 24.3.0](https://github.com/psf/black/compare/24.2.0...24.3.0)
---
 .pre-commit-config.yaml | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml
index b19996b..604b439 100644
--- a/.pre-commit-config.yaml
+++ b/.pre-commit-config.yaml
@@ -9,7 +9,7 @@ repos:
       - id: check-yaml
       - id: debug-statements
   - repo: https://github.com/psf/black
-    rev: 24.2.0
+    rev: 24.3.0
     hooks:
       - id: black
         args: