TEBD merging 3 options

src/qibotn/backends/quimb.py

@@ -24,12 +24,40 @@ def __init__(self, runcard):
             else:
                 raise TypeError("MPS_enabled has an unexpected type")
 
+            global TEBD_enabled
+            global TEBD_option
+            TEBD_enabled = runcard.get("TEBD_enabled")
+            tebd_enabled_value = runcard.get("TEBD_enabled")
+            TEBD_option = runcard.get("TEBD_option")
+            if TEBD_option == True:
+                if tebd_enabled_value is True:
+                    self.tebd_opts = {"dt": 1e-4, "initial_state": "00", "tot_time": 1}
+                elif tebd_enabled_value is False:
+                    self.tebd_opts = None
+                elif isinstance(tebd_enabled_value, dict):
+                    self.tebd_opts = tebd_enabled_value
+            else:
+                if tebd_enabled_value is True:
+                    self.tebd_opts = {
+                        "dt": 1e-4,
+                        "hamiltoninan": "XXZ",
+                        "initial_state": "00",
+                        "tot_time": 1,
+                    }
+                elif tebd_enabled_value is False:
+                    self.tebd_opts = None
+                elif isinstance(tebd_enabled_value, dict):
+                    self.tebd_opts = tebd_enabled_value
+
         else:
             self.MPI_enabled = False
+            self.TEBD_enabled = False
             self.MPS_enabled = False
             self.NCCL_enabled = False
             self.expectation_enabled = False
+            self.TEBD_option = False
             self.mps_opts = None
+            self.tebd_opts = None
 
         self.name = "qibotn"
         self.quimb = quimb
@@ -76,9 +104,17 @@ def execute_circuit(
                 NotImplementedError, "QiboTN quimb backend cannot support expectation"
             )
 
-        state = eval.dense_vector_tn_qu(
-            circuit.to_qasm(), initial_state, self.mps_opts, backend="numpy"
-        )
+        if TEBD_enabled:
+            nqubits = circuit.nqubits
+            if TEBD_option == True:
+                state = eval.tebd_tn_qu(circuit, self.tebd_opts)
+            else:
+                state = eval.tebd_tn_qu_2(circuit, self.tebd_opts)
+
+        else:
+            state = eval.dense_vector_tn_qu(
+                circuit.to_qasm(), initial_state, self.mps_opts, backend="numpy"
+            )
 
         if return_array:
             return state.flatten()

src/qibotn/eval_qu.py

@@ -14,7 +14,6 @@ def init_state_tn(nqubits, init_state_sv):
     """
 
     dims = tuple(2 * np.ones(nqubits, dtype=int))
-
     return qtn.tensor_1d.MatrixProductState.from_dense(init_state_sv, dims)
 
 
@@ -44,3 +43,76 @@ def dense_vector_tn_qu(qasm: str, initial_state, mps_opts, backend="numpy"):
     amplitudes = interim.to_dense(backend=backend)
 
     return amplitudes
+
+
+def tebd_tn_qu(circuit, tebd_opts, mps_opts):
+    """Circuit based TEBD which returns the final evolved state as a dense
+    vector."""
+
+    init_state = tebd_opts["initial_state"]
+    dt = tebd_opts["dt"]
+    T = tebd_opts["tot_time"]
+    nqubits = circuit.nqubits
+
+    initial_state = qtn.MPS_computational_state(init_state)
+
+    import qiskit.qasm2 as qasm
+    from qiskit import QuantumCircuit, transpile
+
+    circ = QuantumCircuit(nqubits)
+    unitary = circuit.unitary()
+
+    qubit_arr = []
+    for i in range(0, nqubits):
+        qubit_arr.append(i)
+    circ.unitary(unitary, qubit_arr)
+
+    transpiled_circ = transpile(circ, basis_gates=["rx", "ry", "rz", "cx"])
+
+    qasm_str = qasm.dumps(transpiled_circ)
+
+    circ_cls = qtn.circuit.CircuitMPS
+    circ_quimb = circ_cls.from_openqasm2_str(
+        qasm_str, psi0=initial_state, gate_opts=mps_opts
+    )
+
+    psi0 = initial_state
+    for i in np.arange(0, T, dt):
+
+        circ_quimb.psi0 = psi0
+        interim = circ_quimb.psi.full_simplify(seq="DRC")
+        psi0 = interim  # - should this be done to update the initial state to the next evolving state
+
+    amplitudes = psi0.to_dense()
+    return amplitudes
+
+
+def tebd_tn_qu_2(circuit, tebd_opts):
+    """Circuit based TEBD which returns the final evolved state as a dense
+    vector."""
+    dt = tebd_opts["dt"]
+    tot_time = tebd_opts["tot_time"]
+    init_state = tebd_opts["initial_state"]
+    nqubits = circuit.nqubits
+
+    initial_state = qtn.MPS_computational_state(init_state)
+
+    i = -1
+    uni = circuit.unitary()
+    h = np.divide((np.log(uni)), -1 * i * dt)
+
+    from qibo import hamiltonians
+
+    ham = hamiltonians.Hamiltonian(nqubits, h)
+    ham_quimb = ham.matrix
+    H = qtn.LocalHam1D(2, H2=ham_quimb)
+
+    tebd = qtn.TEBD(initial_state, H)
+
+    ts = np.arange(0, 1, dt)
+    states = {}
+    for t in tebd.at_times(ts, tol=1e-3):
+        states.update({None: t.to_dense()})
+
+    state = np.array(list(states.values()))[-1]
+    return state

src/qibotn/experiment_with_logging.py

@@ -0,0 +1,53 @@
+# importing necessary packages
+import numpy as np
+import quimb.tensor as qtn
+from qibo import hamiltonians
+
+# tebd opts
+dt = 1e-4
+nqubits = 5
+init_state = "10101"
+tot = 1
+
+# casting ham as crt using td
+ham = hamiltonians.XXZ(nqubits=nqubits, dense=False)
+circuit = ham.circuit(dt=dt)
+
+# openqasm workaround experiments
+# print(circuit.decompose()) # decompose doc says it Returns: Circuit that contains only gates that are supported by OpenQASM and has the same effect as the original circuit.
+
+# build initial state
+psi0 = qtn.MPS_computational_state(init_state)
+
+# extract symb rep terms of ham
+terms_dict = {}
+i = 0
+list_of_terms = ham.terms
+for t in list_of_terms:
+    terms_dict.update({None: t.matrix})
+    i = i + 1
+
+# build quimb ham and tebd object with ts time range
+H = qtn.LocalHam1D(nqubits, H2=terms_dict)
+tebd = qtn.TEBD(psi0, H)
+ts = np.arange(0, tot, dt)
+
+result = {}
+# write evol dense vect at times ts into txt
+# file_path = "tebd_log2.txt"
+# with open(file_path, 'w') as file:
+for t in tebd.at_times(ts, tol=tot):
+    result.update({None: t.to_dense()})
+# file.write("\n"+str(t.to_dense()))
+
+res = list(result.values())
+res = np.array(res[-1])
+print(res)
+"""# extract the final evol dense vect
+with open(file_path, 'r') as file:
+        content = file.read()
+
+state_str = (content.rsplit(']]',2)[-2])+"]]"
+state = np.array(eval(state_str))
+
+print(state)"""

src/qibotn/tebd.py

@@ -0,0 +1,58 @@
+import numpy as np
+import quimb.tensor as qtn
+from qibo.config import raise_error
+
+
+def init_state_tn_tebd(initial_state):
+    """Creates a inital MPS from a binary string."""
+
+    initial_state = qtn.MPS_computational_state(initial_state)
+    return initial_state
+
+
+def tebd_quimb(circuit, tebd_opts):
+    """Symbolic Hamiltonian based TEBD which returns the final evolved state as
+    a dense vector."""
+
+    hamiltonian = tebd_opts["hamiltonian"]
+    dt = tebd_opts["dt"]
+    initial_state = tebd_opts["initial_state"]
+    tot_time = tebd_opts["tot_time"]
+    nqubits = circuit.nqubits
+
+    init_state = init_state_tn_tebd(initial_state)
+    from qibo import hamiltonians
+
+    if hamiltonian == "TFIM":
+        ham = hamiltonians.TFIM(nqubits=nqubits, dense=False)
+    elif hamiltonian == "NIX":
+        ham = hamiltonians.X(nqubits=nqubits, dense=False)
+    elif hamiltonian == "NIY":
+        ham = hamiltonians.Y(nqubits=nqubits, dense=False)
+    elif hamiltonian == "NIZ":
+        ham = hamiltonians.Z(nqubits=nqubits, dense=False)
+    elif hamiltonian == "XXZ":
+        ham = hamiltonians.XXZ(nqubits=nqubits, dense=False)
+    elif hamiltonian == "MC":
+        ham = hamiltonians.MaxCut(nqubits=nqubits, dense=False)
+    else:
+        raise_error(NotImplementedError, "QiboTN does not support custom hamiltonians")
+
+    terms_dict = {}
+    i = 0
+    list_of_terms = ham.terms
+    for t in list_of_terms:
+        terms_dict.update({None: t.matrix})
+        i = i + 1
+
+    H = qtn.LocalHam1D(nqubits, H2=terms_dict)
+
+    tebd = qtn.TEBD(init_state, H)
+    ts = np.arange(0, tot_time, dt)
+
+    states = {}
+    for t in tebd.at_times(ts, tol=1e-3):
+        states.update({None: t.to_dense()})
+
+    state = np.array(list(states.values()))[-1]
+    return state

tests/test_quimb_backend_tebd.py

@@ -0,0 +1,68 @@
+import copy
+import os
+
+import config
+import numpy as np
+import pytest
+import qibo
+from qibo import Circuit, gates
+
+
+def create_init_state(nqubits):
+    init_state = np.ones(nqubits)
+    return init_state
+
+
+def qibo_crt(nqubits, init_state, dt):
+    from numpy import pi
+
+    circ_qibo = Circuit(3)
+    circ_qibo.add(gates.RY(0, theta=pi / 2))
+    circ_qibo.add(gates.RY(1, theta=pi / 4))
+    circ_qibo.add(gates.CNOT(0, 1))
+    state_vec = circ_qibo(init_state).state(numpy=True)
+    return circ_qibo, state_vec
+
+
+@pytest.mark.parametrize(
+    "nqubits, tolerance, is_tebd",
+    [(4, 1e-6, True), (5, 1e-6, False), (6, 1e-3, True), (10, 1e-3, False)],
+)
+def test_eval(nqubits: int, tolerance: float, is_tebd: 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_tebd (bool): True if user selects is TEBD and False if otherwise
+    """
+    # 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())
+
+    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)
+    import eval_qu as ev
+
+    qibo_circ, result_sv = qibo_crt(nqubits, init_state, dt=1e-4)
+
+    # Test quimb
+    if is_tebd:
+        gate_opt = {}
+        gate_opt["dt"] = 1e-4
+        gate_opt["initial_state"] = "101"
+        gate_opt["tot_time"] = 1
+    else:
+        gate_opt = None
+    result_tn = ev.tebd_tn_qu(qibo_circ, gate_opt).flatten()
+
+    assert np.allclose(
+        result_sv, result_tn, atol=tolerance
+    ), "Resulting dense vectors do not match"
+
+
+print(test_eval(nqubits=3, tolerance=1e-6, is_tebd=True))