Initialise MPS dir

scratchpad/qtensor_MPS/constants.py

@@ -0,0 +1,29 @@
+import numpy as np
+
+_hmatrix = (
+    1 / np.sqrt(2) * np.array([[1.0, 1.0], [1.0, -1.0]], dtype=np.complex64)
+)
+_imatrix = np.array([[1.0, 0.0], [0.0, 1.0]], dtype=np.complex64)
+_xmatrix = np.array([[0.0, 1.0], [1.0, 0.0]], dtype=np.complex64)
+_ymatrix = np.array([[0.0, -1j], [1j, 0.0]], dtype=np.complex64)
+_zmatrix = np.array([[1.0, 0.0], [0.0, -1.0]], dtype=np.complex64)
+
+_cnot_matrix = np.array(
+    [
+        [1.0, 0.0, 0.0, 0.0],
+        [0.0, 1.0, 0.0, 0.0],
+        [0.0, 0.0, 0.0, 1.0],
+        [0.0, 0.0, 1.0, 0.0],
+    ]
+)
+_cnot_matrix = np.reshape(_cnot_matrix, newshape=(2, 2, 2, 2))
+
+_swap_matrix = np.array(
+    [
+        [1.0, 0.0, 0.0, 0.0],
+        [0.0, 0.0, 1.0, 0.0],
+        [0.0, 1.0, 0.0, 0.0],
+        [0.0, 0.0, 0.0, 1.0],
+    ]
+)
+_swap_matrix = np.reshape(_swap_matrix, newshape=(2, 2, 2, 2))

scratchpad/qtensor_MPS/mps.py

@@ -0,0 +1,206 @@
+import tensornetwork as tn
+import numpy as np
+
+class MPS:
+    def __init__(self, tensor_name, N, physical_dim = 1) -> None:
+        """
+               (0) physical dim
+                       |
+        (1) bond dim --@-- (2) bond dim
+        """
+        self._N = N
+        self._physical_dim = physical_dim
+        self.name = tensor_name
+
+        if N < 2:
+            raise ValueError("Number of tensors should be >= 2")
+        # Initialise as |0> = [1.0 0.0
+        #                      0.0 0.0]
+        # nodes = [tn.Node(np.array([[1.0], *[[0.0]]*(physical_dim-1)], dtype=np.complex64), name = tensor_name+ "_" +str(0))]
+        # for i in range(N-2):
+        #     node = tn.Node(np.array([[[1.0]], *[[[0.0]]]*(physical_dim-1)], dtype=np.complex64), name = tensor_name + "_" + str(i+1))
+        #     nodes.append(node)
+        # nodes.append(tn.Node(np.array([[1.0], *[[0.0]]*(physical_dim-1)], dtype=np.complex64), name = tensor_name+ "_" +str(0)))
+        nodes = [tn.Node(np.array([[[1.0]], *[[[0.0]]]*(physical_dim -1)],  dtype=np.complex64), name = tensor_name + str(i+1)) for i in range(N-2)]
+        nodes.insert(0, tn.Node(np.array([[1.0], *[[0.0]]*(physical_dim -1)],  dtype=np.complex64), name = tensor_name + str(0)))
+        nodes.append(tn.Node(np.array([[1.0], *[[0.0]]*(physical_dim -1)],  dtype=np.complex64), name = tensor_name + str(N-1)))
+
+        for i in range(1, N-2):
+            tn.connect(nodes[i].get_edge(2), nodes[i+1].get_edge(1))
+
+        if N < 3:
+            tn.connect(nodes[0].get_edge(1), nodes[1].get_edge(1))
+        else:
+            tn.connect(nodes[0].get_edge(1), nodes[1].get_edge(1))
+            tn.connect(nodes[-1].get_edge(1), nodes[-2].get_edge(2))
+
+        self._nodes = nodes
+
+    @staticmethod
+    def construct_mps_from_wavefunction(wavefunction, tensor_name, N, physical_dim = 1) -> 'MPS':
+        """
+        """
+        if wavefunction.size != physical_dim**N:
+            raise ValueError()
+
+        wavefunction = np.reshape(wavefunction, [physical_dim]*N)
+        to_split = tn.Node(wavefunction, axis_names=[str(i) for i in range(N)])
+
+        nodes = []
+        for  i in range(N-1):
+            left_edges = []
+            right_edges = []
+
+            for edge in to_split.get_all_dangling():
+                if edge.name == str(i):
+                    left_edges.append(edge)
+                else:
+                    right_edges.append(edge)
+
+            if nodes:
+                for edge in nodes[-1].get_all_nondangling():
+                    if to_split in edge.get_nodes():
+                        left_edges.append(edge)
+
+            left, right, _ = tn.split_node(to_split, left_edges, right_edges, left_name=tensor_name+str(i))
+
+            nodes.append(left)
+            to_split = right
+        to_split.name = tensor_name + str(N-1)
+        nodes.append(to_split)
+
+        mps = MPS(tensor_name, N, physical_dim)
+        mps._nodes = nodes
+        return mps
+
+    @property
+    def N(self):
+        return self._N
+
+    @property
+    def physical_dim(self):
+        return self._physical_dim
+
+    def get_mps_nodes(self, original) -> list[tn.Node]:
+        if original:
+            return self._nodes
+
+        nodes, edges = tn.copy(self._nodes)
+        return list(nodes.values())
+
+    def get_mps_node(self, index, original) -> tn.Node:
+        return self.get_mps_nodes(original)[index]
+
+    def get_tensors(self, original) -> list[tn.Node]:
+        nodes = []
+
+        if original:
+            nodes = self._nodes
+        else:
+            nodes, edges = tn.copy(self._nodes)
+
+        return list([node.tensor for node in nodes])
+
+    def get_tensor(self, index, original) -> tn.Node:
+        return self.get_tensors(original)[index]
+
+    def is_unitary():
+        pass
+
+    def get_wavefunction(self) -> np.array:
+        nodes = self.get_mps_nodes(False)
+        curr = nodes.pop(0)
+
+        for node in nodes:
+            curr = tn.contract_between(curr, node)
+
+        wavefunction = np.reshape(curr.tensor, newshape=(self._physical_dim ** self._N))
+        return wavefunction
+
+    def apply_single_qubit_gate(self, gate, index) -> None:
+    #     """
+    #     Assumption: Gates are unitary
+
+    #      0
+    #      |
+    #     gate
+    #      |
+    #      1
+
+    #      |
+    #     MPS
+    #      |
+    #     """
+        mps_index_edge = list(self._nodes[index].get_all_dangling())[0]
+        gate_edge = gate[1]
+        temp_node = tn.connect(mps_index_edge, gate_edge)
+
+        new_node = tn.contract(temp_node, name=self._nodes[index].name)
+        self._nodes[index] = new_node
+
+    def apply_two_qubit_gate(self, gate, operating_qubits):
+        """
+            0  1
+            |  |
+            gate
+            |  |
+            2  3
+
+            a   b
+            |   |
+             MPS
+        """
+        # reshape the gate to order-4 tensor
+        # Assumimg operating_qubits are sorted
+        mps_indexA = self.get_mps_node(operating_qubits[0], True).get_all_dangling().pop()
+        mps_indexB = self.get_mps_node(operating_qubits[1], True).get_all_dangling().pop()
+
+        temp_nodesA = tn.connect(mps_indexA, gate.get_edge(2))
+        temp_nodesB = tn.connect(mps_indexB, gate.get_edge(3))
+        left_gate_edge = gate.get_edge(0)
+        right_gate_edge = gate.get_edge(1)
+
+        new_node = tn.contract_between(self._nodes[operating_qubits[0]], self._nodes[operating_qubits[1]])
+        node_gate_edge = tn.flatten_edges_between(new_node, gate)
+        new_node = tn.contract(node_gate_edge)
+
+
+        left_connected_edge = None
+        right_connected_edge = None
+
+        for edge in new_node.get_all_nondangling():
+            index = int(edge.node1.name.split(self.name)[-1])
+
+            if index <= operating_qubits[0]:
+                left_connected_edge = edge
+            else:
+                right_connected_edge = edge
+
+        left_edges = []
+        right_edges = []
+
+        for edge in (left_gate_edge, left_connected_edge):
+            if edge != None:
+                left_edges.append(edge)
+
+        for edge in (right_gate_edge, right_connected_edge):
+            if edge != None:
+                right_edges.append(edge)
+
+        u, s, vdag, _ = tn.split_node_full_svd(
+            new_node,
+            left_edges=left_edges,
+            right_edges=right_edges
+        )
+
+        new_left = u
+        new_right = tn.contract_between(s, vdag)
+
+        new_left.name = self._nodes[operating_qubits[0]].name
+        new_right.name = self._nodes[operating_qubits[1]].name
+
+        self._nodes[operating_qubits[0]] = new_left
+        self._nodes[operating_qubits[1]] = new_right
+
+
+

scratchpad/qtensor_MPS/mps_operation.py

@@ -0,0 +1,15 @@
+import numpy as np
+import tensornetwork as tn
+import numpy as np
+from constants import _xmatrix, _cnot_matrix, _hmatrix
+from copy import deepcopy
+
+def xgate() -> tn.Node:
+    return tn.Node(deepcopy(_xmatrix), name="xgate")
+
+def cnot() -> tn.Node:
+    return tn.Node(deepcopy(_cnot_matrix), name="cnot")
+
+def hgate() -> tn.Node:
+    return tn.Node(deepcopy(_hmatrix), name="hgate")
+

scratchpad/qtensor_MPS/notes.md

@@ -0,0 +1,20 @@
+# Comments:
+- ~Correct the error in one-qubit code : Create function to draw qc/ reshape (2x2)
+- ~Add test for random values
+- ~Add 2 qubit code
+- Add expectation values
+- Implememt calculation of probablities
+   _____  
+  |__A__|
+   | | |
+  _______ 
+  |__A'__|   
+
+- Implement sampling
+- 
+- https://www.nature.com/articles/s41586-023-06096-3
+# Circuits to test
+
+- MPS GHZ implementation
+- How IBM has simulated 127 qubits? (https://arxiv.org/abs/2309.15642)
+- 

scratchpad/qtensor_MPS/test.py

@@ -0,0 +1,85 @@
+import numpy as np
+import tensornetwork as tn
+import xyzpy as xyz
+from mps_operation import xgate, cnot, hgate
+from mps import MPS
+
+def test_from_wavefunction_all_zero_state():
+    wavefunction = np.array([1, 0, 0, 0, 0, 0, 0, 0])
+    mps = MPS.construct_mps_from_wavefunction(wavefunction, 'q', 3, 2)
+
+    # kwargs.setdefault("legend", True)
+    # kwargs.setdefault("compass", True)
+    # kwargs.setdefault("compass_labels", mps.inds)
+    # xyz.visualize_tensor(mps.get_tensors(False), legend=True, compass=True)
+
+    assert isinstance(mps, MPS)
+    assert mps.N == 3
+    assert mps.physical_dim == 2
+    assert np.allclose(mps.get_wavefunction(), wavefunction)
+
+def test_from_wavefunction_random():
+    n = 3
+    wavefunction = np.random.rand(2**n)
+    wavefunction = wavefunction / np.linalg.norm(wavefunction, ord = 2)
+    mps = MPS.construct_mps_from_wavefunction(wavefunction, 'q', n, 2)
+    assert np.allclose(mps.get_wavefunction(), wavefunction)
+#def add random function
+#  
+
+# Entangled
+
+def test_apply_one_qubit_mps_operation_xgate():
+    mps = MPS("q", 2, 2)
+    print("Before applying gate: ", mps.get_wavefunction())
+
+    # nodes = mps.get_tensors(False)
+    # qubits = [node[0] for node in nodes]
+    # q0q1 = 00
+    # On apply x gate |00> -> |10>
+    mps.apply_single_qubit_gate(xgate(), 0)
+
+    print("After applying x gate: ", mps.get_wavefunction())
+
+def test_apply_twoq_cnot_two_qubits():
+    """Tests for correctness of final wavefunction after applying a CNOT
+    to a two-qubit MPS.
+    """
+    # In the following tests, the first qubit is always the control qubit.
+    # Check that CNOT|10> = |11>
+    mps = MPS("q", 2, 2)
+    mps.apply_single_qubit_gate(xgate(), 0)
+
+    print("After applying x gate: ", mps.get_wavefunction())
+
+    mps.apply_two_qubit_gate(cnot(), [0, 1])
+    print("After applying cnot gate: ", mps.get_wavefunction())
+
+
+    correct = np.array([0.0, 0.0, 0.0, 1.0], dtype=np.complex64)
+    curr = mps.get_wavefunction()
+    # assert np.allclose(curr, correct)
+
+def test_apply_gate_for_bell_circuit():
+    mps = MPS("q", 2, 2)
+    mps.apply_single_qubit_gate(hgate(), 0)
+    mps.apply_two_qubit_gate(cnot(), [0,1])
+
+    print("Wavefunction for bell circuit: ", mps.get_wavefunction())
+
+def test_apply_gate_for_ghz_circuit():
+    mps = MPS("q", 3, 2)
+    mps.apply_single_qubit_gate(hgate(), 0)
+    mps.apply_two_qubit_gate(cnot(), [0,1])
+    mps.apply_two_qubit_gate(cnot(), [1,2])
+
+    print("Wavefunction for ghz circuit: ", mps.get_wavefunction())
+
+# test_from_wavefunction_all_zero_state()
+# test_apply_one_qubit_mps_operation_xgate()
+# test_from_wavefunction_random()
+
+#test_apply_twoq_cnot_two_qubits()
+
+# test_apply_gate_for_bell_circuit()
+test_apply_gate_for_ghz_circuit()

scratchpad/qtensor_MPS/test2.py

@@ -0,0 +1,29 @@
+def helper():
+    pass
+
+# Driver Code
+if __name__ == '__main__':
+    t = int(input())
+    for _ in range(t):
+        n, k = [int(i) for i in input().split()]
+        a = [int(i) for i in input().split()]
+        b = [int(i) for i in input().split()]
+
+        for elem in a:
+            if elem not in map:
+                map[elem] = 1
+            map[elem] += 1
+
+        n = map.keys() #distinct elem in a
+
+        cnt = 0 # elements in b not in b
+        temp = []
+
+        for elem in b:
+            if elem not in a:
+                cnt += 1
+
+        if cnt > k:
+            print()
+
+    print(a)