Initialize

tests/config.py

@@ -0,0 +1,12 @@
+from dataclasses import dataclass
+from typing import Optional
+
+
+@dataclass
+class Executor:
+    backend: str
+    platform: Optional[str] = None
+
+
+qibo = Executor(backend="qibojit", platform="numpy")
+quimb = Executor(backend="numpy")

tests/conftest.py

@@ -0,0 +1,66 @@
+"""conftest.py.
+
+Pytest fixtures.
+"""
+
+import sys
+
+import pytest
+
+# backends to be tested
+# TODO: add cutensornet and quimb here as well
+BACKENDS = ["cutensornet"]
+# BACKENDS = ["qmatchatea"]
+
+
+def get_backend(backend_name):
+
+    from qibotn.backends.cutensornet import CuTensorNet
+    from qibotn.backends.qmatchatea import QMatchaTeaBackend
+
+    NAME2BACKEND = {"qmatchatea": QMatchaTeaBackend, "cutensornet": CuTensorNet}
+
+    return NAME2BACKEND[backend_name]()
+
+
+AVAILABLE_BACKENDS = []
+for backend_name in BACKENDS:
+    try:
+        _backend = get_backend(backend_name)
+        AVAILABLE_BACKENDS.append(backend_name)
+    except (ModuleNotFoundError, ImportError):
+        pass
+
+
+def pytest_runtest_setup(item):
+    ALL = {"darwin", "linux"}
+    supported_platforms = ALL.intersection(mark.name for mark in item.iter_markers())
+    plat = sys.platform
+    if supported_platforms and plat not in supported_platforms:  # pragma: no cover
+        # case not covered by workflows
+        pytest.skip(f"Cannot run test on platform {plat}.")
+
+
+@pytest.fixture
+def backend(backend_name):
+    yield get_backend(backend_name)
+
+
+def pytest_runtest_setup(item):
+    ALL = {"darwin", "linux"}
+    supported_platforms = ALL.intersection(mark.name for mark in item.iter_markers())
+    plat = sys.platform
+    if supported_platforms and plat not in supported_platforms:  # pragma: no cover
+        # case not covered by workflows
+        pytest.skip(f"Cannot run test on platform {plat}.")
+
+
+def pytest_configure(config):
+    config.addinivalue_line("markers", "linux: mark test to run only on linux")
+
+
+def pytest_generate_tests(metafunc):
+    module_name = metafunc.module.__name__
+
+    if "backend_name" in metafunc.fixturenames:
+        metafunc.parametrize("backend_name", AVAILABLE_BACKENDS)

tests/test_circuit_execution.py

@@ -0,0 +1,91 @@
+import math
+
+import pytest
+from qibo import Circuit, gates, hamiltonians
+from qibo.symbols import X, Z
+
+from qibotn.backends.qmatchatea import QMatchaTeaBackend
+
+
+def build_observable(nqubits):
+    """Helper function to construct a target observable."""
+    hamiltonian_form = 0
+    for i in range(nqubits):
+        hamiltonian_form += 0.5 * X(i % nqubits) * Z((i + 1) % nqubits)
+
+    hamiltonian = hamiltonians.SymbolicHamiltonian(form=hamiltonian_form)
+    return hamiltonian, hamiltonian_form
+
+
+def build_GHZ(nqubits):
+    """Helper function to construct a layered quantum circuit."""
+    circ = Circuit(nqubits)
+    circ.add(gates.H(0))
+    [circ.add(gates.CNOT(q, q + 1)) for q in range(nqubits - 1)]
+    return circ
+
+
+def construct_targets(nqubits):
+    """Construct strings of 1s and 0s of size `nqubits`."""
+    ones = "1" * nqubits
+    zeros = "0" * nqubits
+    return ones, zeros
+
+
+@pytest.mark.parametrize("nqubits", [2, 10, 40])
+def test_probabilities(backend, nqubits):
+
+    circ = build_GHZ(nqubits=nqubits)
+
+    if isinstance(backend, QMatchaTeaBackend):
+        # unbiased prob
+        out_u = backend.execute_circuit(
+            circuit=circ,
+            prob_type="U",
+            num_samples=1000,
+        ).probabilities()
+
+        math.isclose(out_u[0], 0.5, abs_tol=1e-7)
+        math.isclose(out_u[1], 0.5, abs_tol=1e-7)
+
+        out_g = backend.execute_circuit(
+            circuit=circ,
+            prob_type="G",
+            prob_threshold=1.0,
+        ).probabilities()
+
+        math.isclose(out_g[0], 0.5, abs_tol=1e-7)
+        math.isclose(out_g[1], 0.5, abs_tol=1e-7)
+
+        out_e = backend.execute_circuit(
+            circuit=circ,
+            prob_type="E",
+            prob_threshold=0.2,
+        ).probabilities()
+
+        math.isclose(out_e[0], 0.5, abs_tol=1e-7)
+        math.isclose(out_e[1], 0.5, abs_tol=1e-7)
+
+
+@pytest.mark.parametrize("nqubits", [2, 10, 40])
+@pytest.mark.parametrize("nshots", [100, 1000])
+def test_shots(backend, nqubits, nshots):
+    circ = build_GHZ(nqubits=nqubits)
+    ones, zeros = construct_targets(nqubits)
+
+    # For p = 0.5, sigma = sqrt(nshots * 0.5 * 0.5) = sqrt(nshots)/2.
+    sigma_threshold = 3 * (math.sqrt(nshots) / 2)
+
+    outcome = backend.execute_circuit(circ, nshots=nshots)
+    frequencies = outcome.frequencies()
+
+    shots_ones = frequencies.get(ones, 0)
+    shots_zeros = frequencies.get(zeros, 0)
+
+    # Check that the counts for both outcomes are within the 3-sigma threshold of nshots/2.
+    assert (
+        abs(shots_ones - (nshots / 2)) < sigma_threshold
+    ), f"Count for {ones} deviates too much: {shots_ones} vs expected {nshots/2}"
+    assert (
+        abs(shots_zeros - (nshots / 2)) < sigma_threshold
+    ), f"Count for {zeros} deviates too much: {shots_zeros} vs expected {nshots/2}"

tests/test_cuquantum_cutensor_backend.py

@@ -0,0 +1,199 @@
+import math
+
+import cupy as cp
+import pytest
+import qibo
+from qibo import construct_backend, hamiltonians
+from qibo.models import QFT
+from qibo.symbols import X, Z
+
+ABS_TOL = 1e-7
+
+
+def qibo_qft(nqubits, swaps):
+    circ_qibo = QFT(nqubits, swaps)
+    state_vec = circ_qibo().state(numpy=True)
+    return circ_qibo, state_vec
+
+
+def build_observable(nqubits):
+    """Helper function to construct a target observable."""
+    hamiltonian_form = 0
+    for i in range(nqubits):
+        hamiltonian_form += 0.5 * X(i % nqubits) * Z((i + 1) % nqubits)
+
+    hamiltonian = hamiltonians.SymbolicHamiltonian(form=hamiltonian_form)
+    return hamiltonian, hamiltonian_form
+
+
+def build_observable_dict(nqubits):
+    """Construct a target observable as a dictionary representation.
+
+    Returns a dictionary suitable for `create_hamiltonian_from_dict`.
+    """
+    terms = []
+
+    for i in range(nqubits):
+        term = {
+            "coefficient": 0.5,
+            "operators": [("X", i % nqubits), ("Z", (i + 1) % nqubits)],
+        }
+        terms.append(term)
+
+    return {"terms": terms}
+
+
+@pytest.mark.gpu
+@pytest.mark.parametrize("nqubits", [1, 2, 5, 10])
+def test_eval(nqubits: int, dtype="complex128"):
+    """
+    Args:
+        nqubits (int): Total number of qubits in the system.
+        dtype (str): The data type for precision, 'complex64' for single,
+            'complex128' for double.
+    """
+    # Test qibo
+    qibo.set_backend(backend="numpy")
+    qibo_circ, result_sv = qibo_qft(nqubits, swaps=True)
+    result_sv_cp = cp.asarray(result_sv)
+
+    # Test cutensornet
+    backend = construct_backend(backend="qibotn", platform="cutensornet")
+    # Test with no settings specified. Default is dense vector calculation without MPI or NCCL.
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    print(
+        f"State vector difference: {abs(result_tn.statevector.flatten() - result_sv_cp).max():0.3e}"
+    )
+    assert cp.allclose(
+        result_sv_cp, result_tn.statevector.flatten()
+    ), "Resulting dense vectors do not match"
+
+    # Test with explicit settings specified.
+    comp_set_w_bool = {
+        "MPI_enabled": False,
+        "MPS_enabled": False,
+        "NCCL_enabled": False,
+        "expectation_enabled": False,
+    }
+    backend.configure_tn_simulation(comp_set_w_bool)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    print(
+        f"State vector difference: {abs(result_tn.statevector.flatten() - result_sv_cp).max():0.3e}"
+    )
+    assert cp.allclose(
+        result_sv_cp, result_tn.statevector.flatten()
+    ), "Resulting dense vectors do not match"
+
+
+@pytest.mark.gpu
+@pytest.mark.parametrize("nqubits", [2, 5, 10])
+def test_mps(nqubits: int, dtype="complex128"):
+    """Evaluate MPS with cuQuantum.
+
+    Args:
+        nqubits (int): Total number of qubits in the system.
+        dtype (str): The data type for precision, 'complex64' for single,
+            'complex128' for double.
+    """
+
+    # Test qibo
+    qibo.set_backend(backend="numpy")
+    qibo_circ, result_sv = qibo_qft(nqubits, swaps=True)
+    result_sv_cp = cp.asarray(result_sv)
+
+    # Test cutensornet
+    backend = construct_backend(backend="qibotn", platform="cutensornet")
+    # Test with simple MPS settings specified using bool. Uses the default MPS parameters.
+    comp_set_w_bool = {
+        "MPI_enabled": False,
+        "MPS_enabled": True,
+        "NCCL_enabled": False,
+        "expectation_enabled": False,
+    }
+    backend.configure_tn_simulation(comp_set_w_bool)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    print(
+        f"State vector difference: {abs(result_tn.statevector.flatten() - result_sv_cp).max():0.3e}"
+    )
+    assert cp.allclose(
+        result_tn.statevector.flatten(), result_sv_cp
+    ), "Resulting dense vectors do not match"
+
+    # Test with explicit MPS computation settings specified using Dict. Users able to specify parameters like qr_method etc.
+    comp_set_w_MPS_config_para = {
+        "MPI_enabled": False,
+        "MPS_enabled": {
+            "qr_method": False,
+            "svd_method": {
+                "partition": "UV",
+                "abs_cutoff": 1e-12,
+            },
+        },
+        "NCCL_enabled": False,
+        "expectation_enabled": False,
+    }
+    backend.configure_tn_simulation(comp_set_w_MPS_config_para)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    print(
+        f"State vector difference: {abs(result_tn.statevector.flatten() - result_sv_cp).max():0.3e}"
+    )
+    assert cp.allclose(
+        result_tn.statevector.flatten(), result_sv_cp
+    ), "Resulting dense vectors do not match"
+
+
+@pytest.mark.parametrize("nqubits", [2, 5, 10])
+def test_expectation(nqubits: int, dtype="complex128"):
+
+    # Test qibo
+    qibo_circ, state_vec_qibo = qibo_qft(nqubits, swaps=True)
+    ham, ham_form = build_observable(nqubits)
+    numpy_backend = construct_backend("numpy")
+    exact_expval = numpy_backend.calculate_expectation_state(
+        hamiltonian=ham,
+        state=state_vec_qibo,
+        normalize=False,
+    )
+
+    # Test cutensornet
+    backend = construct_backend(backend="qibotn", platform="cutensornet")
+
+    # Test with simple settings using bool. Uses default Hamilitonian for expectation calculation.
+    comp_set_w_bool = {
+        "MPI_enabled": False,
+        "MPS_enabled": False,
+        "NCCL_enabled": False,
+        "expectation_enabled": True,
+    }
+    backend.configure_tn_simulation(comp_set_w_bool)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    assert math.isclose(
+        exact_expval.item(), result_tn.real.get().item(), abs_tol=ABS_TOL
+    )
+
+    # Test with user defined hamiltonian using "hamiltonians.SymbolicHamiltonian" object.
+    comp_set_w_hamiltonian_obj = {
+        "MPI_enabled": False,
+        "MPS_enabled": False,
+        "NCCL_enabled": False,
+        "expectation_enabled": ham,
+    }
+    backend.configure_tn_simulation(comp_set_w_hamiltonian_obj)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    assert math.isclose(
+        exact_expval.item(), result_tn.real.get().item(), abs_tol=ABS_TOL
+    )
+
+    # Test with user defined hamiltonian using Dictionary object form of hamiltonian.
+    ham_dict = build_observable_dict(nqubits)
+    comp_set_w_hamiltonian_dict = {
+        "MPI_enabled": False,
+        "MPS_enabled": False,
+        "NCCL_enabled": False,
+        "expectation_enabled": ham_dict,
+    }
+    backend.configure_tn_simulation(comp_set_w_hamiltonian_dict)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    assert math.isclose(
+        exact_expval.item(), result_tn.real.get().item(), abs_tol=ABS_TOL
+    )

tests/test_cuquantum_cutensor_mpi_backend.py

@@ -0,0 +1,315 @@
+# mpirun --allow-run-as-root -np 2 python -m pytest --with-mpi test_cuquantum_cutensor_mpi_backend.py
+
+import math
+
+import cupy as cp
+import numpy as np
+import pytest
+import qibo
+from qibo import construct_backend, hamiltonians
+from qibo.models import QFT
+from qibo.symbols import X, Z
+
+ABS_TOL = 1e-7
+
+
+def qibo_qft(nqubits, swaps):
+    circ_qibo = QFT(nqubits, swaps)
+    state_vec = circ_qibo().state(numpy=True)
+    return circ_qibo, state_vec
+
+
+def build_observable(nqubits):
+    """Helper function to construct a target observable."""
+    hamiltonian_form = 0
+    for i in range(nqubits):
+        hamiltonian_form += 0.5 * X(i % nqubits) * Z((i + 1) % nqubits)
+
+    hamiltonian = hamiltonians.SymbolicHamiltonian(form=hamiltonian_form)
+    return hamiltonian, hamiltonian_form
+
+
+def build_observable_dict(nqubits):
+    """Construct a target observable as a dictionary representation.
+
+    Returns a dictionary suitable for `create_hamiltonian_from_dict`.
+    """
+    terms = []
+
+    for i in range(nqubits):
+        term = {
+            "coefficient": 0.5,
+            "operators": [("X", i % nqubits), ("Z", (i + 1) % nqubits)],
+        }
+        terms.append(term)
+
+    return {"terms": terms}
+
+
+@pytest.mark.gpu
+@pytest.mark.mpi
+@pytest.mark.parametrize("nqubits", [1, 2, 5, 7, 10])
+def test_eval_mpi(nqubits: int, dtype="complex128"):
+    """
+    Args:
+        nqubits (int): Total number of qubits in the system.
+        dtype (str): The data type for precision, 'complex64' for single,
+            'complex128' for double.
+    """
+    # Test qibo
+    qibo.set_backend(backend="numpy")
+    qibo_circ, result_sv = qibo_qft(nqubits, swaps=True)
+    result_sv_cp = cp.asarray(result_sv)
+
+    # Test cutensornet
+    backend = construct_backend(backend="qibotn", platform="cutensornet")
+
+    # Test with explicit settings specified.
+    comp_set_w_bool = {
+        "MPI_enabled": True,
+        "MPS_enabled": False,
+        "NCCL_enabled": False,
+        "expectation_enabled": False,
+    }
+    backend.configure_tn_simulation(comp_set_w_bool)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    result_tn_cp = cp.asarray(result_tn.statevector.flatten())
+
+    print(f"State vector difference: {abs(result_tn_cp - result_sv_cp).max():0.3e}")
+
+    if backend.rank == 0:
+
+        assert cp.allclose(
+            result_sv_cp, result_tn_cp
+        ), "Resulting dense vectors do not match"
+    else:
+        assert (
+            isinstance(result_tn_cp, cp.ndarray)
+            and result_tn_cp.size == 1
+            and result_tn_cp.item() == 0
+        ), f"Rank {backend.rank}: result_tn_cp should be scalar/array with 0, got {result_tn_cp}"
+
+
+@pytest.mark.gpu
+@pytest.mark.mpi
+@pytest.mark.parametrize("nqubits", [1, 2, 5, 7, 10])
+def test_expectation_mpi(nqubits: int, dtype="complex128"):
+
+    # Test qibo
+    qibo_circ, state_vec_qibo = qibo_qft(nqubits, swaps=True)
+    ham, ham_form = build_observable(nqubits)
+    numpy_backend = construct_backend("numpy")
+    exact_expval = numpy_backend.calculate_expectation_state(
+        hamiltonian=ham,
+        state=state_vec_qibo,
+        normalize=False,
+    )
+
+    # Test cutensornet
+    backend = construct_backend(backend="qibotn", platform="cutensornet")
+
+    # Test with simple settings using bool. Uses default Hamilitonian for expectation calculation.
+    comp_set_w_bool = {
+        "MPI_enabled": True,
+        "MPS_enabled": False,
+        "NCCL_enabled": False,
+        "expectation_enabled": True,
+    }
+    backend.configure_tn_simulation(comp_set_w_bool)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    if backend.rank == 0:
+        # Compare numerical values
+        assert math.isclose(
+            exact_expval.item(), float(result_tn[0]), abs_tol=ABS_TOL
+        ), f"Rank {backend.rank}: mismatch, expected {exact_expval}, got {result_tn}"
+
+    else:
+        # Rank > 0: must be hardcoded [0] (int)
+        assert (
+            isinstance(result_tn, (np.ndarray, cp.ndarray))
+            and result_tn.size == 1
+            and np.issubdtype(result_tn.dtype, np.integer)
+            and result_tn.item() == 0
+        ), f"Rank {backend.rank}: expected int array [0], got {result_tn}"
+
+    # Test with user defined hamiltonian using "hamiltonians.SymbolicHamiltonian" object.
+    comp_set_w_hamiltonian_obj = {
+        "MPI_enabled": True,
+        "MPS_enabled": False,
+        "NCCL_enabled": False,
+        "expectation_enabled": ham,
+    }
+    backend.configure_tn_simulation(comp_set_w_hamiltonian_obj)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    if backend.rank == 0:
+        # Compare numerical values
+        assert math.isclose(
+            exact_expval.item(), float(result_tn[0]), abs_tol=ABS_TOL
+        ), f"Rank {backend.rank}: mismatch, expected {exact_expval}, got {result_tn}"
+
+    else:
+        # Rank > 0: must be hardcoded [0] (int)
+        assert (
+            isinstance(result_tn, (np.ndarray, cp.ndarray))
+            and result_tn.size == 1
+            and np.issubdtype(result_tn.dtype, np.integer)
+            and result_tn.item() == 0
+        ), f"Rank {backend.rank}: expected int array [0], got {result_tn}"
+
+    # Test with user defined hamiltonian using Dictionary object form of hamiltonian.
+    ham_dict = build_observable_dict(nqubits)
+    comp_set_w_hamiltonian_dict = {
+        "MPI_enabled": True,
+        "MPS_enabled": False,
+        "NCCL_enabled": False,
+        "expectation_enabled": ham_dict,
+    }
+    backend.configure_tn_simulation(comp_set_w_hamiltonian_dict)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    if backend.rank == 0:
+        # Compare numerical values
+        assert math.isclose(
+            exact_expval.item(), float(result_tn[0]), abs_tol=ABS_TOL
+        ), f"Rank {backend.rank}: mismatch, expected {exact_expval}, got {result_tn}"
+
+    else:
+        # Rank > 0: must be hardcoded [0] (int)
+        assert (
+            isinstance(result_tn, (np.ndarray, cp.ndarray))
+            and result_tn.size == 1
+            and np.issubdtype(result_tn.dtype, np.integer)
+            and result_tn.item() == 0
+        ), f"Rank {backend.rank}: expected int array [0], got {result_tn}"
+
+
+@pytest.mark.gpu
+@pytest.mark.mpi
+@pytest.mark.parametrize("nqubits", [1, 2, 5, 7, 10])
+def test_eval_nccl(nqubits: int, dtype="complex128"):
+    """
+    Args:
+        nqubits (int): Total number of qubits in the system.
+        dtype (str): The data type for precision, 'complex64' for single,
+            'complex128' for double.
+    """
+    # Test qibo
+    qibo.set_backend(backend="numpy")
+    qibo_circ, result_sv = qibo_qft(nqubits, swaps=True)
+    result_sv_cp = cp.asarray(result_sv)
+
+    # Test cutensornet
+    backend = construct_backend(backend="qibotn", platform="cutensornet")
+
+    # Test with explicit settings specified.
+    comp_set_w_bool = {
+        "MPI_enabled": False,
+        "MPS_enabled": False,
+        "NCCL_enabled": True,
+        "expectation_enabled": False,
+    }
+    backend.configure_tn_simulation(comp_set_w_bool)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    result_tn_cp = cp.asarray(result_tn.statevector.flatten())
+
+    if backend.rank == 0:
+        assert cp.allclose(
+            result_sv_cp, result_tn_cp
+        ), "Resulting dense vectors do not match"
+    else:
+        assert (
+            isinstance(result_tn_cp, cp.ndarray)
+            and result_tn_cp.size == 1
+            and result_tn_cp.item() == 0
+        ), f"Rank {backend.rank}: result_tn_cp should be scalar/array with 0, got {result_tn_cp}"
+
+
+@pytest.mark.gpu
+@pytest.mark.mpi
+@pytest.mark.parametrize("nqubits", [1, 2, 5, 7, 10])
+def test_expectation_NCCL(nqubits: int, dtype="complex128"):
+
+    # Test qibo
+    qibo_circ, state_vec_qibo = qibo_qft(nqubits, swaps=True)
+    ham, ham_form = build_observable(nqubits)
+    numpy_backend = construct_backend("numpy")
+    exact_expval = numpy_backend.calculate_expectation_state(
+        hamiltonian=ham,
+        state=state_vec_qibo,
+        normalize=False,
+    )
+
+    # Test cutensornet
+    backend = construct_backend(backend="qibotn", platform="cutensornet")
+
+    # Test with simple settings using bool. Uses default Hamilitonian for expectation calculation.
+    comp_set_w_bool = {
+        "MPI_enabled": False,
+        "MPS_enabled": False,
+        "NCCL_enabled": True,
+        "expectation_enabled": True,
+    }
+    backend.configure_tn_simulation(comp_set_w_bool)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    if backend.rank == 0:
+        # Compare numerical values
+        assert math.isclose(
+            exact_expval.item(), float(result_tn[0]), abs_tol=ABS_TOL
+        ), f"Rank {backend.rank}: mismatch, expected {exact_expval}, got {result_tn}"
+
+    else:
+        # Rank > 0: must be hardcoded [0] (int)
+        assert (
+            isinstance(result_tn, (np.ndarray, cp.ndarray))
+            and result_tn.size == 1
+            and np.issubdtype(result_tn.dtype, np.integer)
+            and result_tn.item() == 0
+        ), f"Rank {backend.rank}: expected int array [0], got {result_tn}"
+
+    # Test with user defined hamiltonian using "hamiltonians.SymbolicHamiltonian" object.
+    comp_set_w_hamiltonian_obj = {
+        "MPI_enabled": False,
+        "MPS_enabled": False,
+        "NCCL_enabled": True,
+        "expectation_enabled": ham,
+    }
+    backend.configure_tn_simulation(comp_set_w_hamiltonian_obj)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    if backend.rank == 0:
+        # Compare numerical values
+        assert math.isclose(
+            exact_expval.item(), float(result_tn[0]), abs_tol=ABS_TOL
+        ), f"Rank {backend.rank}: mismatch, expected {exact_expval}, got {result_tn}"
+
+    else:
+        # Rank > 0: must be hardcoded [0] (int)
+        assert (
+            isinstance(result_tn, (np.ndarray, cp.ndarray))
+            and result_tn.size == 1
+            and np.issubdtype(result_tn.dtype, np.integer)
+            and result_tn.item() == 0
+        ), f"Rank {backend.rank}: expected int array [0], got {result_tn}"
+
+    # Test with user defined hamiltonian using Dictionary object form of hamiltonian.
+    ham_dict = build_observable_dict(nqubits)
+    comp_set_w_hamiltonian_dict = {
+        "MPI_enabled": False,
+        "MPS_enabled": False,
+        "NCCL_enabled": True,
+        "expectation_enabled": ham_dict,
+    }
+    backend.configure_tn_simulation(comp_set_w_hamiltonian_dict)
+    result_tn = backend.execute_circuit(circuit=qibo_circ)
+    if backend.rank == 0:
+        # Compare numerical values
+        assert math.isclose(
+            exact_expval.item(), float(result_tn[0]), abs_tol=ABS_TOL
+        ), f"Rank {backend.rank}: mismatch, expected {exact_expval}, got {result_tn}"
+
+    else:
+        # Rank > 0: must be hardcoded [0] (int)
+        assert (
+            isinstance(result_tn, (np.ndarray, cp.ndarray))
+            and result_tn.size == 1
+            and np.issubdtype(result_tn.dtype, np.integer)
+            and result_tn.item() == 0
+        ), f"Rank {backend.rank}: expected int array [0], got {result_tn}"

tests/test_expectation.py

@@ -0,0 +1,47 @@
+import math
+import random
+
+import pytest
+from qibo import Circuit, construct_backend, gates, hamiltonians
+from qibo.symbols import X, Z
+
+
+def build_observable(nqubits):
+    """Helper function to construct a target observable."""
+    hamiltonian_form = 0
+    for i in range(nqubits):
+        hamiltonian_form += 0.5 * X(i % nqubits) * Z((i + 1) % nqubits)
+
+    hamiltonian = hamiltonians.SymbolicHamiltonian(form=hamiltonian_form)
+    return hamiltonian, hamiltonian_form
+
+
+def build_circuit(nqubits, nlayers, seed=42):
+    """Helper function to construct a layered quantum circuit."""
+    random.seed(seed)
+
+    circ = Circuit(nqubits)
+    for _ in range(nlayers):
+        for q in range(nqubits):
+            circ.add(gates.RY(q=q, theta=random.uniform(-math.pi, math.pi)))
+            circ.add(gates.RZ(q=q, theta=random.uniform(-math.pi, math.pi)))
+        [circ.add(gates.CNOT(q % nqubits, (q + 1) % nqubits) for q in range(nqubits))]
+    circ.add(gates.M(*range(nqubits)))
+    return circ
+
+
+@pytest.mark.parametrize("nqubits", [2, 5, 10])
+def test_observable_expval(backend, nqubits):
+    numpy_backend = construct_backend("numpy")
+    ham, ham_form = build_observable(nqubits)
+    circ = build_circuit(nqubits=nqubits, nlayers=1)
+
+    exact_expval = numpy_backend.calculate_expectation_state(
+        hamiltonian=ham,
+        state=circ().state(),
+        normalize=False,
+    )
+
+    tn_expval = backend.expectation(circuit=circ, observable=ham_form)
+
+    assert math.isclose(exact_expval, tn_expval, abs_tol=1e-7)

tests/test_quimb_backend.py

@@ -0,0 +1,66 @@
+import copy
+import os
+
+import config
+import numpy as np
+import pytest
+import qibo
+from qibo.models import QFT
+
+
+def create_init_state(nqubits):
+    init_state = np.random.random(2**nqubits) + 1j * np.random.random(2**nqubits)
+    init_state = init_state / np.sqrt((np.abs(init_state) ** 2).sum())
+    return init_state
+
+
+def qibo_qft(nqubits, init_state, swaps):
+    circ_qibo = QFT(nqubits, swaps)
+    state_vec = circ_qibo(init_state).state(numpy=True)
+    return circ_qibo, state_vec
+
+
+@pytest.mark.parametrize(
+    "nqubits, tolerance, is_mps",
+    [(1, 1e-6, True), (2, 1e-6, False), (5, 1e-3, True), (10, 1e-3, False)],
+)
+def test_eval(nqubits: int, tolerance: float, is_mps: bool):
+    """Evaluate circuit with Quimb backend.
+
+    Args:
+        nqubits (int): Total number of qubits in the system.
+        tolerance (float): Maximum limit allowed for difference in results
+        is_mps (bool): True if state is MPS and False for tensor network structure
+    """
+    # hack quimb to use the correct number of processes
+    # TODO: remove completely, or at least delegate to the backend
+    # implementation
+    os.environ["QUIMB_NUM_PROCS"] = str(os.cpu_count())
+    import qibotn.eval_qu
+
+    init_state = create_init_state(nqubits=nqubits)
+    init_state_tn = copy.deepcopy(init_state)
+
+    # Test qibo
+    qibo.set_backend(backend=config.qibo.backend, platform=config.qibo.platform)
+
+    qibo_circ, result_sv = qibo_qft(nqubits, init_state, swaps=True)
+
+    # Convert to qasm for other backends
+    qasm_circ = qibo_circ.to_qasm()
+
+    # Test quimb
+    if is_mps:
+        gate_opt = {}
+        gate_opt["method"] = "svd"
+        gate_opt["cutoff"] = 1e-6
+        gate_opt["cutoff_mode"] = "abs"
+    else:
+        gate_opt = None
+    result_tn = qibotn.eval_qu.dense_vector_tn_qu(
+        qasm_circ, init_state_tn, gate_opt, backend=config.quimb.backend
+    ).flatten()
+
+    assert np.allclose(
+        result_sv, result_tn, atol=tolerance
+    ), "Resulting dense vectors do not match"