Merge pull request #85 from Zhaoyilunnn/master

qaoa demo

quafu/algorithms/__init__.py

@@ -1,3 +1,5 @@
 """Algorithm module"""
 
-from .hamiltonian import Hamiltonian
+from .hamiltonian import *
+from .ansatz import *
+from .estimator import *

quafu/algorithms/ansatz.py

@@ -13,28 +13,70 @@
 # limitations under the License.
 
 """Ansatz circuits for VQA"""
+from typing import List
+import numpy as np
 
+from quafu.algorithms.hamiltonian import Hamiltonian
 from quafu.circuits.quantum_circuit import QuantumCircuit
 from quafu.synthesis.evolution import ProductFormula
 
 
 class QAOACircuit(QuantumCircuit):
     """QAOA circuit"""
 
-    def __init__(self, pauli: str, num_layers: int = 1):
+    def __init__(self, hamiltonian: Hamiltonian, num_layers: int = 1):
         """Instantiate a QAOAAnsatz"""
-        num_qubits = len(pauli)
-        super().__init__(num_qubits)
+        self._pauli_list = hamiltonian.pauli_list
+        self._coeffs = hamiltonian.coeffs
         self._num_layers = num_layers
         self._evol = ProductFormula()
-        self._build(pauli)
 
-    def _build(self, pauli):
+        # Initialize parameters
+        self._beta = np.zeros(num_layers)
+        self._gamma = np.zeros(num_layers)
+
+        # Build circuit structure
+        num_qubits = len(self._pauli_list[0])
+        super().__init__(num_qubits)
+        self._build()
+
+    @property
+    def parameters(self):
+        """Return complete parameters of the circuit"""
+        paras_list = []
+        for layer in range(self._num_layers):
+            for _ in range(len(self._pauli_list)):
+                paras_list.append(self._gamma[layer])
+            for _ in range(self.num):
+                paras_list.append(self._beta[layer])
+        return paras_list
+
+    def _add_superposition(self):
+        """Apply H gate on all qubits"""
+        for i in range(self.num):
+            self.h(i)
+
+    def _build(self):
         """Construct circuit"""
-        gate_list = self._evol.evol(pauli, 0.0)
-        for g in gate_list:
-            self.add_gate(g)
 
-    # def get_expectations(self):
-    #     """Calculate the expectations of an operator"""
-    #     pass
+        self._add_superposition()
+
+        for layer in range(self._num_layers):
+            # Add H_D layer
+            for pauli_str in self._pauli_list:
+                gate_list = self._evol.evol(pauli_str, self._gamma[layer])
+                for g in gate_list:
+                    self.add_gate(g)
+
+            # Add H_B layer
+            for i in range(self.num):
+                self.rx(i, self._beta[layer])
+
+    def update_params(self, beta: List[float], gamma: List[float]):
+        """Update parameters of QAOA circuit"""
+        # First build parameter list
+        assert len(beta) == self._num_layers and len(gamma) == self._num_layers
+        num_para_gates = len(self.parameterized_gates)
+        assert num_para_gates % self._num_layers == 0
+        self._beta, self._gamma = beta, gamma
+        super().update_params(self.parameters)

quafu/algorithms/estimator.py

@@ -12,11 +12,13 @@
 # See the License for the specific language governing permissions and
 # limitations under the License.
 """Pre-build wrapper to calculate expectation value"""
+import numpy as np
 
 from typing import Optional
 from quafu import QuantumCircuit
 from quafu.simulators.simulator import simulate
 from quafu.tasks.tasks import Task
+from quafu.algorithms.hamiltonian import Hamiltonian
 
 
 class Estimator:
@@ -47,29 +49,39 @@ def __init__(
             self._task.config(backend=self._backend)
             self._task.config(**task_options)
 
-    def _run_real_machine(self, observables):
+    def _run_real_machine(self, observables: Hamiltonian):
         """Submit to quafu service"""
         if not isinstance(self._task, Task):
             raise ValueError("task not set")
-        res, obsexp = self._task.submit(self._circ, observables)
-        return res, obsexp
+        obs = observables.to_legacy_quafu_pauli_list()
+        _, obsexp = self._task.submit(self._circ, obs)
+        return sum(obsexp)
 
-    def _run_simulation(self, observables):
-        """TODO"""
-        sim_res = simulate(self._circ)
-        # TODO
-        # sim_res.calculate_obs()
-        return None, None
+    def _run_simulation(self, observables: Hamiltonian):
+        """Run using quafu simulator"""
+        sim_state = simulate(self._circ, output="state_vector").get_statevector()
+        expectation = np.matmul(
+            np.matmul(sim_state.conj().T, observables.get_matrix()), sim_state
+        ).real
+        return expectation
 
-    def run(self, observables, paras_list=None):
+    def run(self, observables, *params):
         """Calculate estimation for given observables
 
         Args:
             observables: observables to be estimated.
             paras_list: list of parameters of self.circ.
+
+        Returns:
+            Expectation value
         """
-        if paras_list is not None:
-            self._circ.update_params(paras_list)
+        if observables.num_qubits != self._circ.num:
+            raise ValueError(
+                "The number of qubits in the observables does not match the circuit"
+            )
+
+        if params[0] is not None:
+            self._circ.update_params(*params)
 
         if self._backend == "sim":
             return self._run_simulation(observables)

quafu/algorithms/hamiltonian.py

@@ -19,14 +19,10 @@
 from collections.abc import Iterable
 import numpy as np
 from quafu.exceptions import QuafuError
+from quafu.elements.element_gates.matrices import XMatrix, YMatrix, ZMatrix, IdMatrix
 
 
-PAULI_MAT = {
-    "I": np.array([[1, 0], [0, 1]], dtype=complex),
-    "X": np.array([[0, 1], [1, 0]], dtype=complex),
-    "Y": np.array([[0, -1j], [1j, 0]], dtype=complex),
-    "Z": np.array([[1, 0], [0, -1]], dtype=complex)
-}
+PAULI_MAT = {"I": IdMatrix, "X": XMatrix, "Y": YMatrix, "Z": ZMatrix}
 
 
 class Hamiltonian:
@@ -41,6 +37,21 @@ def __init__(self, pauli_str_list: list[str], coeffs: np.ndarray) -> None:
         self._pauli_str_list = pauli_str_list
         self._coeffs = coeffs
 
+    @property
+    def num_qubits(self):
+        """Get the number of qubits of the system"""
+        return len(self.pauli_list[0])
+
+    @property
+    def pauli_list(self):
+        """Get pauli string list"""
+        return self._pauli_str_list
+
+    @property
+    def coeffs(self):
+        """Get coefficients of each pauli string"""
+        return self._coeffs
+
     @staticmethod
     def from_pauli_list(pauli_list: Iterable[tuple[str, complex]]) -> Hamiltonian:
         """
@@ -65,6 +76,17 @@ def from_pauli_list(pauli_list: Iterable[tuple[str, complex]]) -> Hamiltonian:
 
         return Hamiltonian(pauli_str_list, coeffs)
 
+    # TODO(zhaoyilun): delete this in the future
+    def to_legacy_quafu_pauli_list(self):
+        """Transform to legacy quafu pauli list format,
+        this is a temperal function and should be deleted later"""
+        res = []
+        for pauli in self._pauli_str_list:
+            for i, p in enumerate(pauli[::-1]):
+                if p in ["X", "Y", "Z"]:
+                    res.append([p, [i]])
+        return res
+
     def _get_pauli_mat(self, pauli_str: str):
         """Calculate the matrix of a pauli string"""
         mat = None

quafu/synthesis/evolution.py

@@ -22,12 +22,13 @@
 def two_qubit_evol(pauli: str, time: float, cx_structure: str = "chain"):
     """
     Args:
-        pauli: Pauli string
+        pauli: Pauli string (little endian convention)
         time: Evolution time
         cx_structure: TODO
     """
-    qubits = [i for i in range(len(pauli)) if pauli[i] != "I"]
-    labels = np.array([pauli[i] for i in qubits])
+    reversed_pauli = pauli[::-1]
+    qubits = [i for i in range(len(reversed_pauli)) if reversed_pauli[i] != "I"]
+    labels = np.array([reversed_pauli[i] for i in qubits])
     gates = []
 
     if all(labels == "X"):

tests/quafu/algorithms/ansatz_test.py

@@ -1,7 +1,29 @@
+# (C) Copyright 2023 Beijing Academy of Quantum Information Sciences
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+#    http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+import numpy as np
 from quafu.algorithms.ansatz import QAOACircuit
+from quafu.algorithms.hamiltonian import Hamiltonian
 
 
 class TestQAOACircuit:
+    TEST_HAM = Hamiltonian(["IIZZ", "ZZII", "IZZI", "ZIIZ"], np.array([1, 1, 1, 1]))
+
     def test_build(self):
-        qaoa = QAOACircuit("IIZZ")
+        qaoa = QAOACircuit(self.TEST_HAM)
+        print("\n ::: testing ::: \n")
         qaoa.draw_circuit()
+
+    def test_update_params(self):
+        pass

tests/quafu/algorithms/estimator_test.py

@@ -12,9 +12,14 @@
 # See the License for the specific language governing permissions and
 # limitations under the License.
 
+import sys
+import numpy as np
+import pytest
+import math
 from unittest.mock import patch
 from quafu import ExecResult
 from quafu.algorithms.estimator import Estimator
+from quafu.algorithms.hamiltonian import Hamiltonian
 
 from quafu.circuits.quantum_circuit import QuantumCircuit
 from quafu.tasks.tasks import Task
@@ -34,27 +39,42 @@
 class TestEstimator:
     """Test class of Estimator"""
 
-    @patch("quafu.users.userapi.User._load_account_token", autospec=True)
-    @patch("quafu.users.userapi.User.get_available_backends", autospec=True)
-    @patch("quafu.tasks.tasks.Task.send", autospec=True)
-    def test_run(self, mock_send, mock_backends, mock_load_account):
-        """Test Estimator.run"""
-        mock_send.return_value = TEST_EXE_RES
-        mock_backends.return_value = {"ScQ-P10": None}
+    def build_circuit(self):
+        """Build a random circuit for testing purpose"""
         circ = QuantumCircuit(5)
 
         for i in range(5):
             if i % 2 == 0:
                 circ.h(i)
 
+        circ.cnot(0, 1)
         circ.draw_circuit()
         measures = list(range(5))
         circ.measure(measures)
-        test_ising = [["X", [i]] for i in range(5)]
-        test_ising.extend([["ZZ", [i, i + 1]] for i in range(4)])
+        test_ising = Hamiltonian(
+            ["IIIZZ", "ZZIII", "IZZII", "ZIIIZ"], np.array([1, 1, 1, 1])
+        )
+        return circ, test_ising
 
+    @patch("quafu.users.userapi.User._load_account_token", autospec=True)
+    @patch("quafu.users.userapi.User.get_available_backends", autospec=True)
+    @patch("quafu.tasks.tasks.Task.send", autospec=True)
+    def test_run(self, mock_send, mock_backends, mock_load_account):
+        """Test Estimator.run"""
+        mock_send.return_value = TEST_EXE_RES
+        mock_backends.return_value = {"ScQ-P10": None}
+        circ, test_ising = self.build_circuit()
         estimator = Estimator(circ, backend="ScQ-P10")
-        res, obsexp = estimator.run(test_ising)
+        expectation = estimator.run(test_ising, None)
         task = Task()
-        res_org, obsexp_org = task.submit(circ, test_ising)
-        assert res == res_org and obsexp == obsexp_org
+        res_org, obsexp_org = task.submit(circ, test_ising.to_legacy_quafu_pauli_list())
+        assert expectation == sum(obsexp_org)
+
+    @pytest.mark.skipif(
+        sys.platform == "darwin", reason="Avoid error on MacOS arm arch."
+    )
+    def test_run_sim(self):
+        circ, test_ising = self.build_circuit()
+        estimator = Estimator(circ)
+        expectation = estimator.run(test_ising, None)
+        assert math.isclose(expectation, 1.0)