From 89288f848835f3a0992cfcd23720a1889f9b9712 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Mon, 26 Jun 2023 19:13:48 +0200 Subject: [PATCH 01/48] update --- graphix/pattern.py | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/graphix/pattern.py b/graphix/pattern.py index 6c14adcb4..032ec6475 100644 --- a/graphix/pattern.py +++ b/graphix/pattern.py @@ -814,6 +814,11 @@ def connected_edges(self, node, edges): connected: set of tuple set of connected edges """ + # TODO modify that by using the graph nx.graph.edges(node)? and cached get_graph()? + # like in def get_measurement_order_from_flow(self): with self.get_graph() + # FIXME + # BUG + connected = set() for edge in edges: if edge[0] == node: @@ -1011,7 +1016,10 @@ def get_max_degree(self): degree = g.degree() max_degree = max([i for i in dict(degree).values()]) return max_degree - + + # TODO functools.cache() It is called in get measurement order from (g)flow + # + # It is called in get measurement order from (g)flow def get_graph(self): """returns the list of nodes and edges from the command sequence, extracted from 'N' and 'E' commands. From a6da6d3adb9ac5286dfb87c4c7ae4f0bfc063195 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Wed, 28 Jun 2023 15:09:44 +0200 Subject: [PATCH 02/48] up --- graphix/sim/statevec.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 8ebbc5f08..23cff4968 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -258,6 +258,8 @@ def ptrace(self, qargs): rho = np.tensordot(psi, psi.conj(), axes=(qargs, qargs)) # density matrix rho = np.reshape(rho, (2**nqubit_after, 2**nqubit_after)) evals, evecs = np.linalg.eig(rho) # back to statevector + # NOTE works since only one 1 in the eigenvalues corresponding to the state + # TODO use np.eigh since rho is Hermitian? self.psi = np.reshape(evecs[:, np.argmax(evals)], (2,) * nqubit_after) def remove_qubit(self, qarg): From 187e650d0efaeddb8e2e467814b32173aed5b9b9 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Fri, 29 Mar 2024 14:50:10 +0100 Subject: [PATCH 03/48] statevector constructor for input initialization --- .pre-commit-config.yaml | 7 ++ .vscode/settings.json | 11 +++ graphix/ops.py | 20 ++--- graphix/sim/graphix.code-workspace | 7 ++ graphix/sim/statevec.py | 123 ++++++++++++++++++++++++++--- graphix/states.py | 74 +++++++++++++++++ 6 files changed, 223 insertions(+), 19 deletions(-) create mode 100644 .pre-commit-config.yaml create mode 100644 .vscode/settings.json create mode 100644 graphix/sim/graphix.code-workspace create mode 100644 graphix/states.py diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml new file mode 100644 index 000000000..ff5c14609 --- /dev/null +++ b/.pre-commit-config.yaml @@ -0,0 +1,7 @@ +repos: +- repo: https://github.com/psf/black + rev: 22.8.0 + hooks: + - id: black + args: [--line-length=120] + files: ^(graphix|test)/ diff --git a/.vscode/settings.json b/.vscode/settings.json new file mode 100644 index 000000000..0a0133979 --- /dev/null +++ b/.vscode/settings.json @@ -0,0 +1,11 @@ +{ + "workbench.colorCustomizations": { + "activityBar.activeBackground": "#287d9e", + "activityBar.background": "#287d9e", + "activityBar.foreground": "#e7e7e7", + "activityBar.inactiveForeground": "#e7e7e799", + "activityBarBadge.background": "#e599d0", + "activityBarBadge.foreground": "#15202b" + }, + "peacock.color": "#1e5d75" +} \ No newline at end of file diff --git a/graphix/ops.py b/graphix/ops.py index 33d6d507d..6c90d101c 100644 --- a/graphix/ops.py +++ b/graphix/ops.py @@ -7,15 +7,17 @@ import numpy as np - -class States: - plus = np.array([1.0 / np.sqrt(2), 1.0 / np.sqrt(2)]) # plus - minus = np.array([1.0 / np.sqrt(2), -1.0 / np.sqrt(2)]) # minus - zero = np.array([1.0, 0.0]) # zero - one = np.array([0.0, 1.0]) # one - iplus = np.array([1.0 / np.sqrt(2), 1.0j / np.sqrt(2)]) # +1 eigenstate of Pauli Y - iminus = np.array([1.0 / np.sqrt(2), -1.0j / np.sqrt(2)]) # -1 eigenstate of Pauli Y - vec = [plus, minus, zero, one, iplus, iminus] +# TODO modify that + +# Everywhere this is called. use StateVec(State)) +# class States: +# plus = np.array([1.0 / np.sqrt(2), 1.0 / np.sqrt(2)]) # plus +# minus = np.array([1.0 / np.sqrt(2), -1.0 / np.sqrt(2)]) # minus +# zero = np.array([1.0, 0.0]) # zero +# one = np.array([0.0, 1.0]) # one +# iplus = np.array([1.0 / np.sqrt(2), 1.0j / np.sqrt(2)]) # +1 eigenstate of Pauli Y +# iminus = np.array([1.0 / np.sqrt(2), -1.0j / np.sqrt(2)]) # -1 eigenstate of Pauli Y +# vec = [plus, minus, zero, one, iplus, iminus] class Ops: diff --git a/graphix/sim/graphix.code-workspace b/graphix/sim/graphix.code-workspace new file mode 100644 index 000000000..9e68e72b8 --- /dev/null +++ b/graphix/sim/graphix.code-workspace @@ -0,0 +1,7 @@ +{ + "folders": [ + { + "path": ".." + } + ] +} \ No newline at end of file diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 23cff4968..00a763a34 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -1,10 +1,23 @@ from copy import deepcopy +import typing +from typing_extensions import Annotated +from annotated_types import Ge import numpy as np import graphix.sim.base_backend from graphix.clifford import CLIFFORD, CLIFFORD_CONJ, CLIFFORD_MUL from graphix.ops import Ops +import graphix.states +import graphix.pauli +import functools +import warnings + +# Python >= 3.9 +# from collections.abc import Iterable # or use Protocols? +# https://stackoverflow.com/questions/49427944/typehints-for-sized-iterable-in-python +# Python >= 3.8 +# typing.Iterable[T] class StatevectorBackend(graphix.sim.base_backend.Backend): @@ -177,26 +190,116 @@ def meas_op(angle, vop=0, plane="XY", choice=0): [[[[1, 0], [0, 0]], [[0, 0], [1, 0]]], [[[0, 1], [0, 0]], [[0, 0], [0, 1]]]], dtype=np.complex128, ) +PositiveInt = Annotated[int, Ge(0)] # includes 0 class Statevec: - """Simple statevector simulator""" - - def __init__(self, nqubit=1, plus_states=True): + """Statevector object""" + + # TODO at this stage no need for indices just be careful of the ordering in add_nodes + def __init__( + self, + nqubit: typing.Optional[PositiveInt] = None, + state: typing.Union[ + graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[complex] + ] = graphix.states.BasicStates.PLUS, + ): + # always infer nqubit from data + # nqubit = none et pas 1 + # if nqubit is None: on exige iterable fini et on prend son nombre d'élément + # also allow external data. """Initialize statevector Parameters ---------- - nqubit : int, optional: + state : is either + - a single state (:class:`graphix.states.State` object). THen prepares all nodes in that state (tensor product) + - a dictionary mapping the inputs to a :class:`graphix.states.State` object + - an arbitrary :class:`graphix.statevec.Statevec` object (arbitrary input) # TODO work on that since just copy? + nqubit : int, optional: ignored if iterable passed (State, direct data) number of qubits. Defaults to 1. - plus_states : bool, optional + # plus_states : bool, optional whether or not to start all qubits in + state or 0 state. Defaults to + + + Defaults to |+> states and 1 qubit. + If nqubit > 1 and only one state : tensor all of them. Use the tensor method instead of hard code. """ - if plus_states: - self.psi = np.ones((2,) * nqubit) / 2 ** (nqubit / 2) - else: - self.psi = np.zeros((2,) * nqubit) - self.psi[(0,) * nqubit] = 1 + # convert single qubit State object to Statevector + # NOTE make plane, angle attributes of Statevec? + # this will be called by the nqb > 1 case + # instantiate a one + if nqubit == 0: + warnings.warn(f"Called Statevec with 0 qubits. Ignoring the state.") + self.psi = np.array(1, dtype=np.complex128) + + # works only for planar states. Deal with all kind of states? + elif isinstance(state, graphix.states.State): + + if nqubit is None: + raise ValueError("Incorrect value for nqubit.") + + vec = state.get_statevector() + + if nqubit == 1: # or None + self.psi = vec + + # build tensor product |state>^{\otimes nqubit} + # can only be >1 int. + else: + # build tensor product + # comma in tuple is for disambiguation with paranthesed expression + tmp_psi = functools.reduce(np.kron, (vec,) * nqubit) + # reshape + self.psi = tmp_psi.reshape((2,) * nqubit) + + # nqubit is None : on prend la longeur de l'iterable + elif isinstance(state, typing.Iterable): + # iterateur + it = iter(state) + head = next(it) + # type constraint in head doesn't progpagate to all elts + if isinstance(head, graphix.states.State): + # assert isinstance(head, typing.Iterator[graphix.states.State]) + if nqubit is None: + # liste persistante state pour eviter la transience + states = [head] + list(it) + nqubit = len(states) + # sinon on prend nqubit elts + else: # ignore for now + states = [head] + [next(it) for _ in range(nqubit - 1)] + + list_of_sv = [s.get_statevector() for s in states] + tmp_psi = functools.reduce(np.kron, list_of_sv) + # reshape + self.psi = tmp_psi.reshape((2,) * nqubit) + + else: + if nqubit is None: + states = [head] + list(it) + + inferred_size = len(states) + + if inferred_size & (inferred_size - 1) != 0: + raise ValueError(f"Statevector size must be a power of two but is {inferred_size}.") + + nqubit = inferred_size.bit_length() - 1 + + else: # ignore for now + states = [head] + [next(it) for _ in range(2**nqubit - 1)] + + psi = np.array(states) + + if not np.allclose(np.sqrt(np.sum(np.abs(psi) ** 2)), 1): + raise ValueError(f"Statevector must be normalized to one.") + + # just reshape + # NOTE too many conversions to numpy arrays? + self.psi = psi.reshape((2,) * nqubit) + + # if already a valid statevec just copy it. + if isinstance(state, Statevec): + assert nqubit is None or len(state.flatten()) == 2**nqubit + self.psi = state.psi.copy() def __repr__(self): return f"Statevec, data={self.psi}, shape={self.dims()}" diff --git a/graphix/states.py b/graphix/states.py new file mode 100644 index 000000000..edc49bcef --- /dev/null +++ b/graphix/states.py @@ -0,0 +1,74 @@ +""" +quantum states and operators +""" + +import numpy as np +import pydantic +import graphix.pauli +import abc + +# generic class State for all States +class State(abc.ABC): + @abc.abstractmethod + def get_statevector(self) -> np.ndarray: + pass + + +# don't turn it into Statevec here +# Weird not to allow all states? +# Made it inherit from more generic State class. +class PlanarState(pydantic.BaseModel, State): + """Light object used to instantiate backends. + doesn't cover all possible states but this is + covered in :class:`graphix.sim.statevec.Statevec` + and :class:`graphix.sim.densitymatrix.DensityMatrix` + constructors. + + :param plane: One of the three planes (XY, XZ, YZ) + :type plane: :class:`graphix.pauli.Plane` + :param angle: angle IN RADIANS + :type angle: complex + :return: State + :rtype: :class:`graphix.states.State` object + """ + + plane: graphix.pauli.Plane + angle: float + + def __repr__(self): + return f"State objects defined in plane {self.plane} with angle {self.angle}." + + def get_statevector(self) -> np.ndarray: + if self.plane == graphix.pauli.Plane.XY: + return np.array([1, np.exp(1j * self.angle)]) / np.sqrt(2) + + if self.plane == graphix.pauli.Plane.YZ: + return np.array([np.cos(self.angle / 2), 1j * np.sin(self.angle / 2)]) + + if self.plane == graphix.pauli.Plane.XZ: + return np.array([np.cos(self.angle / 2), np.sin(self.angle / 2)]) + # other case never happens since exhaustive + assert False + + +# States namespace for input initialization. +class BasicStates: + ZERO = PlanarState(plane=graphix.pauli.Plane.XZ, angle=0) + ONE = PlanarState(plane=graphix.pauli.Plane.XZ, angle=np.pi) + PLUS = PlanarState(plane=graphix.pauli.Plane.XY, angle=0) + MINUS = PlanarState(plane=graphix.pauli.Plane.XY, angle=np.pi) + PLUS_I = PlanarState(plane=graphix.pauli.Plane.XY, angle=np.pi / 2) + MINUS_I = PlanarState(plane=graphix.pauli.Plane.XY, angle=-np.pi / 2) + + +# Plane.cos.value Plane.cos is an Axis, Axis.value = 0,1,2 (enum) + +# Everywhere this is called. use StateVec(State)) +# class States: +# plus = np.array([1.0 / np.sqrt(2), 1.0 / np.sqrt(2)]) # plus +# minus = np.array([1.0 / np.sqrt(2), -1.0 / np.sqrt(2)]) # minus +# zero = np.array([1.0, 0.0]) # zero +# one = np.array([0.0, 1.0]) # one +# iplus = np.array([1.0 / np.sqrt(2), 1.0j / np.sqrt(2)]) # +1 eigenstate of Pauli Y +# iminus = np.array([1.0 / np.sqrt(2), -1.0j / np.sqrt(2)]) # -1 eigenstate of Pauli Y +# vec = [plus, minus, zero, one, iplus, iminus] From 8d7f18b335cfe9eff40e79e41fa0a9667b135d84 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Fri, 29 Mar 2024 18:48:26 +0100 Subject: [PATCH 04/48] fixed TN files for BasicStates modification --- graphix/sim/tensornet.py | 48 ++++++++++++++++++++-------------- tests/test_statevec_backend.py | 11 +++++--- tests/test_tnsim.py | 5 ++-- 3 files changed, 38 insertions(+), 26 deletions(-) diff --git a/graphix/sim/tensornet.py b/graphix/sim/tensornet.py index 30d942dca..ca3cd6aa1 100644 --- a/graphix/sim/tensornet.py +++ b/graphix/sim/tensornet.py @@ -7,6 +7,10 @@ from graphix.clifford import CLIFFORD, CLIFFORD_CONJ, CLIFFORD_MUL from graphix.ops import Ops, States +from graphix.states import BasicStates +import string +from copy import deepcopy + class TensorNetworkBackend: @@ -259,17 +263,17 @@ def add_qubit(self, index, state="plus"): ind = gen_str() tag = str(index) if state == "plus": - vec = States.plus + vec = BasicStates.PLUS.get_statevector() elif state == "minus": - vec = States.minus + vec = BasicStates.MINUS.get_statevector() elif state == "zero": - vec = States.zero + vec = BasicStates.ZERO.get_statevector() elif state == "one": - vec = States.one + vec = BasicStates.ONE.get_statevector() elif state == "iplus": - vec = States.iplus + vec = BasicStates.PLUS_I.get_statevector() elif state == "iminus": - vec = States.iminus + vec = BasicStates.MINUS_I.get_statevector() else: assert state.shape == (2,), "state must be 2-element np.ndarray" assert np.isclose(np.linalg.norm(state), 1), "state must be normalized" @@ -357,17 +361,17 @@ def measure_single(self, index, basis="Z", bypass_probability_calculation=True, raise Warning("Measurement outcome is chosen but the basis state was given.") proj_vec = basis elif basis == "Z" and result == 0: - proj_vec = States.zero + proj_vec = BasicStates.ZERO.get_statevector() elif basis == "Z" and result == 1: - proj_vec = States.one + proj_vec = BasicStates.ONE.get_statevector() elif basis == "X" and result == 0: - proj_vec = States.plus + proj_vec = BasicStates.PLUS.get_statevector() elif basis == "X" and result == 1: - proj_vec = States.minus + proj_vec = BasicStates.MINUS.get_statevector() elif basis == "Y" and result == 0: - proj_vec = States.iplus + proj_vec = BasicStates.PLUS_I.get_statevector() elif basis == "Y" and result == 1: - proj_vec = States.iminus + proj_vec = BasicStates.MINUS_I.get_statevector() else: raise ValueError("Invalid measurement basis.") else: @@ -414,13 +418,17 @@ def set_graph_state(self, nodes, edges): if node not in ind_dict.keys(): ind = gen_str() self._dangling[str(node)] = ind - self.add_tensor(Tensor(States.plus, [ind], [str(node), "Open"])) + self.add_tensor(Tensor(BasicStates.PLUS.get_statevector(), [ind], [str(node), "Open"])) continue dim_tensor = len(vec_dict[node]) tensor = np.array( [ - outer_product([States.vec[0 + 2 * vec_dict[node][i]] for i in range(dim_tensor)]), - outer_product([States.vec[1 + 2 * vec_dict[node][i]] for i in range(dim_tensor)]), + outer_product( + [BasicStates.VEC[0 + 2 * vec_dict[node][i]].get_statevector() for i in range(dim_tensor)] + ), + outer_product( + [BasicStates.VEC[1 + 2 * vec_dict[node][i]].get_statevector() for i in range(dim_tensor)] + ), ] ) * 2 ** (dim_tensor / 4 - 1.0 / 2) self.add_tensor(Tensor(tensor, ind_dict[node], [str(node), "Open"])) @@ -452,10 +460,10 @@ def get_basis_coefficient(self, basis, normalize=True, indices=None, **kwagrs): node = str(indices[i]) exp = len(indices) - i - 1 if (basis // 2**exp) == 1: - state_out = States.one # project onto |1> + state_out = BasicStates.ONE.get_statevector() # project onto |1> basis -= 2**exp else: - state_out = States.zero # project onto |0> + state_out = BasicStates.ZERO.get_statevector() # project onto |0> tensor = Tensor(state_out, [tn._dangling[node]], [node, f"qubit {i}", "Close"]) # retag old_ind = tn._dangling[node] @@ -707,13 +715,13 @@ def proj_basis(angle, vop, plane, choice): projected state """ if plane == "XY": - vec = States.vec[0 + choice] + vec = BasicStates.VEC[0 + choice].get_statevector() rotU = Ops.Rz(angle) elif plane == "YZ": - vec = States.vec[4 + choice] + vec = BasicStates.VEC[4 + choice].get_statevector() rotU = Ops.Rx(angle) elif plane == "XZ": - vec = States.vec[0 + choice] + vec = States.VEC[0 + choice].get_statevector() rotU = Ops.Ry(-angle) vec = np.matmul(rotU, vec) vec = np.matmul(CLIFFORD[CLIFFORD_CONJ[vop]], vec) diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index 0bf5e52e2..5b7727ba1 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -3,7 +3,7 @@ import numpy as np -from graphix.ops import States +from graphix.states import BasicStates from graphix.sim.statevec import Statevec, meas_op @@ -29,9 +29,12 @@ def test_measurement_into_each_XYZ_basis(self): n = 3 k = 0 # for measurement into |-> returns [[0, 0], ..., [0, 0]] (whose norm is zero) - for state in [States.plus, States.zero, States.one, States.iplus, States.iminus]: - m_op = np.outer(state, state.T.conjugate()) + # NOTE isn't that weird? + for state in [BasicStates.PLUS, BasicStates.ZERO, BasicStates.ONE, BasicStates.PLUS_I, BasicStates.MINUS_I]: + m_op = np.outer(state.get_statevector(), state.get_statevector().T.conjugate()) + # print(m_op) sv = Statevec(nqubit=n) + # print(sv) sv.evolve(m_op, [k]) sv.remove_qubit(k) @@ -41,7 +44,7 @@ def test_measurement_into_each_XYZ_basis(self): def test_measurement_into_minus_state(self): n = 3 k = 0 - m_op = np.outer(States.minus, States.minus.T.conjugate()) + m_op = np.outer(BasicStates.MINUS.get_statevector(), BasicStates.MINUS.get_statevector().T.conjugate()) sv = Statevec(nqubit=n) sv.evolve(m_op, [k]) with self.assertRaises(AssertionError): diff --git a/tests/test_tnsim.py b/tests/test_tnsim.py index f6eb9a960..7ed7bb0b6 100644 --- a/tests/test_tnsim.py +++ b/tests/test_tnsim.py @@ -6,7 +6,8 @@ import tests.random_circuit as rc from graphix.clifford import CLIFFORD -from graphix.ops import Ops, States +from graphix.ops import Ops +from graphix.states import BasicStates from graphix.sim.tensornet import MBQCTensorNet, gen_str from graphix.transpiler import Circuit @@ -25,7 +26,7 @@ def random_op(sites, dtype=np.complex128, seed=0): CZ = Ops.cz -plus = States.plus +plus = BasicStates.PLUS.get_statevector() class TestTN(unittest.TestCase): From cbedf08c45a8f5d0400b4d985de613ba4e5de105 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Tue, 2 Apr 2024 15:10:58 +0200 Subject: [PATCH 05/48] add tests for statevec + some black --- .gitignore | 6 +++ .vscode/settings.json | 7 ++- graphix/sim/statevec.py | 10 ++-- graphix/simulator.py | 5 +- graphix/states.py | 11 +++- tests/test_statevec.py | 110 ++++++++++++++++++++++++++++++++++++++++ 6 files changed, 141 insertions(+), 8 deletions(-) create mode 100644 tests/test_statevec.py diff --git a/.gitignore b/.gitignore index f8203c191..b9c3a65c6 100644 --- a/.gitignore +++ b/.gitignore @@ -14,3 +14,9 @@ graphix/sim/graphix.code-workspace graphix/graphix.code-workspace *~ .vscode/settings.json +graphix/sim/graphix.code-workspace +.vscode/settings.json +.pre-commit-config.yaml +.vscode/settings.json +.pre-commit-config.yaml +.vscode/settings.json diff --git a/.vscode/settings.json b/.vscode/settings.json index 0a0133979..7a4707700 100644 --- a/.vscode/settings.json +++ b/.vscode/settings.json @@ -7,5 +7,10 @@ "activityBarBadge.background": "#e599d0", "activityBarBadge.foreground": "#15202b" }, - "peacock.color": "#1e5d75" + "peacock.color": "#1e5d75", + "python.testing.pytestArgs": [ + "tests" + ], + "python.testing.unittestEnabled": false, + "python.testing.pytestEnabled": true } \ No newline at end of file diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 00a763a34..efabb63a3 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -204,10 +204,7 @@ def __init__( graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[complex] ] = graphix.states.BasicStates.PLUS, ): - # always infer nqubit from data - # nqubit = none et pas 1 - # if nqubit is None: on exige iterable fini et on prend son nombre d'élément - # also allow external data. + """Initialize statevector Parameters @@ -264,9 +261,11 @@ def __init__( # liste persistante state pour eviter la transience states = [head] + list(it) nqubit = len(states) + self.Nqubit = nqubit # sinon on prend nqubit elts else: # ignore for now states = [head] + [next(it) for _ in range(nqubit - 1)] + self.Nqubit = nqubit list_of_sv = [s.get_statevector() for s in states] tmp_psi = functools.reduce(np.kron, list_of_sv) @@ -295,11 +294,14 @@ def __init__( # just reshape # NOTE too many conversions to numpy arrays? self.psi = psi.reshape((2,) * nqubit) + # for in all cases + self.Nqubit = nqubit # if already a valid statevec just copy it. if isinstance(state, Statevec): assert nqubit is None or len(state.flatten()) == 2**nqubit self.psi = state.psi.copy() + self.Nqubit = state.Nqubit def __repr__(self): return f"Statevec, data={self.psi}, shape={self.dims()}" diff --git a/graphix/simulator.py b/graphix/simulator.py index 2681d2077..616197a4a 100644 --- a/graphix/simulator.py +++ b/graphix/simulator.py @@ -44,7 +44,6 @@ def __init__(self, pattern, backend="statevector", noise_model=None, **kwargs): elif backend == "densitymatrix": if noise_model is None: self.noise_model = None - # no noise: no need to compute probabilities self.backend = DensityMatrixBackend(pattern, **kwargs) warnings.warn( "Simulating using densitymatrix backend with no noise. To add noise to the simulation, give an object of `graphix.noise_models.Noisemodel` to `noise_model` keyword argument." @@ -52,6 +51,7 @@ def __init__(self, pattern, backend="statevector", noise_model=None, **kwargs): if noise_model is not None: self.set_noise_model(noise_model) # if noise: have to compute the probabilities + # NOTE : could remove, pr_calc defaults to True now. self.backend = DensityMatrixBackend(pattern, pr_calc=True, **kwargs) elif backend in {"tensornetwork", "mps"} and noise_model is None: self.noise_model = None @@ -84,7 +84,8 @@ def run(self): the output quantum state, in the representation depending on the backend used. """ - + # use add_nodes or write a new method? + # self.backend.initialize_inputs(self.pattern.input_nodes, option, ...) self.backend.add_nodes(self.pattern.input_nodes) if self.noise_model is None: for cmd in self.pattern: diff --git a/graphix/states.py b/graphix/states.py index edc49bcef..9dfdadfa7 100644 --- a/graphix/states.py +++ b/graphix/states.py @@ -9,6 +9,12 @@ # generic class State for all States class State(abc.ABC): + """Abstract base class for states objects. + Only requirement for concrete classes is to have + a get_statevector() method that returns the statevector + representation of the state + """ + @abc.abstractmethod def get_statevector(self) -> np.ndarray: pass @@ -36,7 +42,7 @@ class PlanarState(pydantic.BaseModel, State): angle: float def __repr__(self): - return f"State objects defined in plane {self.plane} with angle {self.angle}." + return f"PlanarState object defined in plane {self.plane} with angle {self.angle}." def get_statevector(self) -> np.ndarray: if self.plane == graphix.pauli.Plane.XY: @@ -59,6 +65,9 @@ class BasicStates: MINUS = PlanarState(plane=graphix.pauli.Plane.XY, angle=np.pi) PLUS_I = PlanarState(plane=graphix.pauli.Plane.XY, angle=np.pi / 2) MINUS_I = PlanarState(plane=graphix.pauli.Plane.XY, angle=-np.pi / 2) + # remove that in the end + # need in TN backend + VEC = [PLUS, MINUS, ZERO, ONE, PLUS_I, MINUS_I] # Plane.cos.value Plane.cos is an Axis, Axis.value = 0,1,2 (enum) diff --git a/tests/test_statevec.py b/tests/test_statevec.py new file mode 100644 index 000000000..d7223db2a --- /dev/null +++ b/tests/test_statevec.py @@ -0,0 +1,110 @@ +import unittest +import numpy as np + +from graphix.states import BasicStates, PlanarState +import graphix.pauli +import graphix.random_objects as randobj +from graphix.sim.statevec import Statevec +import functools + + +class TestStatevec(unittest.TestCase): + """Test for Statevec class. Particularly new constructor.""" + + def setUp(self): + # set up the random numbers + self.rng = np.random.default_rng() # seed=422 + + # Errors: types, size, + + # test injitializing one qubit in plus state + def test_default_success(self): + vec = Statevec(nqubit=1) + np.testing.assert_allclose(vec.psi, np.array([1, 1] / np.sqrt(2))) + assert vec.Nqubit == 1 + + def test_basicstates_success(self): + # minus + vec = Statevec(nqubit=1, state=BasicStates.MINUS) + np.testing.assert_allclose(vec.psi, np.array([1, -1] / np.sqrt(2))) + assert vec.Nqubit == 1 + # zero + vec = Statevec(nqubit=1, state=BasicStates.ZERO) + np.testing.assert_allclose(vec.psi, np.array([1, 0]), rtol=0, atol=1e-15) + assert vec.Nqubit == 1 + # one + vec = Statevec(nqubit=1, state=BasicStates.ONE) + np.testing.assert_allclose(vec.psi, np.array([0, 1]), rtol=0, atol=1e-15) + assert vec.Nqubit == 1 + # plus_i + vec = Statevec(nqubit=1, state=BasicStates.PLUS_I) + np.testing.assert_allclose(vec.psi, np.array([1, 1j] / np.sqrt(2))) + assert vec.Nqubit == 1 + # minus_i + vec = Statevec(nqubit=1, state=BasicStates.MINUS_I) + np.testing.assert_allclose(vec.psi, np.array([1, -1j] / np.sqrt(2))) + assert vec.Nqubit == 1 + + # even more tests? + def test_default_tensor_success(self): + nqb = self.rng.integers(2, 5) + # print(f"nqb is {nqb}") + vec = Statevec(nqubit=nqb) + np.testing.assert_allclose(vec.psi, np.ones(((2,) * nqb)) / (np.sqrt(2)) ** nqb) + assert vec.Nqubit == nqb + + vec = Statevec(nqubit=nqb, state=BasicStates.MINUS_I) + sv_list = [BasicStates.MINUS_I.get_statevector() for _ in range(nqb)] + sv = functools.reduce(np.kron, sv_list) + np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) + assert vec.Nqubit == nqb + + rand_angle = self.rng.random() * 2 * np.pi + rand_plane = self.rng.choice(np.array([i for i in graphix.pauli.Plane])) + state = PlanarState(plane=rand_plane, angle=rand_angle) + vec = Statevec(nqubit=nqb, state=state) + sv_list = [state.get_statevector() for _ in range(nqb)] + sv = functools.reduce(np.kron, sv_list) + np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) + assert vec.Nqubit == nqb + + def test_data_success(self): + nqb = self.rng.integers(2, 5) + l = 2 ** nqb + rand_vec = self.rng.random(l) + 1j * self.rng.random(l) + rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) + vec = Statevec(state=rand_vec) + np.testing.assert_allclose(vec.psi, rand_vec.reshape((2,) * nqb)) + assert vec.Nqubit == nqb + + # fail: incorrect len + def test_data_dim_fail(self): + l = 5 + rand_vec = self.rng.random(l) + 1j * self.rng.random(l) + rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) + with self.assertRaises(ValueError): + vec = Statevec(state=rand_vec) + + # fail: not normalized + def test_data_norm_fail(self): + nqb = self.rng.integers(2, 5) + l = 2 ** nqb + rand_vec = self.rng.random(l) + 1j * self.rng.random(l) + with self.assertRaises(ValueError): + vec = Statevec(state=rand_vec) + + # fail: no nqubit provided + def test_default_fail(self): + with self.assertRaises(ValueError): + vec = Statevec() + + def test_copy(self): + nqb = self.rng.integers(2, 5) + l = 2 ** nqb + rand_vec = self.rng.random(l) + 1j * self.rng.random(l) + rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) + test_vec = Statevec(state=rand_vec) + vec = Statevec(state=test_vec) + + np.testing.assert_allclose(vec.psi, test_vec.psi) + assert vec.Nqubit == test_vec.Nqubit From 818844d12ba003433d20f340c887901cf85443ea Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Wed, 3 Apr 2024 16:27:13 +0200 Subject: [PATCH 06/48] modifying Backends and simulator, add tests --- graphix/sim/statevec.py | 30 ++++++++++++---- graphix/simulator.py | 2 +- tests/test_statevec.py | 21 ++++++++--- tests/test_statevec_backend.py | 64 +++++++++++++++++++++++++++++++--- 4 files changed, 99 insertions(+), 18 deletions(-) diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index efabb63a3..871c559d7 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -23,7 +23,15 @@ class StatevectorBackend(graphix.sim.base_backend.Backend): """MBQC simulator with statevector method.""" - def __init__(self, pattern, max_qubit_num=20, pr_calc=True): + def __init__( + self, + pattern, + input_state: typing.Union[ + graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[complex] + ] = graphix.states.BasicStates.PLUS, + max_qubit_num=20, + pr_calc=True, + ): """ Parameters ----------- @@ -49,9 +57,12 @@ def __init__(self, pattern, max_qubit_num=20, pr_calc=True): self.to_trace_loc = [] self.max_qubit_num = max_qubit_num if pattern.max_space() > max_qubit_num: - raise ValueError("Pattern.max_space is larger than max_qubit_num. Increase max_qubit_num and try again") + raise ValueError("Pattern.max_space is larger than max_qubit_num. Increase max_qubit_num and try again.") super().__init__(pr_calc) + # initialize input qubits to desired init_state + self.add_nodes(pattern.input_nodes, input_state) + def qubit_dim(self): """Returns the qubit number in the internal statevector @@ -61,7 +72,7 @@ def qubit_dim(self): """ return len(self.state.dims()) - def add_nodes(self, nodes): + def add_nodes(self, nodes, input_state): """add new qubit to internal statevector and assign the corresponding node number to list self.node_index. @@ -73,7 +84,7 @@ def add_nodes(self, nodes): if not self.state: self.state = Statevec(nqubit=0) n = len(nodes) - sv_to_add = Statevec(nqubit=n) + sv_to_add = Statevec(nqubit=n, state=input_state) self.state.tensor(sv_to_add) self.node_index.extend(nodes) self.Nqubit += n @@ -262,7 +273,7 @@ def __init__( states = [head] + list(it) nqubit = len(states) self.Nqubit = nqubit - # sinon on prend nqubit elts + # else take nqubit elts else: # ignore for now states = [head] + [next(it) for _ in range(nqubit - 1)] self.Nqubit = nqubit @@ -437,8 +448,13 @@ def tensor(self, other): """ psi_self = self.psi.flatten() psi_other = other.psi.flatten() - total_num = len(self.dims()) + len(other.dims()) - self.psi = np.kron(psi_self, psi_other).reshape((2,) * total_num) + + # NOTE on tensor form not vector + # deprecated + # total_num = len(self.dims()) + len(other.dims()) + + self.Nqubit += other.Nqubit + self.psi = np.kron(psi_self, psi_other).reshape((2,) * self.Nqubit) def CNOT(self, qubits): """apply CNOT diff --git a/graphix/simulator.py b/graphix/simulator.py index 616197a4a..28bd29257 100644 --- a/graphix/simulator.py +++ b/graphix/simulator.py @@ -86,7 +86,7 @@ def run(self): """ # use add_nodes or write a new method? # self.backend.initialize_inputs(self.pattern.input_nodes, option, ...) - self.backend.add_nodes(self.pattern.input_nodes) + self.backend.add_nodes(self.pattern.input_nodes, state=state) if self.noise_model is None: for cmd in self.pattern: if cmd[0] == "N": diff --git a/tests/test_statevec.py b/tests/test_statevec.py index d7223db2a..749c74841 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -59,6 +59,7 @@ def test_default_tensor_success(self): np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) assert vec.Nqubit == nqb + # tensor of same state rand_angle = self.rng.random() * 2 * np.pi rand_plane = self.rng.choice(np.array([i for i in graphix.pauli.Plane])) state = PlanarState(plane=rand_plane, angle=rand_angle) @@ -68,12 +69,22 @@ def test_default_tensor_success(self): np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) assert vec.Nqubit == nqb + # tensor of different states + rand_angles = self.rng.random(nqb) * 2 * np.pi + rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) + states = [PlanarState(plane = i, angle = j) for i, j in zip(rand_planes, rand_angles)] + vec = Statevec(nqubit=nqb, state=states) + sv_list = [state.get_statevector() for state in states] + sv = functools.reduce(np.kron, sv_list) + np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) + assert vec.Nqubit == nqb + def test_data_success(self): nqb = self.rng.integers(2, 5) l = 2 ** nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) - vec = Statevec(state=rand_vec) + vec = Statevec(state = rand_vec) np.testing.assert_allclose(vec.psi, rand_vec.reshape((2,) * nqb)) assert vec.Nqubit == nqb @@ -83,7 +94,7 @@ def test_data_dim_fail(self): rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) with self.assertRaises(ValueError): - vec = Statevec(state=rand_vec) + vec = Statevec(state = rand_vec) # fail: not normalized def test_data_norm_fail(self): @@ -91,7 +102,7 @@ def test_data_norm_fail(self): l = 2 ** nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) with self.assertRaises(ValueError): - vec = Statevec(state=rand_vec) + vec = Statevec(state = rand_vec) # fail: no nqubit provided def test_default_fail(self): @@ -103,8 +114,8 @@ def test_copy(self): l = 2 ** nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) - test_vec = Statevec(state=rand_vec) - vec = Statevec(state=test_vec) + test_vec = Statevec(state = rand_vec) + vec = Statevec(state = test_vec) np.testing.assert_allclose(vec.psi, test_vec.psi) assert vec.Nqubit == test_vec.Nqubit diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index 5b7727ba1..5a5375e7b 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -2,10 +2,10 @@ from copy import deepcopy import numpy as np - -from graphix.states import BasicStates -from graphix.sim.statevec import Statevec, meas_op - +from graphix import Circuit +from graphix.states import BasicStates, PlanarState +from graphix.sim.statevec import Statevec, meas_op, StatevectorBackend +import graphix.pauli class TestStatevec(unittest.TestCase): def test_remove_one_qubit(self): @@ -25,11 +25,12 @@ def test_remove_one_qubit(self): np.testing.assert_almost_equal(np.abs(sv.psi.flatten().dot(sv2.psi.flatten().conj())), 1) + #TODO This is a weird test! def test_measurement_into_each_XYZ_basis(self): n = 3 k = 0 # for measurement into |-> returns [[0, 0], ..., [0, 0]] (whose norm is zero) - # NOTE isn't that weird? + # NOTE weird choice (MINUS is orthogonal to PLUS so zero) for state in [BasicStates.PLUS, BasicStates.ZERO, BasicStates.ONE, BasicStates.PLUS_I, BasicStates.MINUS_I]: m_op = np.outer(state.get_statevector(), state.get_statevector().T.conjugate()) # print(m_op) @@ -50,6 +51,59 @@ def test_measurement_into_minus_state(self): with self.assertRaises(AssertionError): sv.remove_qubit(k) +class TestStatevecNew(unittest.TestCase): + def setUp(self): + # set up the random numbers + self.rng = np.random.default_rng() # seed=422 + + circ = Circuit(1) + circ.h(0) + self.hadamardpattern = circ.transpile() + + # test initialization only + def test_init_success(self): + + # plus state (default) + backend = StatevectorBackend(self.hadamardpattern) + vec = Statevec(nqubit=1) + np.testing.assert_allclose(vec.psi, backend.state.psi) + assert backend.state.Nqubit == 1 + + # minus state + backend = StatevectorBackend(self.hadamardpattern, input_state = BasicStates.MINUS) + vec = Statevec(nqubit=1, state = BasicStates.MINUS) + np.testing.assert_allclose(vec.psi, backend.state.psi) + assert backend.state.Nqubit == 1 + + # random planar state + rand_angle = self.rng.random() * 2 * np.pi + rand_plane = self.rng.choice(np.array([i for i in graphix.pauli.Plane])) + state = PlanarState(plane = rand_plane, angle = rand_angle) + backend = StatevectorBackend(self.hadamardpattern, input_state = state) + vec = Statevec(nqubit=1, state = state) + np.testing.assert_allclose(vec.psi, backend.state.psi) + assert backend.state.Nqubit == 1 + + # incorrect number of dimensions + # only one input node,, two states provided + # doesn't fail! just takes the first qubit! + # Discard second qubit so can be whatever + + rand_angle = self.rng.random(2) * 2 * np.pi + rand_plane = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), 2) + + state = PlanarState(plane = rand_plane[0], angle = rand_angle[0]) + state2 = PlanarState(plane = rand_plane[1], angle = rand_angle[1]) + #with self.assertRaises(ValueError): + backend = StatevectorBackend(self.hadamardpattern, input_state = [state, state2]) + vec = Statevec(nqubit=1, state = state) + np.testing.assert_allclose(vec.psi, backend.state.psi) + assert backend.state.Nqubit == 1 + + + + + if __name__ == "__main__": unittest.main() From 46bdfa53aa9fa6e0af61cd5610a9fadc2346d6d0 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Wed, 3 Apr 2024 16:59:52 +0200 Subject: [PATCH 07/48] updates --- graphix/sim/statevec.py | 2 +- graphix/simulator.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 871c559d7..aa7975266 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -72,7 +72,7 @@ def qubit_dim(self): """ return len(self.state.dims()) - def add_nodes(self, nodes, input_state): + def add_nodes(self, nodes, input_state=graphix.states.BasicStates.PLUS): """add new qubit to internal statevector and assign the corresponding node number to list self.node_index. diff --git a/graphix/simulator.py b/graphix/simulator.py index 28bd29257..d606d9a60 100644 --- a/graphix/simulator.py +++ b/graphix/simulator.py @@ -86,7 +86,7 @@ def run(self): """ # use add_nodes or write a new method? # self.backend.initialize_inputs(self.pattern.input_nodes, option, ...) - self.backend.add_nodes(self.pattern.input_nodes, state=state) + # self.backend.add_nodes(self.pattern.input_nodes, input_state=state) if self.noise_model is None: for cmd in self.pattern: if cmd[0] == "N": From 86ccb49e446c346672dcb9c1ab627897f0c52fbd Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Wed, 3 Apr 2024 17:09:24 +0200 Subject: [PATCH 08/48] cancelled modification in tensor. Indeed, the attribute Statevec.Nqubit leads to problems that have to be fixed --- graphix/sim/statevec.py | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index aa7975266..4e0bcc3ef 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -239,6 +239,7 @@ def __init__( if nqubit == 0: warnings.warn(f"Called Statevec with 0 qubits. Ignoring the state.") self.psi = np.array(1, dtype=np.complex128) + self.Nqubit = 0 # works only for planar states. Deal with all kind of states? elif isinstance(state, graphix.states.State): @@ -451,10 +452,11 @@ def tensor(self, other): # NOTE on tensor form not vector # deprecated - # total_num = len(self.dims()) + len(other.dims()) - - self.Nqubit += other.Nqubit - self.psi = np.kron(psi_self, psi_other).reshape((2,) * self.Nqubit) + total_num = len(self.dims()) + len(other.dims()) + print(total_num) + # self.Nqubit += other.Nqubit + self.psi = np.kron(psi_self, psi_other).reshape((2,) * total_num) + self.Nqubit = len(self.dims()) def CNOT(self, qubits): """apply CNOT From c2487a4c48f9e6a2f3d68cc889863a50142227c4 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Wed, 3 Apr 2024 17:13:04 +0200 Subject: [PATCH 09/48] Update statevec.py --- graphix/sim/statevec.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 4e0bcc3ef..e2349629e 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -87,7 +87,7 @@ def add_nodes(self, nodes, input_state=graphix.states.BasicStates.PLUS): sv_to_add = Statevec(nqubit=n, state=input_state) self.state.tensor(sv_to_add) self.node_index.extend(nodes) - self.Nqubit += n + self.Nqubit += self.state.Nqubit # n def entangle_nodes(self, edge): """Apply CZ gate to two connected nodes From 882cac7ca5bb07882be96154008562209edfdc7a Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Wed, 3 Apr 2024 19:43:59 +0200 Subject: [PATCH 10/48] commented all statevec.Nqubit attributes some test fixed, some other not (related to DM and TN backend) --- graphix/sim/statevec.py | 14 ++++++------ tests/test_statevec.py | 41 ++++++++++++++++++++++++---------- tests/test_statevec_backend.py | 12 ++++++---- 3 files changed, 44 insertions(+), 23 deletions(-) diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index e2349629e..c9eaddb52 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -87,7 +87,7 @@ def add_nodes(self, nodes, input_state=graphix.states.BasicStates.PLUS): sv_to_add = Statevec(nqubit=n, state=input_state) self.state.tensor(sv_to_add) self.node_index.extend(nodes) - self.Nqubit += self.state.Nqubit # n + self.Nqubit += n def entangle_nodes(self, edge): """Apply CZ gate to two connected nodes @@ -239,7 +239,7 @@ def __init__( if nqubit == 0: warnings.warn(f"Called Statevec with 0 qubits. Ignoring the state.") self.psi = np.array(1, dtype=np.complex128) - self.Nqubit = 0 + # self.Nqubit = 0 # works only for planar states. Deal with all kind of states? elif isinstance(state, graphix.states.State): @@ -273,11 +273,11 @@ def __init__( # liste persistante state pour eviter la transience states = [head] + list(it) nqubit = len(states) - self.Nqubit = nqubit + # self.Nqubit = nqubit # else take nqubit elts else: # ignore for now states = [head] + [next(it) for _ in range(nqubit - 1)] - self.Nqubit = nqubit + # self.Nqubit = nqubit list_of_sv = [s.get_statevector() for s in states] tmp_psi = functools.reduce(np.kron, list_of_sv) @@ -307,13 +307,13 @@ def __init__( # NOTE too many conversions to numpy arrays? self.psi = psi.reshape((2,) * nqubit) # for in all cases - self.Nqubit = nqubit + # self.Nqubit = nqubit # if already a valid statevec just copy it. if isinstance(state, Statevec): assert nqubit is None or len(state.flatten()) == 2**nqubit self.psi = state.psi.copy() - self.Nqubit = state.Nqubit + # self.Nqubit = state.Nqubit def __repr__(self): return f"Statevec, data={self.psi}, shape={self.dims()}" @@ -456,7 +456,7 @@ def tensor(self, other): print(total_num) # self.Nqubit += other.Nqubit self.psi = np.kron(psi_self, psi_other).reshape((2,) * total_num) - self.Nqubit = len(self.dims()) + # self.Nqubit = len(self.dims()) def CNOT(self, qubits): """apply CNOT diff --git a/tests/test_statevec.py b/tests/test_statevec.py index 749c74841..d2bfcfaac 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -21,29 +21,39 @@ def setUp(self): def test_default_success(self): vec = Statevec(nqubit=1) np.testing.assert_allclose(vec.psi, np.array([1, 1] / np.sqrt(2))) - assert vec.Nqubit == 1 + # assert vec.Nqubit == 1 + assert len(vec.dims()) == 1 def test_basicstates_success(self): # minus vec = Statevec(nqubit=1, state=BasicStates.MINUS) np.testing.assert_allclose(vec.psi, np.array([1, -1] / np.sqrt(2))) - assert vec.Nqubit == 1 + # assert vec.Nqubit == 1 + assert len(vec.dims()) == 1 + # zero vec = Statevec(nqubit=1, state=BasicStates.ZERO) np.testing.assert_allclose(vec.psi, np.array([1, 0]), rtol=0, atol=1e-15) - assert vec.Nqubit == 1 + # assert vec.Nqubit == 1 + assert len(vec.dims()) == 1 + # one vec = Statevec(nqubit=1, state=BasicStates.ONE) np.testing.assert_allclose(vec.psi, np.array([0, 1]), rtol=0, atol=1e-15) - assert vec.Nqubit == 1 + # assert vec.Nqubit == 1 + assert len(vec.dims()) == 1 + # plus_i vec = Statevec(nqubit=1, state=BasicStates.PLUS_I) np.testing.assert_allclose(vec.psi, np.array([1, 1j] / np.sqrt(2))) - assert vec.Nqubit == 1 + # assert vec.Nqubit == 1 + assert len(vec.dims()) == 1 + # minus_i vec = Statevec(nqubit=1, state=BasicStates.MINUS_I) np.testing.assert_allclose(vec.psi, np.array([1, -1j] / np.sqrt(2))) - assert vec.Nqubit == 1 + #assert vec.Nqubit == 1 + assert len(vec.dims()) == 1 # even more tests? def test_default_tensor_success(self): @@ -51,13 +61,15 @@ def test_default_tensor_success(self): # print(f"nqb is {nqb}") vec = Statevec(nqubit=nqb) np.testing.assert_allclose(vec.psi, np.ones(((2,) * nqb)) / (np.sqrt(2)) ** nqb) - assert vec.Nqubit == nqb + # assert vec.Nqubit == nqb + assert len(vec.dims()) == nqb vec = Statevec(nqubit=nqb, state=BasicStates.MINUS_I) sv_list = [BasicStates.MINUS_I.get_statevector() for _ in range(nqb)] sv = functools.reduce(np.kron, sv_list) np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) - assert vec.Nqubit == nqb + # assert vec.Nqubit == nqb + assert len(vec.dims()) == nqb # tensor of same state rand_angle = self.rng.random() * 2 * np.pi @@ -67,7 +79,8 @@ def test_default_tensor_success(self): sv_list = [state.get_statevector() for _ in range(nqb)] sv = functools.reduce(np.kron, sv_list) np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) - assert vec.Nqubit == nqb + # assert vec.Nqubit == nqb + assert len(vec.dims()) == nqb # tensor of different states rand_angles = self.rng.random(nqb) * 2 * np.pi @@ -77,7 +90,8 @@ def test_default_tensor_success(self): sv_list = [state.get_statevector() for state in states] sv = functools.reduce(np.kron, sv_list) np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) - assert vec.Nqubit == nqb + # assert vec.Nqubit == nqb + assert len(vec.dims()) == nqb def test_data_success(self): nqb = self.rng.integers(2, 5) @@ -86,7 +100,9 @@ def test_data_success(self): rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) vec = Statevec(state = rand_vec) np.testing.assert_allclose(vec.psi, rand_vec.reshape((2,) * nqb)) - assert vec.Nqubit == nqb + # assert vec.Nqubit == nqb + assert len(vec.dims()) == nqb + # fail: incorrect len def test_data_dim_fail(self): @@ -118,4 +134,5 @@ def test_copy(self): vec = Statevec(state = test_vec) np.testing.assert_allclose(vec.psi, test_vec.psi) - assert vec.Nqubit == test_vec.Nqubit + # assert vec.Nqubit == test_vec.Nqubit + assert len(vec.dims()) == len(test_vec.dims()) diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index 5a5375e7b..38e807cf8 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -67,13 +67,15 @@ def test_init_success(self): backend = StatevectorBackend(self.hadamardpattern) vec = Statevec(nqubit=1) np.testing.assert_allclose(vec.psi, backend.state.psi) - assert backend.state.Nqubit == 1 + # assert backend.state.Nqubit == 1 + assert len(backend.state.dims()) == 1 # minus state backend = StatevectorBackend(self.hadamardpattern, input_state = BasicStates.MINUS) vec = Statevec(nqubit=1, state = BasicStates.MINUS) np.testing.assert_allclose(vec.psi, backend.state.psi) - assert backend.state.Nqubit == 1 + # assert backend.state.Nqubit == 1 + assert len(backend.state.dims()) == 1 # random planar state rand_angle = self.rng.random() * 2 * np.pi @@ -82,7 +84,8 @@ def test_init_success(self): backend = StatevectorBackend(self.hadamardpattern, input_state = state) vec = Statevec(nqubit=1, state = state) np.testing.assert_allclose(vec.psi, backend.state.psi) - assert backend.state.Nqubit == 1 + # assert backend.state.Nqubit == 1 + assert len(backend.state.dims()) == 1 # incorrect number of dimensions # only one input node,, two states provided @@ -98,7 +101,8 @@ def test_init_success(self): backend = StatevectorBackend(self.hadamardpattern, input_state = [state, state2]) vec = Statevec(nqubit=1, state = state) np.testing.assert_allclose(vec.psi, backend.state.psi) - assert backend.state.Nqubit == 1 + # assert backend.state.Nqubit == 1 + assert len(backend.state.dims()) == 1 From 341536927f9f4f2eaf26d727814053f38d687906 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Thu, 4 Apr 2024 18:33:25 +0200 Subject: [PATCH 11/48] update statevec.__init__. DM ongoing --- graphix/sim/density_matrix.py | 151 +++++++++++++++++++++++++++------ graphix/sim/statevec.py | 98 ++++++++------------- tests/test_density_matrix.py | 3 +- tests/test_statevec.py | 33 ++++++- tests/test_statevec_backend.py | 20 +++-- 5 files changed, 201 insertions(+), 104 deletions(-) diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index e139c88ec..3cacae401 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -7,19 +7,40 @@ import numpy as np -import graphix.sim.base_backend +from graphix.linalg_validations import check_square, check_hermitian, check_unit_trace from graphix.channels import KrausChannel -from graphix.clifford import CLIFFORD -from graphix.linalg_validations import check_hermitian, check_square, check_unit_trace from graphix.ops import Ops -from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, meas_op +from graphix.clifford import CLIFFORD +from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, meas_op, Statevec +import graphix.sim.base_backend +import graphix.states + +import typing +from typing_extensions import Annotated +from annotated_types import Ge +import functools + + +PositiveInt = Annotated[int, Ge(0)] # includes 0 class DensityMatrix: """DensityMatrix object.""" - def __init__(self, data=None, plus_state=True, nqubit=1): + def __init__( + self, + nqubit: typing.Optional[PositiveInt] = None, + state: typing.Union[ + graphix.states.State, + "DensityMatrix", + Statevec, + typing.Iterable[graphix.states.State], + typing.Iterable[complex], + ] = graphix.states.BasicStates.PLUS, + ): + """ + rewrite! Parameters ---------- data : DensityMatrix, list, tuple, np.ndarray or None @@ -27,33 +48,107 @@ def __init__(self, data=None, plus_state=True, nqubit=1): nqubit : int Number of qubits. Default is 1. If both `data` and `nqubit` are specified, `nqubit` is ignored. """ - if data is None: - assert nqubit >= 0 - self.Nqubit = nqubit - if plus_state: - self.rho = np.ones((2**nqubit, 2**nqubit)) / 2**nqubit - else: - self.rho = np.zeros((2**nqubit, 2**nqubit)) - self.rho[0, 0] = 1.0 - else: - if isinstance(data, DensityMatrix): - data = data.rho - elif isinstance(data, (list, tuple)): - data = np.asarray(data, dtype=complex) - elif isinstance(data, np.ndarray): - pass - else: - raise TypeError("data must be DensityMatrix, list, tuple, or np.ndarray.") + if nqubit == 0: + self.rho = np.array(1, dtype=np.complex128) + + elif isinstance(state, graphix.states.State): + + if nqubit is None: + raise ValueError("Incorrect value for nqubit.") - assert check_square(data) - self.Nqubit = len(data).bit_length() - 1 + # or directly get_dm from get_statevec? Can make it inherit from abc? + vec = state.get_statevector() + dm = np.outer(vec, vec.conj()) + # these vecs are normalized + if nqubit == 1: # or None + self.rho = dm + + # build tensor product |state>^{\otimes nqubit} + # can only be >1 int. + else: + # build tensor product + # comma in tuple is for disambiguation with paranthesed expression + self.rho = functools.reduce(np.kron, (dm,) * nqubit) + # no reshape + + # nqubit is None : on prend la longeur de l'iterable + elif isinstance(state, typing.Iterable): + # iterateur + it = iter(state) + head = next(it) + # type constraint in head doesn't progpagate to all elts + if isinstance(head, graphix.states.State): + # assert isinstance(head, typing.Iterator[graphix.states.State]) + if nqubit is None: + # liste persistante state pour eviter la transience + states = [head] + list(it) + nqubit = len(states) + # self.Nqubit = nqubit + # else take nqubit elts + else: # ignore for now + states = [head] + [next(it) for _ in range(nqubit - 1)] + + # self.Nqubit = nqubit + + list_of_sv = [s.get_statevector() for s in states] + list_of_dm = [np.outer(sv, sv.conj()) for sv in list_of_sv] + self.rho = functools.reduce(np.kron, list_of_dm) + # no reshape + # self.psi = tmp_psi.reshape((2,) * nqubit) + + else: # now 2**n by 2**n matrices also iterable? + if nqubit is None: + # BUG what to do with that? + states = [head] + list(it) + # need a shape so just np.ndarray? + inferred_shape = state.shape + + # this is done in check_square for matrix do it on shape + if len(inferred_shape) != 2: + raise ValueError(f"The object has {len(inferred_shape)} axes but must have 2 to be a matrix.") + if inferred_shape[0] != inferred_shape[1]: + raise ValueError(f"Matrix must be square but has different dimensions {inferred_shape}.") + inferred_size = inferred_shape[0] + if inferred_size & (inferred_size - 1) != 0: + raise ValueError(f"Matrix size must be a power of two but is {inferred_size}.") + + nqubit = inferred_size.bit_length() - 1 + + else: # ignore for now + # BUG what to do with that? + states = state[: 2**nqubit, : 2**nqubit].copy() + # [head] + [next(it) for _ in range(2**nqubit - 1)] + + # psi = np.array(states) + + dm = np.array([np.outer(s, s.conj()) for s in states]) + + # hermicity and trace checked later. Or do it from statevec? + # from pure states, should be ok by default so just when input data? + assert check_hermitian(self.rho) + assert check_unit_trace(self.rho) + + # TODO now: statevec and densitymatrix + # if already a valid statevec transform it to DensityMatrix. + if isinstance(state, Statevec): + if nqubit is not None or len(state.flatten()) != 2**nqubit: + raise ValueError( + f"Inconsistent parameters between nqubit = {nqubit} and the size of the provided statevector = {len(state.flatten())}." + ) + vect = state.psi.copy() + self.rho = np.outer(vect, vect.conj()) + + # if DensityMatrix, just copy it + elif isinstance(state, DensityMatrix): + if nqubit is not None or state.dims() != (2**nqubit, 2**nqubit): + raise ValueError( + f"Inconsistent parameters between nqubit = {nqubit} and the shape of the provided density matrix = {state.dims()}." + ) - self.rho = data - assert check_hermitian(self.rho) - assert check_unit_trace(self.rho) + self.rho = state.rho.copy() def __repr__(self): - return f"DensityMatrix, data={self.rho}, shape={self.dims()}" + return f"DensityMatrix object , with density matrix {self.rho} and shape{self.dims()}." def evolve_single(self, op, i): """Single-qubit operation. diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index c9eaddb52..ed3844ad7 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -236,87 +236,59 @@ def __init__( # NOTE make plane, angle attributes of Statevec? # this will be called by the nqb > 1 case # instantiate a one - if nqubit == 0: - warnings.warn(f"Called Statevec with 0 qubits. Ignoring the state.") - self.psi = np.array(1, dtype=np.complex128) - # self.Nqubit = 0 - - # works only for planar states. Deal with all kind of states? - elif isinstance(state, graphix.states.State): + if isinstance(state, Statevec): + # assert nqubit is None or len(state.flatten()) == 2**nqubit + if nqubit is not None and len(state.flatten()) != 2**nqubit: + raise ValueError( + f"Inconsistent parameters between nqubit = {nqubit} and the inferred number of qubit = {len(state.flatten())}." + ) + self.psi = state.psi.copy() + return + if isinstance(state, graphix.states.State): if nqubit is None: raise ValueError("Incorrect value for nqubit.") + input_list = [state] * nqubit + elif isinstance(state, typing.Iterable): + input_list = list(state) + else: + raise TypeError("Incorrect type for input state") - vec = state.get_statevector() - - if nqubit == 1: # or None - self.psi = vec - - # build tensor product |state>^{\otimes nqubit} - # can only be >1 int. - else: - # build tensor product - # comma in tuple is for disambiguation with paranthesed expression - tmp_psi = functools.reduce(np.kron, (vec,) * nqubit) - # reshape - self.psi = tmp_psi.reshape((2,) * nqubit) + if len(input_list) == 0: + if nqubit is not None and nqubit != 0: + raise ValueError("nqubit is not null but input state is empty.") - # nqubit is None : on prend la longeur de l'iterable - elif isinstance(state, typing.Iterable): - # iterateur - it = iter(state) - head = next(it) - # type constraint in head doesn't progpagate to all elts - if isinstance(head, graphix.states.State): - # assert isinstance(head, typing.Iterator[graphix.states.State]) + # warnings.warn(f"Called Statevec with 0 qubits. Ignoring the state.") + self.psi = np.array(1, dtype=np.complex128) + # self.Nqubit = 0 + else: + if isinstance(input_list[0], graphix.states.State): if nqubit is None: - # liste persistante state pour eviter la transience - states = [head] + list(it) - nqubit = len(states) - # self.Nqubit = nqubit - # else take nqubit elts - else: # ignore for now - states = [head] + [next(it) for _ in range(nqubit - 1)] - # self.Nqubit = nqubit - - list_of_sv = [s.get_statevector() for s in states] + nqubit = len(input_list) + elif nqubit != len(input_list): + raise ValueError("Mismatch between nqubit and length of input state") + list_of_sv = [s.get_statevector() for s in input_list] tmp_psi = functools.reduce(np.kron, list_of_sv) # reshape self.psi = tmp_psi.reshape((2,) * nqubit) - else: if nqubit is None: - states = [head] + list(it) - - inferred_size = len(states) - - if inferred_size & (inferred_size - 1) != 0: - raise ValueError(f"Statevector size must be a power of two but is {inferred_size}.") - - nqubit = inferred_size.bit_length() - 1 - - else: # ignore for now - states = [head] + [next(it) for _ in range(2**nqubit - 1)] - - psi = np.array(states) - + length = len(input_list) + if length & (length - 1): + raise ValueError("Length is not a power of two") + nqubit = length.bit_length() - 1 + elif nqubit != len(input_list).bit_length() - 1: + raise ValueError("Mismatch between nqubit and length of input state") + psi = np.array(input_list) if not np.allclose(np.sqrt(np.sum(np.abs(psi) ** 2)), 1): - raise ValueError(f"Statevector must be normalized to one.") - + raise ValueError("Input state is not normalized") # just reshape # NOTE too many conversions to numpy arrays? self.psi = psi.reshape((2,) * nqubit) - # for in all cases - # self.Nqubit = nqubit - - # if already a valid statevec just copy it. - if isinstance(state, Statevec): - assert nqubit is None or len(state.flatten()) == 2**nqubit - self.psi = state.psi.copy() # self.Nqubit = state.Nqubit def __repr__(self): - return f"Statevec, data={self.psi}, shape={self.dims()}" + return f"Statevec object with statevector {self.psi} and length {self.dims()}." def evolve_single(self, op, i): """Single-qubit operation diff --git a/tests/test_density_matrix.py b/tests/test_density_matrix.py index 22df7fad6..e9bb8bb84 100644 --- a/tests/test_density_matrix.py +++ b/tests/test_density_matrix.py @@ -863,7 +863,8 @@ def test_correct_byproduct(self): rho = backend.state.rho backend = StatevectorBackend(pattern) - backend.add_nodes([0, 1, 2]) + # node 0 initialized in Backend + backend.add_nodes([1, 2]) backend.entangle_nodes((0, 1)) backend.entangle_nodes((1, 2)) backend.measure(backend.pattern[-4]) diff --git a/tests/test_statevec.py b/tests/test_statevec.py index d2bfcfaac..d3ac543e6 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -4,10 +4,12 @@ from graphix.states import BasicStates, PlanarState import graphix.pauli import graphix.random_objects as randobj -from graphix.sim.statevec import Statevec +from graphix.sim.statevec import Statevec, NotNormalizedError import functools + + class TestStatevec(unittest.TestCase): """Test for Statevec class. Particularly new constructor.""" @@ -58,7 +60,7 @@ def test_basicstates_success(self): # even more tests? def test_default_tensor_success(self): nqb = self.rng.integers(2, 5) - # print(f"nqb is {nqb}") + print(f"nqb is {nqb}") vec = Statevec(nqubit=nqb) np.testing.assert_allclose(vec.psi, np.ones(((2,) * nqb)) / (np.sqrt(2)) ** nqb) # assert vec.Nqubit == nqb @@ -112,6 +114,16 @@ def test_data_dim_fail(self): with self.assertRaises(ValueError): vec = Statevec(state = rand_vec) + # with less qubit than number of qubits inferred from a correct state vect + # returns a truncated statevec that is hence not normalized + # NOTE weird behaviour?? + def test_data_dim_fail_mismatch(self): + nqb = 3 + rand_vec = self.rng.random(2 ** nqb) + 1j * self.rng.random(2 ** nqb) + rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) + with self.assertRaises(ValueError): + vec = Statevec(nqubit = 2, state = rand_vec) + # fail: not normalized def test_data_norm_fail(self): nqb = self.rng.integers(2, 5) @@ -120,19 +132,32 @@ def test_data_norm_fail(self): with self.assertRaises(ValueError): vec = Statevec(state = rand_vec) - # fail: no nqubit provided + # fail: no nqubit provided. State defaults to PLUS def test_default_fail(self): with self.assertRaises(ValueError): vec = Statevec() - def test_copy(self): + # try copying Statevec input + def test_copy_success(self): nqb = self.rng.integers(2, 5) l = 2 ** nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) test_vec = Statevec(state = rand_vec) + # try to copy it vec = Statevec(state = test_vec) np.testing.assert_allclose(vec.psi, test_vec.psi) # assert vec.Nqubit == test_vec.Nqubit assert len(vec.dims()) == len(test_vec.dims()) + + # try calling with incorrect number of qubits compared to inferred one + def test_copy_fail(self): + nqb = self.rng.integers(2, 5) + l = 2 ** nqb + rand_vec = self.rng.random(l) + 1j * self.rng.random(l) + rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) + test_vec = Statevec(state = rand_vec) + + with self.assertRaises(ValueError): + vec = Statevec(nqubit = l - 1, state = test_vec) diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index 38e807cf8..3f45d1bf6 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -87,8 +87,12 @@ def test_init_success(self): # assert backend.state.Nqubit == 1 assert len(backend.state.dims()) == 1 - # incorrect number of dimensions - # only one input node,, two states provided + + + + def test_init_fail(self): + # incorrect number of dimensions for State input + # only one input node, two states provided # doesn't fail! just takes the first qubit! # Discard second qubit so can be whatever @@ -97,12 +101,12 @@ def test_init_success(self): state = PlanarState(plane = rand_plane[0], angle = rand_angle[0]) state2 = PlanarState(plane = rand_plane[1], angle = rand_angle[1]) - #with self.assertRaises(ValueError): - backend = StatevectorBackend(self.hadamardpattern, input_state = [state, state2]) - vec = Statevec(nqubit=1, state = state) - np.testing.assert_allclose(vec.psi, backend.state.psi) - # assert backend.state.Nqubit == 1 - assert len(backend.state.dims()) == 1 + with self.assertRaises(ValueError): + StatevectorBackend(self.hadamardpattern, input_state = [state, state2]) + # vec = Statevec(nqubit=1, state = state) + # np.testing.assert_allclose(vec.psi, backend.state.psi) + # # assert backend.state.Nqubit == 1 + # assert len(backend.state.dims()) == 1 From ab6781dfdde5ac9d22714781ba49f2fe3f62bde3 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Thu, 4 Apr 2024 19:12:32 +0200 Subject: [PATCH 12/48] Constructors for statevec and DM --- graphix/sim/density_matrix.py | 120 +++++++--------------------------- tests/test_density_matrix.py | 22 ++----- tests/test_statevec.py | 2 +- 3 files changed, 32 insertions(+), 112 deletions(-) diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index 3cacae401..b9317456b 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -7,7 +7,7 @@ import numpy as np -from graphix.linalg_validations import check_square, check_hermitian, check_unit_trace +from graphix.linalg_validations import check_square, check_hermitian, check_unit_trace, check_psd from graphix.channels import KrausChannel from graphix.ops import Ops from graphix.clifford import CLIFFORD @@ -30,12 +30,13 @@ class DensityMatrix: def __init__( self, nqubit: typing.Optional[PositiveInt] = None, - state: typing.Union[ + data: typing.Union[ graphix.states.State, "DensityMatrix", Statevec, typing.Iterable[graphix.states.State], typing.Iterable[complex], + typing.Iterable[typing.Iterable[complex]], ] = graphix.states.BasicStates.PLUS, ): @@ -48,104 +49,31 @@ def __init__( nqubit : int Number of qubits. Default is 1. If both `data` and `nqubit` are specified, `nqubit` is ignored. """ - if nqubit == 0: - self.rho = np.array(1, dtype=np.complex128) - elif isinstance(state, graphix.states.State): - - if nqubit is None: - raise ValueError("Incorrect value for nqubit.") - - # or directly get_dm from get_statevec? Can make it inherit from abc? - vec = state.get_statevector() - dm = np.outer(vec, vec.conj()) - # these vecs are normalized - if nqubit == 1: # or None - self.rho = dm - - # build tensor product |state>^{\otimes nqubit} - # can only be >1 int. - else: - # build tensor product - # comma in tuple is for disambiguation with paranthesed expression - self.rho = functools.reduce(np.kron, (dm,) * nqubit) - # no reshape - - # nqubit is None : on prend la longeur de l'iterable - elif isinstance(state, typing.Iterable): - # iterateur - it = iter(state) - head = next(it) - # type constraint in head doesn't progpagate to all elts - if isinstance(head, graphix.states.State): - # assert isinstance(head, typing.Iterator[graphix.states.State]) - if nqubit is None: - # liste persistante state pour eviter la transience - states = [head] + list(it) - nqubit = len(states) - # self.Nqubit = nqubit - # else take nqubit elts - else: # ignore for now - states = [head] + [next(it) for _ in range(nqubit - 1)] - - # self.Nqubit = nqubit - - list_of_sv = [s.get_statevector() for s in states] - list_of_dm = [np.outer(sv, sv.conj()) for sv in list_of_sv] - self.rho = functools.reduce(np.kron, list_of_dm) - # no reshape - # self.psi = tmp_psi.reshape((2,) * nqubit) - - else: # now 2**n by 2**n matrices also iterable? - if nqubit is None: - # BUG what to do with that? - states = [head] + list(it) - # need a shape so just np.ndarray? - inferred_shape = state.shape - - # this is done in check_square for matrix do it on shape - if len(inferred_shape) != 2: - raise ValueError(f"The object has {len(inferred_shape)} axes but must have 2 to be a matrix.") - if inferred_shape[0] != inferred_shape[1]: - raise ValueError(f"Matrix must be square but has different dimensions {inferred_shape}.") - inferred_size = inferred_shape[0] - if inferred_size & (inferred_size - 1) != 0: - raise ValueError(f"Matrix size must be a power of two but is {inferred_size}.") - - nqubit = inferred_size.bit_length() - 1 - - else: # ignore for now - # BUG what to do with that? - states = state[: 2**nqubit, : 2**nqubit].copy() - # [head] + [next(it) for _ in range(2**nqubit - 1)] - - # psi = np.array(states) - - dm = np.array([np.outer(s, s.conj()) for s in states]) - - # hermicity and trace checked later. Or do it from statevec? - # from pure states, should be ok by default so just when input data? - assert check_hermitian(self.rho) - assert check_unit_trace(self.rho) - - # TODO now: statevec and densitymatrix - # if already a valid statevec transform it to DensityMatrix. - if isinstance(state, Statevec): - if nqubit is not None or len(state.flatten()) != 2**nqubit: - raise ValueError( - f"Inconsistent parameters between nqubit = {nqubit} and the size of the provided statevector = {len(state.flatten())}." - ) - vect = state.psi.copy() - self.rho = np.outer(vect, vect.conj()) - - # if DensityMatrix, just copy it - elif isinstance(state, DensityMatrix): - if nqubit is not None or state.dims() != (2**nqubit, 2**nqubit): + def check_size_consistency(mat): + if nqubit is not None and mat.shape != (2**nqubit, 2**nqubit): raise ValueError( - f"Inconsistent parameters between nqubit = {nqubit} and the shape of the provided density matrix = {state.dims()}." + f"Inconsistent parameters between nqubit = {nqubit} and the shape of the provided density matrix = {mat.shape}." ) - self.rho = state.rho.copy() + if isinstance(data, DensityMatrix): + check_size_consistency(data) + self.rho = data.rho.copy() + self.Nqubit = data.Nqubit + return + if isinstance(data, typing.Iterable): + input_list = list(data) + if len(input_list) != 0: + if isinstance(input_list[0], typing.Iterable) and isinstance(input_list[0][0], complex): + self.rho = np.array(input_list) + check_size_consistency(self.rho) + assert check_unit_trace(self.rho) + assert check_psd(self.rho) + self.Nqubit = self.rho.shape[0].bit_length() - 1 + return + statevec = Statevec(nqubit, data) + self.rho = np.outer(statevec.psi, statevec.psi.conj()) + self.Nqubit = len(statevec.dims()) def __repr__(self): return f"DensityMatrix object , with density matrix {self.rho} and shape{self.dims()}." diff --git a/tests/test_density_matrix.py b/tests/test_density_matrix.py index e9bb8bb84..a8268b661 100644 --- a/tests/test_density_matrix.py +++ b/tests/test_density_matrix.py @@ -34,7 +34,7 @@ def test_init_with_invalid_data_fail(self): # check with hermitian dm but not unit trace with self.assertRaises(ValueError): - DensityMatrix(data=randobj.rand_herm(2 ** np.random.randint(2, 5))) + DensityMatrix(randobj.rand_herm(2 ** np.random.randint(2, 5))) # check with non hermitian dm but unit trace with self.assertRaises(ValueError): @@ -80,13 +80,9 @@ def test_init_without_data_success(self): def test_init_with_data_success(self): # don't use rand_dm here since want to check for n in range(3): - data = randobj.rand_herm(2**n) - - data /= np.trace(data) - dm = DensityMatrix(data=data) + dm = randobj.rand_dm(2**n) assert dm.Nqubit == n assert dm.rho.shape == (2**n, 2**n) - assert np.allclose(dm.rho, data) def test_evolve_single_fail(self): dm = DensityMatrix(nqubit=2) @@ -425,7 +421,7 @@ def test_evolve_fail(self): def test_normalize(self): # tmp = np.random.rand(4, 4) + 1j * np.random.rand(4, 4) - data = randobj.rand_herm(2 ** np.random.randint(2, 4)) + data = randobj.rand_dm(2 ** np.random.randint(2, 4), dm_dtype=False) dm = DensityMatrix(data / data.trace()) dm.normalize() @@ -479,10 +475,8 @@ def test_ptrace(self): def test_apply_dephasing_channel(self): # check on single qubit first # # create random density matrix - # data = randobj.rand_herm(2 ** np.random.randint(2, 4)) - data = randobj.rand_herm(2) - data /= np.trace(data) - dm = DensityMatrix(data=data) + # data = randobj.rand_dm(2 ** np.random.randint(2, 4)) + dm = randobj.rand_dm(2) # copy of initial dm rho_test = deepcopy(dm.rho) @@ -563,10 +557,8 @@ def test_apply_dephasing_channel(self): def test_apply_depolarising_channel(self): # check on single qubit first # # create random density matrix - # data = randobj.rand_herm(2 ** np.random.randint(2, 4)) - data = randobj.rand_herm(2) - data /= np.trace(data) - dm = DensityMatrix(data=data) + # data = randobj.rand_dm(2 ** np.random.randint(2, 4)) + dm = randobj.rand_dm(2) # copy of initial dm rho_test = deepcopy(dm.rho) diff --git a/tests/test_statevec.py b/tests/test_statevec.py index d3ac543e6..e4a8629f0 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -4,7 +4,7 @@ from graphix.states import BasicStates, PlanarState import graphix.pauli import graphix.random_objects as randobj -from graphix.sim.statevec import Statevec, NotNormalizedError +from graphix.sim.statevec import Statevec import functools From 2c0c2ee8c2c5f3cdd40f42381feae119876cbf36 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Fri, 5 Apr 2024 15:51:49 +0200 Subject: [PATCH 13/48] Fix tests for input states --- graphix/sim/density_matrix.py | 53 ++++++++++++++++--------------- graphix/sim/statevec.py | 57 ++++++++++++++++++---------------- graphix/sim/tensornet.py | 3 ++ graphix/types.py | 10 ++++++ tests/test_density_matrix.py | 31 +++++++++--------- tests/test_statevec.py | 39 ++++++++++++----------- tests/test_statevec_backend.py | 4 +-- 7 files changed, 106 insertions(+), 91 deletions(-) create mode 100644 graphix/types.py diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index b9317456b..b343972fd 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -4,8 +4,12 @@ """ from copy import deepcopy +import typing +import functools +import numbers import numpy as np +import pydantic from graphix.linalg_validations import check_square, check_hermitian, check_unit_trace, check_psd from graphix.channels import KrausChannel @@ -14,14 +18,16 @@ from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, meas_op, Statevec import graphix.sim.base_backend import graphix.states +import graphix.types -import typing -from typing_extensions import Annotated -from annotated_types import Ge -import functools - - -PositiveInt = Annotated[int, Ge(0)] # includes 0 +Data = typing.Union[ + graphix.states.State, + "DensityMatrix", + Statevec, + typing.Iterable[graphix.states.State], + typing.Iterable[numbers.Number], + typing.Iterable[typing.Iterable[numbers.Number]], +] class DensityMatrix: @@ -29,15 +35,8 @@ class DensityMatrix: def __init__( self, - nqubit: typing.Optional[PositiveInt] = None, - data: typing.Union[ - graphix.states.State, - "DensityMatrix", - Statevec, - typing.Iterable[graphix.states.State], - typing.Iterable[complex], - typing.Iterable[typing.Iterable[complex]], - ] = graphix.states.BasicStates.PLUS, + data: Data = graphix.states.BasicStates.PLUS, + nqubit: typing.Optional[graphix.types.PositiveInt] = None, ): """ @@ -47,8 +46,9 @@ def __init__( data : DensityMatrix, list, tuple, np.ndarray or None Density matrix of shape (2**nqubits, 2**nqubits). nqubit : int - Number of qubits. Default is 1. If both `data` and `nqubit` are specified, `nqubit` is ignored. + Number of qubits. Default is 1. If both `data` and `nqubit` are specified, consistency is checked. """ + pydantic.TypeAdapter(typing.Optional[graphix.types.PositiveInt]).validate_python(nqubit) def check_size_consistency(mat): if nqubit is not None and mat.shape != (2**nqubit, 2**nqubit): @@ -64,14 +64,15 @@ def check_size_consistency(mat): if isinstance(data, typing.Iterable): input_list = list(data) if len(input_list) != 0: - if isinstance(input_list[0], typing.Iterable) and isinstance(input_list[0][0], complex): + if isinstance(input_list[0], typing.Iterable) and isinstance(input_list[0][0], numbers.Number): self.rho = np.array(input_list) + assert check_square(self.rho) check_size_consistency(self.rho) assert check_unit_trace(self.rho) assert check_psd(self.rho) self.Nqubit = self.rho.shape[0].bit_length() - 1 return - statevec = Statevec(nqubit, data) + statevec = Statevec(data, nqubit) self.rho = np.outer(statevec.psi, statevec.psi.conj()) self.Nqubit = len(statevec.dims()) @@ -306,7 +307,7 @@ def apply_channel(self, channel: KrausChannel, qargs): class DensityMatrixBackend(graphix.sim.base_backend.Backend): """MBQC simulator with density matrix method.""" - def __init__(self, pattern, max_qubit_num=12, pr_calc=True): + def __init__(self, pattern, max_qubit_num=12, pr_calc=True, input_state: Data = graphix.states.BasicStates.PLUS): """ Parameters ---------- @@ -331,7 +332,10 @@ def __init__(self, pattern, max_qubit_num=12, pr_calc=True): raise ValueError("Pattern.max_space is larger than max_qubit_num. Increase max_qubit_num and try again.") super().__init__(pr_calc) - def add_nodes(self, nodes, qubit_to_add=None): + # initialize input qubits to desired init_state + self.add_nodes(pattern.input_nodes, input_state) + + def add_nodes(self, nodes, input_state: Data = graphix.states.BasicStates.PLUS): """add new qubit to the internal density matrix and asign the corresponding node number to list self.node_index. @@ -345,12 +349,7 @@ def add_nodes(self, nodes, qubit_to_add=None): if not self.state: self.state = DensityMatrix(nqubit=0) n = len(nodes) - if qubit_to_add is None: - dm_to_add = DensityMatrix(nqubit=n) - else: - assert isinstance(qubit_to_add, DensityMatrix) - assert qubit_to_add.nqubit == 1 - dm_to_add = qubit_to_add + dm_to_add = DensityMatrix(nqubit=n, data=input_state) self.state.tensor(dm_to_add) self.node_index.extend(nodes) self.Nqubit += n diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index ed3844ad7..a1cd0a3d5 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -1,17 +1,18 @@ from copy import deepcopy +import numbers import typing -from typing_extensions import Annotated -from annotated_types import Ge import numpy as np +import functools +import pydantic +import warnings import graphix.sim.base_backend from graphix.clifford import CLIFFORD, CLIFFORD_CONJ, CLIFFORD_MUL from graphix.ops import Ops import graphix.states import graphix.pauli -import functools -import warnings +import graphix.types # Python >= 3.9 # from collections.abc import Iterable # or use Protocols? @@ -27,7 +28,7 @@ def __init__( self, pattern, input_state: typing.Union[ - graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[complex] + graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[numbers.Number] ] = graphix.states.BasicStates.PLUS, max_qubit_num=20, pr_calc=True, @@ -84,7 +85,7 @@ def add_nodes(self, nodes, input_state=graphix.states.BasicStates.PLUS): if not self.state: self.state = Statevec(nqubit=0) n = len(nodes) - sv_to_add = Statevec(nqubit=n, state=input_state) + sv_to_add = Statevec(nqubit=n, data=input_state) self.state.tensor(sv_to_add) self.node_index.extend(nodes) self.Nqubit += n @@ -201,7 +202,6 @@ def meas_op(angle, vop=0, plane="XY", choice=0): [[[[1, 0], [0, 0]], [[0, 0], [1, 0]]], [[[0, 1], [0, 0]], [[0, 0], [0, 1]]]], dtype=np.complex128, ) -PositiveInt = Annotated[int, Ge(0)] # includes 0 class Statevec: @@ -210,17 +210,17 @@ class Statevec: # TODO at this stage no need for indices just be careful of the ordering in add_nodes def __init__( self, - nqubit: typing.Optional[PositiveInt] = None, - state: typing.Union[ - graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[complex] + data: typing.Union[ + graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[numbers.Number] ] = graphix.states.BasicStates.PLUS, + nqubit: typing.Optional[graphix.types.PositiveInt] = None, ): """Initialize statevector Parameters ---------- - state : is either + data : is either - a single state (:class:`graphix.states.State` object). THen prepares all nodes in that state (tensor product) - a dictionary mapping the inputs to a :class:`graphix.states.State` object - an arbitrary :class:`graphix.statevec.Statevec` object (arbitrary input) # TODO work on that since just copy? @@ -232,27 +232,25 @@ def __init__( Defaults to |+> states and 1 qubit. If nqubit > 1 and only one state : tensor all of them. Use the tensor method instead of hard code. """ - # convert single qubit State object to Statevector - # NOTE make plane, angle attributes of Statevec? - # this will be called by the nqb > 1 case - # instantiate a one - if isinstance(state, Statevec): + pydantic.TypeAdapter(typing.Optional[graphix.types.PositiveInt]).validate_python(nqubit) + + if isinstance(data, Statevec): # assert nqubit is None or len(state.flatten()) == 2**nqubit - if nqubit is not None and len(state.flatten()) != 2**nqubit: + if nqubit is not None and len(data.flatten()) != 2**nqubit: raise ValueError( - f"Inconsistent parameters between nqubit = {nqubit} and the inferred number of qubit = {len(state.flatten())}." + f"Inconsistent parameters between nqubit = {nqubit} and the inferred number of qubit = {len(data.flatten())}." ) - self.psi = state.psi.copy() + self.psi = data.psi.copy() return - if isinstance(state, graphix.states.State): + if isinstance(data, graphix.states.State): if nqubit is None: - raise ValueError("Incorrect value for nqubit.") - input_list = [state] * nqubit - elif isinstance(state, typing.Iterable): - input_list = list(state) + nqubit = 1 + input_list = [data] * nqubit + elif isinstance(data, typing.Iterable): + input_list = list(data) else: - raise TypeError("Incorrect type for input state") + raise TypeError(f"Incorrect type for data: {type(data)}") if len(input_list) == 0: if nqubit is not None and nqubit != 0: @@ -263,6 +261,7 @@ def __init__( # self.Nqubit = 0 else: if isinstance(input_list[0], graphix.states.State): + graphix.types.check_list_elements(input_list, graphix.states.State) if nqubit is None: nqubit = len(input_list) elif nqubit != len(input_list): @@ -271,7 +270,8 @@ def __init__( tmp_psi = functools.reduce(np.kron, list_of_sv) # reshape self.psi = tmp_psi.reshape((2,) * nqubit) - else: + elif isinstance(input_list[0], numbers.Number): + graphix.types.check_list_elements(input_list, numbers.Number) if nqubit is None: length = len(input_list) if length & (length - 1): @@ -285,6 +285,10 @@ def __init__( # just reshape # NOTE too many conversions to numpy arrays? self.psi = psi.reshape((2,) * nqubit) + else: + raise TypeError( + f"First element of data has type {type(input_list[0])} whereas Number or State is expected" + ) # self.Nqubit = state.Nqubit def __repr__(self): @@ -425,7 +429,6 @@ def tensor(self, other): # NOTE on tensor form not vector # deprecated total_num = len(self.dims()) + len(other.dims()) - print(total_num) # self.Nqubit += other.Nqubit self.psi = np.kron(psi_self, psi_other).reshape((2,) * total_num) # self.Nqubit = len(self.dims()) diff --git a/graphix/sim/tensornet.py b/graphix/sim/tensornet.py index ca3cd6aa1..d0475962c 100644 --- a/graphix/sim/tensornet.py +++ b/graphix/sim/tensornet.py @@ -70,6 +70,9 @@ def __init__(self, pattern, graph_prep="auto", **kwargs): self._decomposed_cz = _get_decomposed_cz() self._isolated_nodes = pattern.get_isolated_nodes() + # initialize input qubits to desired init_state + self.add_nodes(pattern.input_nodes) + def add_nodes(self, nodes): """Add nodes to the network diff --git a/graphix/types.py b/graphix/types.py new file mode 100644 index 000000000..165c40e7e --- /dev/null +++ b/graphix/types.py @@ -0,0 +1,10 @@ +import annotated_types +import typing_extensions + +PositiveInt = typing_extensions.Annotated[int, annotated_types.Ge(0)] # includes 0 + + +def check_list_elements(l, ty): + for index, item in enumerate(l): + if not isinstance(item, ty): + raise TypeError(f"data[{index}] has type {type(item)} whereas {ty} is expected") diff --git a/tests/test_density_matrix.py b/tests/test_density_matrix.py index a8268b661..69364c32e 100644 --- a/tests/test_density_matrix.py +++ b/tests/test_density_matrix.py @@ -3,6 +3,7 @@ from copy import deepcopy import numpy as np +import pydantic import graphix.random_objects as randobj from graphix import Circuit @@ -10,17 +11,17 @@ from graphix.ops import Ops from graphix.sim.density_matrix import DensityMatrix, DensityMatrixBackend from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, Statevec, StatevectorBackend - +import graphix.states class TestDensityMatrix(unittest.TestCase): """Test for DensityMatrix class.""" def test_init_without_data_fail(self): - with self.assertRaises(AssertionError): + with self.assertRaises(pydantic.ValidationError): DensityMatrix(nqubit=-2) - with self.assertRaises(TypeError): + with self.assertRaises(pydantic.ValidationError): DensityMatrix(nqubit="hello") - with self.assertRaises(TypeError): + with self.assertRaises(pydantic.ValidationError): DensityMatrix(nqubit=[]) def test_init_with_invalid_data_fail(self): @@ -29,7 +30,7 @@ def test_init_with_invalid_data_fail(self): with self.assertRaises(TypeError): DensityMatrix(1) # deprecated data shape (these test might be unnecessary) - with self.assertRaises(ValueError): + with self.assertRaises(TypeError): DensityMatrix([1, 2, [3]]) # check with hermitian dm but not unit trace @@ -52,7 +53,7 @@ def test_init_with_invalid_data_fail(self): DensityMatrix(data=np.random.rand(3, 2)) # check higher dimensional matrix - with self.assertRaises(ValueError): + with self.assertRaises(TypeError): DensityMatrix(data=np.random.rand(2, 2, 3)) # check square and hermitian but with incorrect dimension (non-qubit type) @@ -70,7 +71,7 @@ def test_init_without_data_success(self): assert dm.rho.shape == (2**n, 2**n) assert np.allclose(dm.rho, expected_density_matrix) - dm = DensityMatrix(plus_state=False, nqubit=n) + dm = DensityMatrix(data=graphix.states.BasicStates.ZERO, nqubit=n) expected_density_matrix = np.zeros((2**n, 2**n)) expected_density_matrix[0, 0] = 1 assert dm.Nqubit == n @@ -778,16 +779,15 @@ def test_init_success(self): backend = DensityMatrixBackend(pattern) assert backend.pattern == pattern assert backend.results == pattern.results - assert backend.state is None - assert backend.node_index == [] - assert backend.Nqubit == 0 + assert backend.node_index == [0] + assert backend.Nqubit == 1 assert backend.max_qubit_num == 12 def test_add_nodes(self): circ = Circuit(1) pattern = circ.transpile().pattern backend = DensityMatrixBackend(pattern) - backend.add_nodes([0, 1]) + backend.add_nodes([1]) expected_matrix = np.array([0.25] * 16).reshape(4, 4) np.testing.assert_allclose(backend.state.rho, expected_matrix) @@ -795,7 +795,7 @@ def test_entangle_nodes(self): circ = Circuit(1) pattern = circ.transpile().pattern backend = DensityMatrixBackend(pattern) - backend.add_nodes([0, 1]) + backend.add_nodes([1]) backend.entangle_nodes((0, 1)) expected_matrix = np.array([[1, 1, 1, -1], [1, 1, 1, -1], [1, 1, 1, -1], [-1, -1, -1, 1]]) / 4 np.testing.assert_allclose(backend.state.rho, expected_matrix) @@ -809,7 +809,7 @@ def test_measure(self): pattern = circ.transpile().pattern backend = DensityMatrixBackend(pattern) - backend.add_nodes([0, 1, 2]) + backend.add_nodes([1, 2]) backend.entangle_nodes((0, 1)) backend.entangle_nodes((1, 2)) backend.measure(backend.pattern[-4]) @@ -825,7 +825,7 @@ def test_measure_pr_calc(self): pattern = circ.transpile().pattern backend = DensityMatrixBackend(pattern, pr_calc=True) - backend.add_nodes([0, 1, 2]) + backend.add_nodes([1, 2]) backend.entangle_nodes((0, 1)) backend.entangle_nodes((1, 2)) backend.measure(backend.pattern[-4]) @@ -844,7 +844,8 @@ def test_correct_byproduct(self): pattern = circ.transpile().pattern backend = DensityMatrixBackend(pattern) - backend.add_nodes([0, 1, 2]) + # node 0 initialized in Backend + backend.add_nodes([1, 2]) backend.entangle_nodes((0, 1)) backend.entangle_nodes((1, 2)) backend.measure(backend.pattern[-4]) diff --git a/tests/test_statevec.py b/tests/test_statevec.py index e4a8629f0..c0793c561 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -28,31 +28,31 @@ def test_default_success(self): def test_basicstates_success(self): # minus - vec = Statevec(nqubit=1, state=BasicStates.MINUS) + vec = Statevec(nqubit=1, data=BasicStates.MINUS) np.testing.assert_allclose(vec.psi, np.array([1, -1] / np.sqrt(2))) # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # zero - vec = Statevec(nqubit=1, state=BasicStates.ZERO) + vec = Statevec(nqubit=1, data=BasicStates.ZERO) np.testing.assert_allclose(vec.psi, np.array([1, 0]), rtol=0, atol=1e-15) # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # one - vec = Statevec(nqubit=1, state=BasicStates.ONE) + vec = Statevec(nqubit=1, data=BasicStates.ONE) np.testing.assert_allclose(vec.psi, np.array([0, 1]), rtol=0, atol=1e-15) # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # plus_i - vec = Statevec(nqubit=1, state=BasicStates.PLUS_I) + vec = Statevec(nqubit=1, data=BasicStates.PLUS_I) np.testing.assert_allclose(vec.psi, np.array([1, 1j] / np.sqrt(2))) # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # minus_i - vec = Statevec(nqubit=1, state=BasicStates.MINUS_I) + vec = Statevec(nqubit=1, data=BasicStates.MINUS_I) np.testing.assert_allclose(vec.psi, np.array([1, -1j] / np.sqrt(2))) #assert vec.Nqubit == 1 assert len(vec.dims()) == 1 @@ -66,7 +66,7 @@ def test_default_tensor_success(self): # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb - vec = Statevec(nqubit=nqb, state=BasicStates.MINUS_I) + vec = Statevec(nqubit=nqb, data=BasicStates.MINUS_I) sv_list = [BasicStates.MINUS_I.get_statevector() for _ in range(nqb)] sv = functools.reduce(np.kron, sv_list) np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) @@ -77,7 +77,7 @@ def test_default_tensor_success(self): rand_angle = self.rng.random() * 2 * np.pi rand_plane = self.rng.choice(np.array([i for i in graphix.pauli.Plane])) state = PlanarState(plane=rand_plane, angle=rand_angle) - vec = Statevec(nqubit=nqb, state=state) + vec = Statevec(nqubit=nqb, data=state) sv_list = [state.get_statevector() for _ in range(nqb)] sv = functools.reduce(np.kron, sv_list) np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) @@ -88,7 +88,7 @@ def test_default_tensor_success(self): rand_angles = self.rng.random(nqb) * 2 * np.pi rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) states = [PlanarState(plane = i, angle = j) for i, j in zip(rand_planes, rand_angles)] - vec = Statevec(nqubit=nqb, state=states) + vec = Statevec(nqubit=nqb, data=states) sv_list = [state.get_statevector() for state in states] sv = functools.reduce(np.kron, sv_list) np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) @@ -100,7 +100,7 @@ def test_data_success(self): l = 2 ** nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) - vec = Statevec(state = rand_vec) + vec = Statevec(data=rand_vec) np.testing.assert_allclose(vec.psi, rand_vec.reshape((2,) * nqb)) # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb @@ -112,7 +112,7 @@ def test_data_dim_fail(self): rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) with self.assertRaises(ValueError): - vec = Statevec(state = rand_vec) + vec = Statevec(data=rand_vec) # with less qubit than number of qubits inferred from a correct state vect # returns a truncated statevec that is hence not normalized @@ -122,7 +122,7 @@ def test_data_dim_fail_mismatch(self): rand_vec = self.rng.random(2 ** nqb) + 1j * self.rng.random(2 ** nqb) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) with self.assertRaises(ValueError): - vec = Statevec(nqubit = 2, state = rand_vec) + vec = Statevec(nqubit = 2, data=rand_vec) # fail: not normalized def test_data_norm_fail(self): @@ -130,12 +130,11 @@ def test_data_norm_fail(self): l = 2 ** nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) with self.assertRaises(ValueError): - vec = Statevec(state = rand_vec) + vec = Statevec(data=rand_vec) - # fail: no nqubit provided. State defaults to PLUS - def test_default_fail(self): - with self.assertRaises(ValueError): - vec = Statevec() + def test_defaults_to_one(self): + vec = Statevec() + assert len(vec.dims()) == 1 # try copying Statevec input def test_copy_success(self): @@ -143,9 +142,9 @@ def test_copy_success(self): l = 2 ** nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) - test_vec = Statevec(state = rand_vec) + test_vec = Statevec(data=rand_vec) # try to copy it - vec = Statevec(state = test_vec) + vec = Statevec(data=test_vec) np.testing.assert_allclose(vec.psi, test_vec.psi) # assert vec.Nqubit == test_vec.Nqubit @@ -157,7 +156,7 @@ def test_copy_fail(self): l = 2 ** nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) - test_vec = Statevec(state = rand_vec) + test_vec = Statevec(data=rand_vec) with self.assertRaises(ValueError): - vec = Statevec(nqubit = l - 1, state = test_vec) + vec = Statevec(nqubit=l - 1, data=test_vec) diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index 3f45d1bf6..e7ab78598 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -72,7 +72,7 @@ def test_init_success(self): # minus state backend = StatevectorBackend(self.hadamardpattern, input_state = BasicStates.MINUS) - vec = Statevec(nqubit=1, state = BasicStates.MINUS) + vec = Statevec(nqubit=1, data=BasicStates.MINUS) np.testing.assert_allclose(vec.psi, backend.state.psi) # assert backend.state.Nqubit == 1 assert len(backend.state.dims()) == 1 @@ -82,7 +82,7 @@ def test_init_success(self): rand_plane = self.rng.choice(np.array([i for i in graphix.pauli.Plane])) state = PlanarState(plane = rand_plane, angle = rand_angle) backend = StatevectorBackend(self.hadamardpattern, input_state = state) - vec = Statevec(nqubit=1, state = state) + vec = Statevec(nqubit=1, data=state) np.testing.assert_allclose(vec.psi, backend.state.psi) # assert backend.state.Nqubit == 1 assert len(backend.state.dims()) == 1 From 506f9f6b24b5b02f3138d2405c45e0eb936e9726 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Mon, 29 Apr 2024 17:03:29 +0200 Subject: [PATCH 14/48] implem arbitrary inputs for circuit statevector --- graphix/linalg_validations.py | 2 +- graphix/ops.py | 2 + graphix/pattern.py | 6 +- graphix/random_objects.py | 29 +-- graphix/sim/density_matrix.py | 30 ++- graphix/sim/statevec.py | 10 +- graphix/states.py | 6 +- graphix/transpiler.py | 11 +- tests/test_density_matrix.py | 329 ++++++++++++++++++++++++--------- tests/test_pattern.py | 43 ++++- tests/test_statevec.py | 54 +++--- tests/test_statevec_backend.py | 46 +++-- 12 files changed, 392 insertions(+), 176 deletions(-) diff --git a/graphix/linalg_validations.py b/graphix/linalg_validations.py index f452504c8..e4fa101c1 100644 --- a/graphix/linalg_validations.py +++ b/graphix/linalg_validations.py @@ -35,7 +35,7 @@ def check_psd(matrix: np.ndarray, tol: float = 1e-15) -> bool: evals = np.linalg.eigvalsh(matrix) if not all(evals >= -tol): - raise ValueError("The matrix is not positive semi-definite.") + raise ValueError("The matrix {matrix} is not positive semi-definite.") return True diff --git a/graphix/ops.py b/graphix/ops.py index 6c90d101c..d8aa4b9ae 100644 --- a/graphix/ops.py +++ b/graphix/ops.py @@ -23,6 +23,8 @@ class Ops: """Basic single- and two-qubits operators""" + # class attributes. DOn' need to instantiate. Inherited by all class members + x = np.array([[0, 1], [1, 0]]) y = np.array([[0, -1j], [1j, 0]]) z = np.array([[1, 0], [0, -1]]) diff --git a/graphix/pattern.py b/graphix/pattern.py index 032ec6475..fcbd80c5b 100644 --- a/graphix/pattern.py +++ b/graphix/pattern.py @@ -818,7 +818,7 @@ def connected_edges(self, node, edges): # like in def get_measurement_order_from_flow(self): with self.get_graph() # FIXME # BUG - + connected = set() for edge in edges: if edge[0] == node: @@ -1016,9 +1016,9 @@ def get_max_degree(self): degree = g.degree() max_degree = max([i for i in dict(degree).values()]) return max_degree - + # TODO functools.cache() It is called in get measurement order from (g)flow - # + # # It is called in get measurement order from (g)flow def get_graph(self): """returns the list of nodes and edges from the command sequence, diff --git a/graphix/random_objects.py b/graphix/random_objects.py index b583c2849..2ed2e61bc 100644 --- a/graphix/random_objects.py +++ b/graphix/random_objects.py @@ -1,5 +1,5 @@ import numpy as np -import numpy.typing as npt +from typing import Optional import scipy.linalg from scipy.stats import unitary_group @@ -29,27 +29,28 @@ def rand_unit(l: int): UNITS = np.array([1, 1j]) -def rand_dm(dim: int, rank: int = None, dm_dtype=True) -> DensityMatrix: - """Returns a "density matrix" as a DensityMatrix object ie a positive-semidefinite (hence Hermitian) matrix with unit trace - Note, not a proper DM since its dim can be something else than a power of 2. - The rank is random between 1 (pure) and dim if not specified - Thanks to Ulysse Chabaud. +def rand_dm(dim: int, rank: Optional[int] = None, dm_dtype=True) -> DensityMatrix | np.ndarray: + """Utility to generate random density matrices (positive semi-definite matrices with unit trace). + Returns either a :class:`graphix.sim.density_matrix.DensityMatrix` or a :class:`np.ndarray` depending on the parameter `dm_dtype`. :param dim: Linear dimension of the (square) matrix :type dim: int - :param rank: If rank not specified then random between 1 and matrix dimension - If rank is one : then pure state else mixed state. Defaults to None + :param rank: Rank of the density matrix (1 = pure state). If not specified then sent to dim (maximal rank). + Defaults to None :type rank: int, optional - :param dm_dtype: If True returns a :class:`graphix.sim.density_matrix.DensityMatrix` or a numpy.ndarray if False. Defaults to True. + :param dm_dtype: If `True` returns a :class:`graphix.sim.density_matrix.DensityMatrix` object. If `False`returns a :class:`np.ndarray` :type dm_dtype: bool, optional - :return: Random density matrix as a :class:`graphix.sim.density_matrix.DensityMatrix` object or a numpy.ndarray. - :rtype: :class:`graphix.sim.density_matrix.DensityMatrix` or numpy.ndarray. - + :return: the density matrix in the specified format. + :rtype: DensityMatrix | np.ndarray + .. note:: + Thanks to Ulysse Chabaud. + .. warning:: + Note that setting `dm_dtype=False` allows to generate "density matrices" inconsistent with qubits i.e. with dimensions not being powers of 2. """ # if not provided, use a random value. if rank is None: - rank = np.random.randint(1, dim + 1) + rank = dim # np.random.randint(1, dim + 1) evals = np.random.rand(rank) @@ -68,7 +69,7 @@ def rand_dm(dim: int, rank: int = None, dm_dtype=True) -> DensityMatrix: return dm -def rand_gauss_cpx_mat(dim: int, sig: float = 1 / np.sqrt(2)) -> npt.NDArray: +def rand_gauss_cpx_mat(dim: int, sig: float = 1 / np.sqrt(2)) -> np.ndarray: """ Returns a square array of standard normal complex random variates. Code from QuTiP: https://qutip.org/docs/4.0.2/modules/qutip/random_objects.html diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index b343972fd..cb5384106 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -10,6 +10,7 @@ import numpy as np import pydantic +import typing from graphix.linalg_validations import check_square, check_hermitian, check_unit_trace, check_psd from graphix.channels import KrausChannel @@ -35,7 +36,7 @@ class DensityMatrix: def __init__( self, - data: Data = graphix.states.BasicStates.PLUS, + data: typing.Optional[Data] = graphix.states.BasicStates.PLUS, nqubit: typing.Optional[graphix.types.PositiveInt] = None, ): @@ -58,26 +59,33 @@ def check_size_consistency(mat): if isinstance(data, DensityMatrix): check_size_consistency(data) + # safe: https://numpy.org/doc/stable/reference/generated/numpy.ndarray.copy.html self.rho = data.rho.copy() self.Nqubit = data.Nqubit return if isinstance(data, typing.Iterable): input_list = list(data) if len(input_list) != 0: - if isinstance(input_list[0], typing.Iterable) and isinstance(input_list[0][0], numbers.Number): - self.rho = np.array(input_list) - assert check_square(self.rho) - check_size_consistency(self.rho) - assert check_unit_trace(self.rho) - assert check_psd(self.rho) - self.Nqubit = self.rho.shape[0].bit_length() - 1 - return + # do try except else? + # needed since Object are iterable but not subscribable! + try: + if isinstance(input_list[0], typing.Iterable) and isinstance(input_list[0][0], numbers.Number): + self.rho = np.array(input_list) + assert check_square(self.rho) + check_size_consistency(self.rho) + assert check_unit_trace(self.rho) + assert check_psd(self.rho) + self.Nqubit = self.rho.shape[0].bit_length() - 1 + return + except TypeError: + pass statevec = Statevec(data, nqubit) + # NOTE sthis works since np.outer flattens the inputs! self.rho = np.outer(statevec.psi, statevec.psi.conj()) self.Nqubit = len(statevec.dims()) def __repr__(self): - return f"DensityMatrix object , with density matrix {self.rho} and shape{self.dims()}." + return f"DensityMatrix object, with density matrix {self.rho} and shape {self.dims()}." def evolve_single(self, op, i): """Single-qubit operation. @@ -172,6 +180,8 @@ def expectation_single(self, op, i): return np.trace(st1.rho) + # TODO + # @property def dims(self): return self.rho.shape diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index a1cd0a3d5..c31090047 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -203,6 +203,10 @@ def meas_op(angle, vop=0, plane="XY", choice=0): dtype=np.complex128, ) +SV_Data = typing.Union[ + graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[numbers.Number] +] + class Statevec: """Statevector object""" @@ -210,9 +214,7 @@ class Statevec: # TODO at this stage no need for indices just be careful of the ordering in add_nodes def __init__( self, - data: typing.Union[ - graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[numbers.Number] - ] = graphix.states.BasicStates.PLUS, + data: typing.Optional[SV_Data] = graphix.states.BasicStates.PLUS, nqubit: typing.Optional[graphix.types.PositiveInt] = None, ): @@ -265,7 +267,7 @@ def __init__( if nqubit is None: nqubit = len(input_list) elif nqubit != len(input_list): - raise ValueError("Mismatch between nqubit and length of input state") + raise ValueError("Mismatch between nqubit and length of input state.") list_of_sv = [s.get_statevector() for s in input_list] tmp_psi = functools.reduce(np.kron, list_of_sv) # reshape diff --git a/graphix/states.py b/graphix/states.py index 9dfdadfa7..81ecd4877 100644 --- a/graphix/states.py +++ b/graphix/states.py @@ -9,7 +9,7 @@ # generic class State for all States class State(abc.ABC): - """Abstract base class for states objects. + """Abstract base class for single qubit states objects. Only requirement for concrete classes is to have a get_statevector() method that returns the statevector representation of the state @@ -19,6 +19,10 @@ class State(abc.ABC): def get_statevector(self) -> np.ndarray: pass + def get_densitymatrix(self) -> np.ndarray: + # return DM in 2**n x 2**n dim (2x2 here) + return np.outer(self.get_statevector(), self.get_statevector().conj()) + # don't turn it into Statevec here # Weird not to allow all states? diff --git a/graphix/transpiler.py b/graphix/transpiler.py index 0e1de7583..ff87cf52e 100644 --- a/graphix/transpiler.py +++ b/graphix/transpiler.py @@ -16,7 +16,10 @@ import graphix.sim.base_backend from graphix.ops import Ops from graphix.pattern import Pattern -from graphix.sim.statevec import Statevec +import graphix.sim.base_backend +from graphix.sim.statevec import Statevec, SV_Data +from graphix.sim.density_matrix import DensityMatrix, Data +import graphix.pauli @dataclasses.dataclass @@ -1359,8 +1362,8 @@ def _sort_outputs(self, pattern: Pattern, output_nodes: Sequence[int]): elif cmd[1] in old_out: cmd[1] = output_nodes[old_out.index(cmd[1])] - def simulate_statevector(self, input_state: Optional[Statevec] = None) -> SimulateResult: - """Run statevector simultion of the gate sequence, using graphix.Statevec + def simulate_statevector(self, input_state: Optional[SV_Data] = None): + """Run statevector simulation of the gate sequence, using graphix.Statevec Parameters ---------- @@ -1375,7 +1378,7 @@ def simulate_statevector(self, input_state: Optional[Statevec] = None) -> Simula if input_state is None: state = Statevec(nqubit=self.width) else: - state = input_state + state = Statevec(nqubit=self.width, data=input_state) classical_measures = [] diff --git a/tests/test_density_matrix.py b/tests/test_density_matrix.py index 69364c32e..7a53410f5 100644 --- a/tests/test_density_matrix.py +++ b/tests/test_density_matrix.py @@ -1,9 +1,10 @@ import random import unittest +import pydantic from copy import deepcopy import numpy as np -import pydantic +import pytest import graphix.random_objects as randobj from graphix import Circuit @@ -12,87 +13,186 @@ from graphix.sim.density_matrix import DensityMatrix, DensityMatrixBackend from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, Statevec, StatevectorBackend import graphix.states +import graphix.pauli +import tests.random_circuit +import functools + class TestDensityMatrix(unittest.TestCase): """Test for DensityMatrix class.""" + def setUp(self): + # set up the random numbers + self.rng = np.random.default_rng() + def test_init_without_data_fail(self): - with self.assertRaises(pydantic.ValidationError): + with pytest.raises(pydantic.ValidationError): DensityMatrix(nqubit=-2) - with self.assertRaises(pydantic.ValidationError): + with pytest.raises(pydantic.ValidationError): DensityMatrix(nqubit="hello") - with self.assertRaises(pydantic.ValidationError): + with pytest.raises(pydantic.ValidationError): DensityMatrix(nqubit=[]) def test_init_with_invalid_data_fail(self): - with self.assertRaises(TypeError): + with pytest.raises(TypeError): DensityMatrix("hello") - with self.assertRaises(TypeError): + with pytest.raises(TypeError): DensityMatrix(1) # deprecated data shape (these test might be unnecessary) - with self.assertRaises(TypeError): + with pytest.raises(TypeError): DensityMatrix([1, 2, [3]]) # check with hermitian dm but not unit trace - with self.assertRaises(ValueError): + with pytest.raises(ValueError): DensityMatrix(randobj.rand_herm(2 ** np.random.randint(2, 5))) # check with non hermitian dm but unit trace - with self.assertRaises(ValueError): + with pytest.raises(ValueError): l = 2 ** np.random.randint(2, 5) tmp = np.random.rand(l, l) + 1j * np.random.rand(l, l) DensityMatrix(data=tmp / np.trace(tmp)) # check with non hermitian dm and not unit trace - with self.assertRaises(ValueError): + with pytest.raises(ValueError): l = 2 ** np.random.randint(2, 5) # np.random.randint(2, 20) DensityMatrix(data=np.random.rand(l, l) + 1j * np.random.rand(l, l)) # check not square matrix - with self.assertRaises(ValueError): + with pytest.raises(ValueError): # l = 2 ** np.random.randint(2, 5) # np.random.randint(2, 20) DensityMatrix(data=np.random.rand(3, 2)) # check higher dimensional matrix - with self.assertRaises(TypeError): + with pytest.raises(TypeError): DensityMatrix(data=np.random.rand(2, 2, 3)) # check square and hermitian but with incorrect dimension (non-qubit type) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): # not really a dm since not PSD but ok. data = randobj.rand_herm(5) data /= np.trace(data) DensityMatrix(data=data) def test_init_without_data_success(self): - for n in range(3): - dm = DensityMatrix(nqubit=n) - expected_density_matrix = np.outer(np.ones((2,) * n), np.ones((2,) * n)) / 2**n - assert dm.Nqubit == n - assert dm.rho.shape == (2**n, 2**n) - assert np.allclose(dm.rho, expected_density_matrix) + n = np.random.randint(2, 5) - dm = DensityMatrix(data=graphix.states.BasicStates.ZERO, nqubit=n) - expected_density_matrix = np.zeros((2**n, 2**n)) - expected_density_matrix[0, 0] = 1 - assert dm.Nqubit == n - assert dm.rho.shape == (2**n, 2**n) - assert np.allclose(dm.rho, expected_density_matrix) + dm = DensityMatrix(nqubit=n) + expected_density_matrix = np.outer(np.ones((2,) * n), np.ones((2,) * n)) / 2**n + assert dm.Nqubit == n + assert dm.rho.shape == (2**n, 2**n) + assert np.allclose(dm.rho, expected_density_matrix) + + dm = DensityMatrix(data=graphix.states.BasicStates.ZERO, nqubit=n) + expected_density_matrix = np.zeros((2**n, 2**n)) + expected_density_matrix[0, 0] = 1 + assert dm.Nqubit == n + assert dm.rho.shape == (2**n, 2**n) + assert np.allclose(dm.rho, expected_density_matrix) def test_init_with_data_success(self): - # don't use rand_dm here since want to check - for n in range(3): - dm = randobj.rand_dm(2**n) - assert dm.Nqubit == n - assert dm.rho.shape == (2**n, 2**n) + # explicitely use dm_dtype=False to check the constructor. + # since rand_dm relies on the constructor data validation. + n = np.random.randint(2, 5) + data = randobj.rand_dm(2**n) + dm = DensityMatrix(data=data) + assert dm.Nqubit == n + assert dm.rho.shape == (2**n, 2**n) + + def test_init_with_state_sucess(self): + # both "numerical" statevec and Statevec object + # relies on Statevec constructor validation + + nqb = self.rng.integers(2, 5) + print(f"nqb is {nqb}") + rand_angles = self.rng.random(nqb) * 2 * np.pi + rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) + states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] + vec = Statevec(data=states) + # flattens input! + expected_dm = np.outer(vec.psi, vec.psi.conj()) + + # input with a State object + dm = DensityMatrix(data=states) + assert dm.dims() == (2**nqb, 2**nqb) + assert np.allclose(dm.rho, expected_dm) + + def test_init_with_state_fail(self): + + nqb = 2 + rand_angles = self.rng.random(nqb) * 2 * np.pi + rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) + states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] + + with pytest.raises(ValueError): + dm = DensityMatrix(nqubit=1, data=states) + + with pytest.raises(ValueError): + dm = DensityMatrix(nqubit=3, data=states) + + def test_init_with_statevec_sucess(self): + # both "numerical" statevec and Statevec object + # relies on Statevec constructor validation + + nqb = self.rng.integers(2, 5) + rand_angles = self.rng.random(nqb) * 2 * np.pi + rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) + states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] + vec = Statevec(data=states) + # flattens input! + expected_dm = np.outer(vec.psi, vec.psi.conj()) + + # input with a Statevec object + dm = DensityMatrix(data=vec) + assert dm.dims() == (2**nqb, 2**nqb) + assert np.allclose(dm.rho, expected_dm) + + sv_list = [state.get_statevector() for state in states] + sv = functools.reduce(np.kron, sv_list) + + # input with a statevector DATA (not Statevec object) + dm2 = DensityMatrix(data=sv) + + print("dims", dm.dims()) + assert dm2.dims() == (2**nqb, 2**nqb) + assert np.allclose(dm2.rho, dm.rho) + assert np.allclose(dm2.rho, expected_dm) + + def test_init_with_densitymatrix_sucess(self): + # both "numerical" densitymatrix and DensityMatrix object + + nqb = self.rng.integers(2, 5) + rand_angles = self.rng.random(nqb) * 2 * np.pi + rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) + print("planes", rand_planes) + states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] + vec = Statevec(data=states) + expected_dm = np.outer(vec.psi, vec.psi.conj()) + + # input with a huge density matrix + dm_list = [state.get_densitymatrix() for state in states] + num_dm = functools.reduce(np.kron, dm_list) + + dm = DensityMatrix(data=num_dm) + + assert dm.dims() == (2**nqb, 2**nqb) + assert np.allclose(dm.rho, expected_dm) + + # check copying + dm2 = DensityMatrix(dm) + assert dm2.dims() == (2**nqb, 2**nqb) + assert np.allclose(dm2.rho, expected_dm) + assert np.allclose(dm2.rho, dm.rho) + + # assert dm.Nqubit == n + # assert dm.rho.shape == (2**n, 2**n) def test_evolve_single_fail(self): dm = DensityMatrix(nqubit=2) # generate random 4 x 4 unitary matrix op = randobj.rand_unit(4) - with self.assertRaises(AssertionError): + with pytest.raises(AssertionError): dm.evolve_single(op, 2) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.evolve_single(op, 1) def test_evolve_single_success(self): @@ -116,16 +216,16 @@ def test_expectation_single_fail(self): # generate random 4 x 4 unitary matrix op = randobj.rand_unit(4) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.expectation_single(op, 2) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.expectation_single(op, 1) # wrong qubit indices op = randobj.rand_unit(2) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.expectation_single(op, -3) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.expectation_single(op, nqb + 3) def test_expectation_single_success(self): @@ -152,13 +252,13 @@ def test_expectation_single_success(self): psi1 = psi1.reshape(2**nqb) # watch out ordering. Expval unitary is cpx so psi1 on the right to match DM. - np.testing.assert_allclose(np.dot(psi.conjugate(), psi1), dm.expectation_single(op, target_qubit)) + assert np.allclose(np.dot(psi.conjugate(), psi1), dm.expectation_single(op, target_qubit)) def test_tensor_fail(self): dm = DensityMatrix(nqubit=1) - with self.assertRaises(TypeError): + with pytest.raises(TypeError): dm.tensor("hello") - with self.assertRaises(TypeError): + with pytest.raises(TypeError): dm.tensor(1) def test_tensor_without_data_success(self): @@ -183,15 +283,15 @@ def test_tensor_with_data_success(self): def test_cnot_fail(self): dm = DensityMatrix(nqubit=2) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.cnot((1, 1)) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.cnot((-1, 1)) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.cnot((1, -1)) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.cnot((1, 2)) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.cnot((2, 1)) def test_cnot_success(self): @@ -232,19 +332,19 @@ def test_cnot_success(self): psi = np.tensordot(CNOT_TENSOR, psi, ((2, 3), edge)) psi = np.moveaxis(psi, (0, 1), edge) expected_matrix = np.outer(psi, psi.conj()) - np.testing.assert_allclose(rho, expected_matrix) + assert np.allclose(rho, expected_matrix) def test_swap_fail(self): dm = DensityMatrix(nqubit=2) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.swap((1, 1)) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.swap((-1, 1)) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.swap((1, -1)) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.swap((1, 2)) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.swap((2, 1)) def test_swap_success(self): @@ -270,11 +370,11 @@ def test_swap_success(self): def test_entangle_fail(self): dm = DensityMatrix(nqubit=3) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.entangle((1, 1)) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.entangle(((1, 3))) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.entangle((0, 1, 2)) def test_entangle_success(self): @@ -307,7 +407,7 @@ def test_entangle_success(self): psi = np.tensordot(CZ_TENSOR, psi, ((2, 3), edge)) psi = np.moveaxis(psi, (0, 1), edge) expected_matrix = np.outer(psi, psi.conj()) - np.testing.assert_allclose(rho, expected_matrix) + assert np.allclose(rho, expected_matrix) def test_evolve_success(self): # single-qubit gate @@ -331,7 +431,7 @@ def test_evolve_success(self): dm.evolve(op, [i]) dm_single.evolve_single(op, i) - np.testing.assert_allclose(dm.rho, dm_single.rho) + assert np.allclose(dm.rho, dm_single.rho) # 2-qubit gate @@ -357,7 +457,7 @@ def test_evolve_success(self): psi = np.tensordot(op.reshape((2,) * 2 * N_qubits_op), psi, ((2, 3), edge)) psi = np.moveaxis(psi, (0, 1), edge) expected_matrix = np.outer(psi, psi.conj()) - np.testing.assert_allclose(rho, expected_matrix) + assert np.allclose(rho, expected_matrix) # 3-qubit gate N_qubits = np.random.randint(3, 5) @@ -382,7 +482,7 @@ def test_evolve_success(self): psi = np.tensordot(op.reshape((2,) * 2 * N_qubits_op), psi, ((3, 4, 5), targets)) psi = np.moveaxis(psi, (0, 1, 2), targets) expected_matrix = np.outer(psi, psi.conj()) - np.testing.assert_allclose(rho, expected_matrix) + assert np.allclose(rho, expected_matrix) def test_evolve_fail(self): # test on 3-qubit gate just in case. @@ -395,27 +495,27 @@ def test_evolve_fail(self): dm = DensityMatrix(nqubit=N_qubits) # dimension mismatch - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.evolve(op, (1, 1)) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.evolve(op, (0, 1, 2, 3)) # incorrect range - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.evolve(op, (-1, 0, 1)) # repeated index - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.evolve(op, (0, 1, 1)) # check not square matrix - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.evolve(np.random.rand(2, 3), (0, 1)) # check higher dimensional matrix - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.evolve(np.random.rand(2, 2, 3), (0, 1)) # check square but with incorrect dimension (non-qubit type) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): dm.evolve(np.random.rand(5, 5), (0, 1)) # TODO the test for normalization is done at initialization with data. Now check that all operations conserve the norm. @@ -430,10 +530,10 @@ def test_normalize(self): def test_ptrace_fail(self): dm = DensityMatrix(nqubit=0) - with self.assertRaises(AssertionError): + with pytest.raises(AssertionError): dm.ptrace((0,)) dm = DensityMatrix(nqubit=2) - with self.assertRaises(AssertionError): + with pytest.raises(AssertionError): dm.ptrace((2,)) def test_ptrace(self): @@ -502,8 +602,8 @@ def test_apply_dephasing_channel(self): + np.sqrt(prob) ** 2 * Ops.z @ rho_test @ Ops.z.conj().T ) - np.testing.assert_allclose(expected_dm.trace(), 1.0) - np.testing.assert_allclose(dm.rho, expected_dm) + assert np.allclose(expected_dm.trace(), 1.0) + assert np.allclose(dm.rho, expected_dm) N_qubits = np.random.randint(2, 5) @@ -552,8 +652,8 @@ def test_apply_dephasing_channel(self): ) ** 2 * np.outer(psi_evolvedb, psi_evolvedb.conj()) # compare - np.testing.assert_allclose(expected_dm.trace(), 1.0) - np.testing.assert_allclose(dm.rho, expected_dm) + assert np.allclose(expected_dm.trace(), 1.0) + assert np.allclose(dm.rho, expected_dm) def test_apply_depolarising_channel(self): # check on single qubit first @@ -586,8 +686,8 @@ def test_apply_depolarising_channel(self): + np.sqrt(prob / 3.0) ** 2 * Ops.z @ rho_test @ Ops.z.conj().T ) - np.testing.assert_allclose(expected_dm.trace(), 1.0) - np.testing.assert_allclose(dm.rho, expected_dm) + assert np.allclose(expected_dm.trace(), 1.0) + assert np.allclose(dm.rho, expected_dm) # chek against statevector backend by hand for now. # create random density matrix @@ -652,8 +752,8 @@ def test_apply_depolarising_channel(self): ) # compare - np.testing.assert_allclose(expected_dm.trace(), 1.0) - np.testing.assert_allclose(dm.rho, expected_dm) + assert np.allclose(expected_dm.trace(), 1.0) + assert np.allclose(dm.rho, expected_dm) def test_apply_random_channel_one_qubit(self): """ @@ -704,8 +804,8 @@ def test_apply_random_channel_one_qubit(self): expected_dm += elem["coef"] * np.conj(elem["coef"]) * np.outer(psi_evolved, np.conj(psi_evolved)) # compare - np.testing.assert_allclose(expected_dm.trace(), 1.0) - np.testing.assert_allclose(dm.rho, expected_dm) + assert np.allclose(expected_dm.trace(), 1.0) + assert np.allclose(dm.rho, expected_dm) def test_apply_random_channel_two_qubits(self): """ @@ -745,8 +845,8 @@ def test_apply_random_channel_two_qubits(self): psi_evolved = np.moveaxis(psi_evolved, (0, 1), qubits) expected_dm += elem["coef"] * np.conj(elem["coef"]) * np.outer(psi_evolved, np.conj(psi_evolved)) - np.testing.assert_allclose(expected_dm.trace(), 1.0) - np.testing.assert_allclose(dm.rho, expected_dm) + assert np.allclose(expected_dm.trace(), 1.0) + assert np.allclose(dm.rho, expected_dm) def test_apply_channel_fail(self): """ @@ -761,21 +861,77 @@ def test_apply_channel_fail(self): # build DensityMatrix dm = DensityMatrix(data=np.outer(psi, psi.conj())) - with self.assertRaises(TypeError): + with pytest.raises(TypeError): dm.apply_channel("a", [i]) -class DensityMatrixBackendTest(unittest.TestCase): +class TestDensityMatrixBackend(unittest.TestCase): """Test for DensityMatrixBackend class.""" + def setUp(self): + # set up the random numbers + self.rng = np.random.default_rng() # seed=422 + + circ = Circuit(1) + circ.h(0) + self.hadamardpattern = circ.transpile() + + self.nqb = self.rng.integers(2, 5) + # just want to test the initialization + self.depth = 1 + rand_circ = tests.random_circuit.get_rand_circuit(self.nqb, self.depth) + self.randpattern = rand_circ.transpile() + # print(self.randpattern, self.nqb, self.depth) + + # test initialization only + def test_init_success(self): + + # plus state (default) + backend = DensityMatrixBackend(self.hadamardpattern) + dm = DensityMatrix(nqubit=1) + assert np.allclose(dm.rho, backend.state.rho) + # assert backend.state.Nqubit == 1 + assert backend.state.dims() == (2, 2) + + # minus state + backend = DensityMatrixBackend(self.randpattern, input_state=graphix.states.BasicStates.MINUS) + dm = DensityMatrix(nqubit=self.nqb, data=graphix.states.BasicStates.MINUS) + assert np.allclose(dm.rho, backend.state.rho) + # assert backend.state.Nqubit == 1 + assert backend.state.dims() == (2**self.nqb, 2**self.nqb) + + rand_angles = self.rng.random(self.nqb) * 2 * np.pi + rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), self.nqb) + states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] + + expected_dm = DensityMatrix(data=states).rho + + backend = DensityMatrixBackend(self.randpattern, input_state=states) + dm = backend.state + + assert dm.dims() == (2**self.nqb, 2**self.nqb) + assert np.allclose(dm.rho, expected_dm) + assert backend.Nqubit == self.nqb + def test_init_fail(self): - with self.assertRaises(TypeError): + + rand_angles = self.rng.random(self.nqb + 1) * 2 * np.pi + rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), self.nqb + 1) + states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] + + # test init from State Iterable with incorrect size + with pytest.raises(ValueError): + backend = DensityMatrixBackend(self.randpattern, input_state=states) + + # don't provide required pattern argument + with pytest.raises(TypeError): DensityMatrixBackend() + def test_init_success(self): circ = Circuit(1) circ.rx(0, np.pi / 2) - pattern = circ.transpile().pattern + pattern = circ.transpile() backend = DensityMatrixBackend(pattern) assert backend.pattern == pattern assert backend.results == pattern.results @@ -783,13 +939,14 @@ def test_init_success(self): assert backend.Nqubit == 1 assert backend.max_qubit_num == 12 + def test_add_nodes(self): circ = Circuit(1) pattern = circ.transpile().pattern backend = DensityMatrixBackend(pattern) backend.add_nodes([1]) expected_matrix = np.array([0.25] * 16).reshape(4, 4) - np.testing.assert_allclose(backend.state.rho, expected_matrix) + assert np.allclose(backend.state.rho, expected_matrix) def test_entangle_nodes(self): circ = Circuit(1) @@ -798,10 +955,10 @@ def test_entangle_nodes(self): backend.add_nodes([1]) backend.entangle_nodes((0, 1)) expected_matrix = np.array([[1, 1, 1, -1], [1, 1, 1, -1], [1, 1, 1, -1], [-1, -1, -1, 1]]) / 4 - np.testing.assert_allclose(backend.state.rho, expected_matrix) + assert np.allclose(backend.state.rho, expected_matrix) backend.entangle_nodes((0, 1)) - np.testing.assert_allclose(backend.state.rho, np.array([0.25] * 16).reshape(4, 4)) + assert np.allclose(backend.state.rho, np.array([0.25] * 16).reshape(4, 4)) def test_measure(self): circ = Circuit(1) @@ -867,7 +1024,7 @@ def test_correct_byproduct(self): backend.finalize() psi = backend.state.psi - np.testing.assert_allclose(rho, np.outer(psi, psi.conj())) + assert np.allclose(rho, np.outer(psi, psi.conj())) if __name__ == "__main__": diff --git a/tests/test_pattern.py b/tests/test_pattern.py index 29caf7ccc..44f5289d9 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -6,8 +6,13 @@ import tests.random_circuit as rc from graphix.pattern import CommandNode, Pattern -from graphix.simulator import PatternSimulator from graphix.transpiler import Circuit +import graphix.states +from graphix.simulator import PatternSimulator +from graphix.sim.statevec import Statevec +from graphix.sim.density_matrix import DensityMatrix +import graphix.ops + SEED = 42 rc.set_seed(SEED) @@ -491,5 +496,41 @@ def assert_equal_edge(edge, ref): return False +# for testing with arbitrary inputs +# SV and DM backend +class TestPatternSim(unittest.TestCase): + def setUp(self): + # set up the random numbers + self.rng = np.random.default_rng() # seed=422 + + self.circ = Circuit(1) + self.circ.h(0) + self.hadamardpattern = self.circ.transpile() + + # self.nqb = self.rng.integers(2, 5) + # # just want to test the initialization + # self.depth = 1 + # rand_circ = tests.random_circuit.get_rand_circuit(self.nqb, self.depth) + # self.randpattern = rand_circ.transpile() + + def test_SV_sim(self): + nqb = 1 + rand_angles = self.rng.random(nqb) * 2 * np.pi + rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) + states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] + + out = self.hadamardpattern.simulate_pattern(backend="statevector", input_state=states) + + ref = Statevec(states) + ref.evolve_single(graphix.ops.Ops.h, 0) + + assert np.allclose(out.psi, ref.psi) + + out_circ = self.circ.simulate_statevector(input_state = states) + + assert np.allclose(out_circ.psi, ref.psi) + assert np.allclose(out.psi, out_circ.psi) + + if __name__ == "__main__": unittest.main() diff --git a/tests/test_statevec.py b/tests/test_statevec.py index c0793c561..efb7bf8f0 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -1,4 +1,5 @@ import unittest +import pytest import numpy as np from graphix.states import BasicStates, PlanarState @@ -8,8 +9,6 @@ import functools - - class TestStatevec(unittest.TestCase): """Test for Statevec class. Particularly new constructor.""" @@ -22,39 +21,39 @@ def setUp(self): # test injitializing one qubit in plus state def test_default_success(self): vec = Statevec(nqubit=1) - np.testing.assert_allclose(vec.psi, np.array([1, 1] / np.sqrt(2))) + assert np.allclose(vec.psi, np.array([1, 1] / np.sqrt(2))) # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 def test_basicstates_success(self): # minus vec = Statevec(nqubit=1, data=BasicStates.MINUS) - np.testing.assert_allclose(vec.psi, np.array([1, -1] / np.sqrt(2))) + assert np.allclose(vec.psi, np.array([1, -1] / np.sqrt(2))) # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # zero vec = Statevec(nqubit=1, data=BasicStates.ZERO) - np.testing.assert_allclose(vec.psi, np.array([1, 0]), rtol=0, atol=1e-15) + assert np.allclose(vec.psi, np.array([1, 0]), rtol=0, atol=1e-15) # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # one vec = Statevec(nqubit=1, data=BasicStates.ONE) - np.testing.assert_allclose(vec.psi, np.array([0, 1]), rtol=0, atol=1e-15) + assert np.allclose(vec.psi, np.array([0, 1]), rtol=0, atol=1e-15) # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # plus_i vec = Statevec(nqubit=1, data=BasicStates.PLUS_I) - np.testing.assert_allclose(vec.psi, np.array([1, 1j] / np.sqrt(2))) + assert np.allclose(vec.psi, np.array([1, 1j] / np.sqrt(2))) # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # minus_i vec = Statevec(nqubit=1, data=BasicStates.MINUS_I) - np.testing.assert_allclose(vec.psi, np.array([1, -1j] / np.sqrt(2))) - #assert vec.Nqubit == 1 + assert np.allclose(vec.psi, np.array([1, -1j] / np.sqrt(2))) + # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # even more tests? @@ -62,14 +61,14 @@ def test_default_tensor_success(self): nqb = self.rng.integers(2, 5) print(f"nqb is {nqb}") vec = Statevec(nqubit=nqb) - np.testing.assert_allclose(vec.psi, np.ones(((2,) * nqb)) / (np.sqrt(2)) ** nqb) + assert np.allclose(vec.psi, np.ones(((2,) * nqb)) / (np.sqrt(2)) ** nqb) # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb vec = Statevec(nqubit=nqb, data=BasicStates.MINUS_I) sv_list = [BasicStates.MINUS_I.get_statevector() for _ in range(nqb)] sv = functools.reduce(np.kron, sv_list) - np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) + assert np.allclose(vec.psi, sv.reshape((2,) * nqb)) # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb @@ -80,38 +79,37 @@ def test_default_tensor_success(self): vec = Statevec(nqubit=nqb, data=state) sv_list = [state.get_statevector() for _ in range(nqb)] sv = functools.reduce(np.kron, sv_list) - np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) + assert np.allclose(vec.psi, sv.reshape((2,) * nqb)) # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb # tensor of different states rand_angles = self.rng.random(nqb) * 2 * np.pi rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) - states = [PlanarState(plane = i, angle = j) for i, j in zip(rand_planes, rand_angles)] + states = [PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] vec = Statevec(nqubit=nqb, data=states) sv_list = [state.get_statevector() for state in states] sv = functools.reduce(np.kron, sv_list) - np.testing.assert_allclose(vec.psi, sv.reshape((2,) * nqb)) + assert np.allclose(vec.psi, sv.reshape((2,) * nqb)) # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb def test_data_success(self): nqb = self.rng.integers(2, 5) - l = 2 ** nqb + l = 2**nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) vec = Statevec(data=rand_vec) - np.testing.assert_allclose(vec.psi, rand_vec.reshape((2,) * nqb)) + assert np.allclose(vec.psi, rand_vec.reshape((2,) * nqb)) # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb - # fail: incorrect len def test_data_dim_fail(self): l = 5 rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): vec = Statevec(data=rand_vec) # with less qubit than number of qubits inferred from a correct state vect @@ -119,17 +117,17 @@ def test_data_dim_fail(self): # NOTE weird behaviour?? def test_data_dim_fail_mismatch(self): nqb = 3 - rand_vec = self.rng.random(2 ** nqb) + 1j * self.rng.random(2 ** nqb) + rand_vec = self.rng.random(2**nqb) + 1j * self.rng.random(2**nqb) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) - with self.assertRaises(ValueError): - vec = Statevec(nqubit = 2, data=rand_vec) + with pytest.raises(ValueError): + vec = Statevec(nqubit=2, data=rand_vec) # fail: not normalized def test_data_norm_fail(self): nqb = self.rng.integers(2, 5) - l = 2 ** nqb + l = 2**nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): vec = Statevec(data=rand_vec) def test_defaults_to_one(self): @@ -139,24 +137,24 @@ def test_defaults_to_one(self): # try copying Statevec input def test_copy_success(self): nqb = self.rng.integers(2, 5) - l = 2 ** nqb + l = 2**nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) test_vec = Statevec(data=rand_vec) # try to copy it vec = Statevec(data=test_vec) - np.testing.assert_allclose(vec.psi, test_vec.psi) + assert np.allclose(vec.psi, test_vec.psi) # assert vec.Nqubit == test_vec.Nqubit - assert len(vec.dims()) == len(test_vec.dims()) + assert len(vec.dims()) == len(test_vec.dims()) # try calling with incorrect number of qubits compared to inferred one def test_copy_fail(self): nqb = self.rng.integers(2, 5) - l = 2 ** nqb + l = 2**nqb rand_vec = self.rng.random(l) + 1j * self.rng.random(l) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) test_vec = Statevec(data=rand_vec) - with self.assertRaises(ValueError): + with pytest.raises(ValueError): vec = Statevec(nqubit=l - 1, data=test_vec) diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index e7ab78598..5ff03954d 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -1,4 +1,6 @@ import unittest +import pytest + from copy import deepcopy import numpy as np @@ -7,7 +9,8 @@ from graphix.sim.statevec import Statevec, meas_op, StatevectorBackend import graphix.pauli -class TestStatevec(unittest.TestCase): + +class TestStatevecBackend(unittest.TestCase): def test_remove_one_qubit(self): n = 10 k = 3 @@ -25,7 +28,7 @@ def test_remove_one_qubit(self): np.testing.assert_almost_equal(np.abs(sv.psi.flatten().dot(sv2.psi.flatten().conj())), 1) - #TODO This is a weird test! + # TODO This is a weird test! def test_measurement_into_each_XYZ_basis(self): n = 3 k = 0 @@ -48,10 +51,14 @@ def test_measurement_into_minus_state(self): m_op = np.outer(BasicStates.MINUS.get_statevector(), BasicStates.MINUS.get_statevector().T.conjugate()) sv = Statevec(nqubit=n) sv.evolve(m_op, [k]) - with self.assertRaises(AssertionError): + with pytest.raises(AssertionError): sv.remove_qubit(k) + class TestStatevecNew(unittest.TestCase): + + # more tests not really needed since redundant with Statevec constructor tests + def setUp(self): # set up the random numbers self.rng = np.random.default_rng() # seed=422 @@ -66,30 +73,29 @@ def test_init_success(self): # plus state (default) backend = StatevectorBackend(self.hadamardpattern) vec = Statevec(nqubit=1) - np.testing.assert_allclose(vec.psi, backend.state.psi) + assert np.allclose(vec.psi, backend.state.psi) # assert backend.state.Nqubit == 1 assert len(backend.state.dims()) == 1 - # minus state - backend = StatevectorBackend(self.hadamardpattern, input_state = BasicStates.MINUS) + # minus state + backend = StatevectorBackend(self.hadamardpattern, input_state=BasicStates.MINUS) vec = Statevec(nqubit=1, data=BasicStates.MINUS) - np.testing.assert_allclose(vec.psi, backend.state.psi) + assert np.allclose(vec.psi, backend.state.psi) # assert backend.state.Nqubit == 1 assert len(backend.state.dims()) == 1 # random planar state rand_angle = self.rng.random() * 2 * np.pi rand_plane = self.rng.choice(np.array([i for i in graphix.pauli.Plane])) - state = PlanarState(plane = rand_plane, angle = rand_angle) - backend = StatevectorBackend(self.hadamardpattern, input_state = state) + state = PlanarState(plane=rand_plane, angle=rand_angle) + backend = StatevectorBackend(self.hadamardpattern, input_state=state) vec = Statevec(nqubit=1, data=state) - np.testing.assert_allclose(vec.psi, backend.state.psi) + assert np.allclose(vec.psi, backend.state.psi) # assert backend.state.Nqubit == 1 assert len(backend.state.dims()) == 1 + # data input and Statevec input - - def test_init_fail(self): # incorrect number of dimensions for State input # only one input node, two states provided @@ -99,18 +105,10 @@ def test_init_fail(self): rand_angle = self.rng.random(2) * 2 * np.pi rand_plane = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), 2) - state = PlanarState(plane = rand_plane[0], angle = rand_angle[0]) - state2 = PlanarState(plane = rand_plane[1], angle = rand_angle[1]) - with self.assertRaises(ValueError): - StatevectorBackend(self.hadamardpattern, input_state = [state, state2]) - # vec = Statevec(nqubit=1, state = state) - # np.testing.assert_allclose(vec.psi, backend.state.psi) - # # assert backend.state.Nqubit == 1 - # assert len(backend.state.dims()) == 1 - - - - + state = PlanarState(plane=rand_plane[0], angle=rand_angle[0]) + state2 = PlanarState(plane=rand_plane[1], angle=rand_angle[1]) + with pytest.raises(ValueError): + StatevectorBackend(self.hadamardpattern, input_state=[state, state2]) if __name__ == "__main__": From 10169b406e05e55e72efd3ca7902edc3296a0b01 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Mon, 29 Apr 2024 17:24:09 +0200 Subject: [PATCH 15/48] update test_pattern.py Add test to check simulation at the circuit level is the same as at pattern level in the statevec backend; DM to come --- tests/test_pattern.py | 35 +++++++++++++++++------------------ 1 file changed, 17 insertions(+), 18 deletions(-) diff --git a/tests/test_pattern.py b/tests/test_pattern.py index 44f5289d9..efc508031 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -12,6 +12,7 @@ from graphix.sim.statevec import Statevec from graphix.sim.density_matrix import DensityMatrix import graphix.ops +import tests.random_circuit SEED = 42 @@ -507,30 +508,28 @@ def setUp(self): self.circ.h(0) self.hadamardpattern = self.circ.transpile() - # self.nqb = self.rng.integers(2, 5) - # # just want to test the initialization - # self.depth = 1 - # rand_circ = tests.random_circuit.get_rand_circuit(self.nqb, self.depth) - # self.randpattern = rand_circ.transpile() + self.nqb = self.rng.integers(2, 5) + # just want to test the initialization + self.depth = 2 + self.rand_circ = tests.random_circuit.get_rand_circuit(self.nqb, self.depth) + self.randpattern = self.rand_circ.transpile() - def test_SV_sim(self): - nqb = 1 - rand_angles = self.rng.random(nqb) * 2 * np.pi - rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) + def test_sv_sim(self): + rand_angles = self.rng.random(self.nqb) * 2 * np.pi + rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), self.nqb) states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] - out = self.hadamardpattern.simulate_pattern(backend="statevector", input_state=states) + out = self.randpattern.simulate_pattern(backend = "statevector", input_state = states) - ref = Statevec(states) - ref.evolve_single(graphix.ops.Ops.h, 0) + out_circ = self.rand_circ.simulate_statevector(input_state = states) - assert np.allclose(out.psi, ref.psi) - - out_circ = self.circ.simulate_statevector(input_state = states) - - assert np.allclose(out_circ.psi, ref.psi) - assert np.allclose(out.psi, out_circ.psi) + # MBQC is up to a global phase! + np.testing.assert_almost_equal(np.abs(np.dot(out.psi.flatten().conjugate(), out_circ.psi.flatten())), 1) + # assert np.allclose(out.psi, out_circ.psi) + + def test_dm_sim(self): + pass if __name__ == "__main__": unittest.main() From db6d54e3ef256502c06b02d67d82c4e26d8fdae2 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Fri, 3 May 2024 14:12:17 +0200 Subject: [PATCH 16/48] Circuit measure handling --- tests/test_density_matrix.py | 4 ++-- tests/test_pattern.py | 4 ++-- tests/test_statevec_backend.py | 2 +- 3 files changed, 5 insertions(+), 5 deletions(-) diff --git a/tests/test_density_matrix.py b/tests/test_density_matrix.py index 7a53410f5..51022ae3a 100644 --- a/tests/test_density_matrix.py +++ b/tests/test_density_matrix.py @@ -880,7 +880,7 @@ def setUp(self): # just want to test the initialization self.depth = 1 rand_circ = tests.random_circuit.get_rand_circuit(self.nqb, self.depth) - self.randpattern = rand_circ.transpile() + self.randpattern = rand_circ.transpile().pattern # print(self.randpattern, self.nqb, self.depth) # test initialization only @@ -931,7 +931,7 @@ def test_init_fail(self): def test_init_success(self): circ = Circuit(1) circ.rx(0, np.pi / 2) - pattern = circ.transpile() + pattern = circ.transpile().pattern backend = DensityMatrixBackend(pattern) assert backend.pattern == pattern assert backend.results == pattern.results diff --git a/tests/test_pattern.py b/tests/test_pattern.py index efc508031..96a88da61 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -512,7 +512,7 @@ def setUp(self): # just want to test the initialization self.depth = 2 self.rand_circ = tests.random_circuit.get_rand_circuit(self.nqb, self.depth) - self.randpattern = self.rand_circ.transpile() + self.randpattern = self.rand_circ.transpile().pattern def test_sv_sim(self): rand_angles = self.rng.random(self.nqb) * 2 * np.pi @@ -521,7 +521,7 @@ def test_sv_sim(self): out = self.randpattern.simulate_pattern(backend = "statevector", input_state = states) - out_circ = self.rand_circ.simulate_statevector(input_state = states) + out_circ = self.rand_circ.simulate_statevector(input_state = states).statevec # MBQC is up to a global phase! np.testing.assert_almost_equal(np.abs(np.dot(out.psi.flatten().conjugate(), out_circ.psi.flatten())), 1) diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index 5ff03954d..ebdc3935c 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -65,7 +65,7 @@ def setUp(self): circ = Circuit(1) circ.h(0) - self.hadamardpattern = circ.transpile() + self.hadamardpattern = circ.transpile().pattern # test initialization only def test_init_success(self): From 25aa33840f08bdf66199427b6c45083841d4171b Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Fri, 3 May 2024 14:14:06 +0200 Subject: [PATCH 17/48] Some simplifications --- graphix/linalg_validations.py | 4 +--- graphix/ops.py | 12 ------------ graphix/sim/density_matrix.py | 5 ++--- graphix/sim/statevec.py | 8 ++------ graphix/sim/tensornet.py | 3 +-- graphix/transpiler.py | 2 +- graphix/types.py | 2 +- tests/test_density_matrix.py | 6 +++--- tests/test_statevec.py | 16 +--------------- 9 files changed, 12 insertions(+), 46 deletions(-) diff --git a/graphix/linalg_validations.py b/graphix/linalg_validations.py index e4fa101c1..b3cd0a525 100644 --- a/graphix/linalg_validations.py +++ b/graphix/linalg_validations.py @@ -21,13 +21,11 @@ def check_square(matrix: np.ndarray) -> bool: def check_psd(matrix: np.ndarray, tol: float = 1e-15) -> bool: """ check if a density matrix is positive semidefinite by diagonalizing. - After check_square and check_hermitian (osef) so that it already is square with power of 2 dimension. - Parameters ---------- matrix : np.ndarray - matrix to check. Normally already square and 2**n x 2**n + matrix to check tol : float tolerance on the small negatives. Default 1e-15. """ diff --git a/graphix/ops.py b/graphix/ops.py index d8aa4b9ae..b74366333 100644 --- a/graphix/ops.py +++ b/graphix/ops.py @@ -7,18 +7,6 @@ import numpy as np -# TODO modify that - -# Everywhere this is called. use StateVec(State)) -# class States: -# plus = np.array([1.0 / np.sqrt(2), 1.0 / np.sqrt(2)]) # plus -# minus = np.array([1.0 / np.sqrt(2), -1.0 / np.sqrt(2)]) # minus -# zero = np.array([1.0, 0.0]) # zero -# one = np.array([0.0, 1.0]) # one -# iplus = np.array([1.0 / np.sqrt(2), 1.0j / np.sqrt(2)]) # +1 eigenstate of Pauli Y -# iminus = np.array([1.0 / np.sqrt(2), -1.0j / np.sqrt(2)]) # -1 eigenstate of Pauli Y -# vec = [plus, minus, zero, one, iplus, iminus] - class Ops: """Basic single- and two-qubits operators""" diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index cb5384106..6abbb67ee 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -9,7 +9,6 @@ import numbers import numpy as np -import pydantic import typing from graphix.linalg_validations import check_square, check_hermitian, check_unit_trace, check_psd @@ -37,7 +36,7 @@ class DensityMatrix: def __init__( self, data: typing.Optional[Data] = graphix.states.BasicStates.PLUS, - nqubit: typing.Optional[graphix.types.PositiveInt] = None, + nqubit: typing.Optional[graphix.types.PositiveOrNullInt] = None, ): """ @@ -49,7 +48,7 @@ def __init__( nqubit : int Number of qubits. Default is 1. If both `data` and `nqubit` are specified, consistency is checked. """ - pydantic.TypeAdapter(typing.Optional[graphix.types.PositiveInt]).validate_python(nqubit) + assert nqubit is None or isinstance(nqubit, numbers.Integral) and nqubit >= 0 def check_size_consistency(mat): if nqubit is not None and mat.shape != (2**nqubit, 2**nqubit): diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index c31090047..cc42b06fd 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -215,7 +215,7 @@ class Statevec: def __init__( self, data: typing.Optional[SV_Data] = graphix.states.BasicStates.PLUS, - nqubit: typing.Optional[graphix.types.PositiveInt] = None, + nqubit: typing.Optional[graphix.types.PositiveOrNullInt] = None, ): """Initialize statevector @@ -234,7 +234,7 @@ def __init__( Defaults to |+> states and 1 qubit. If nqubit > 1 and only one state : tensor all of them. Use the tensor method instead of hard code. """ - pydantic.TypeAdapter(typing.Optional[graphix.types.PositiveInt]).validate_python(nqubit) + assert nqubit is None or isinstance(nqubit, numbers.Integral) and nqubit >= 0 if isinstance(data, Statevec): # assert nqubit is None or len(state.flatten()) == 2**nqubit @@ -428,12 +428,8 @@ def tensor(self, other): psi_self = self.psi.flatten() psi_other = other.psi.flatten() - # NOTE on tensor form not vector - # deprecated total_num = len(self.dims()) + len(other.dims()) - # self.Nqubit += other.Nqubit self.psi = np.kron(psi_self, psi_other).reshape((2,) * total_num) - # self.Nqubit = len(self.dims()) def CNOT(self, qubits): """apply CNOT diff --git a/graphix/sim/tensornet.py b/graphix/sim/tensornet.py index d0475962c..a40c1b263 100644 --- a/graphix/sim/tensornet.py +++ b/graphix/sim/tensornet.py @@ -6,13 +6,12 @@ from quimb.tensor import Tensor, TensorNetwork from graphix.clifford import CLIFFORD, CLIFFORD_CONJ, CLIFFORD_MUL -from graphix.ops import Ops, States +from graphix.ops import Ops from graphix.states import BasicStates import string from copy import deepcopy - class TensorNetworkBackend: """Tensor Network Simulator for MBQC diff --git a/graphix/transpiler.py b/graphix/transpiler.py index ff87cf52e..315fc857f 100644 --- a/graphix/transpiler.py +++ b/graphix/transpiler.py @@ -1362,7 +1362,7 @@ def _sort_outputs(self, pattern: Pattern, output_nodes: Sequence[int]): elif cmd[1] in old_out: cmd[1] = output_nodes[old_out.index(cmd[1])] - def simulate_statevector(self, input_state: Optional[SV_Data] = None): + def simulate_statevector(self, input_state: Optional[SV_Data] = None) -> SimulateResult: """Run statevector simulation of the gate sequence, using graphix.Statevec Parameters diff --git a/graphix/types.py b/graphix/types.py index 165c40e7e..c761768d9 100644 --- a/graphix/types.py +++ b/graphix/types.py @@ -1,7 +1,7 @@ import annotated_types import typing_extensions -PositiveInt = typing_extensions.Annotated[int, annotated_types.Ge(0)] # includes 0 +PositiveOrNullInt = typing_extensions.Annotated[int, annotated_types.Ge(0)] # includes 0 def check_list_elements(l, ty): diff --git a/tests/test_density_matrix.py b/tests/test_density_matrix.py index 51022ae3a..049117f80 100644 --- a/tests/test_density_matrix.py +++ b/tests/test_density_matrix.py @@ -26,11 +26,11 @@ def setUp(self): self.rng = np.random.default_rng() def test_init_without_data_fail(self): - with pytest.raises(pydantic.ValidationError): + with pytest.raises(AssertionError): DensityMatrix(nqubit=-2) - with pytest.raises(pydantic.ValidationError): + with pytest.raises(AssertionError): DensityMatrix(nqubit="hello") - with pytest.raises(pydantic.ValidationError): + with pytest.raises(AssertionError): DensityMatrix(nqubit=[]) def test_init_with_invalid_data_fail(self): diff --git a/tests/test_statevec.py b/tests/test_statevec.py index efb7bf8f0..e822f5efa 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -22,38 +22,32 @@ def setUp(self): def test_default_success(self): vec = Statevec(nqubit=1) assert np.allclose(vec.psi, np.array([1, 1] / np.sqrt(2))) - # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 def test_basicstates_success(self): # minus vec = Statevec(nqubit=1, data=BasicStates.MINUS) assert np.allclose(vec.psi, np.array([1, -1] / np.sqrt(2))) - # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # zero vec = Statevec(nqubit=1, data=BasicStates.ZERO) assert np.allclose(vec.psi, np.array([1, 0]), rtol=0, atol=1e-15) - # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # one vec = Statevec(nqubit=1, data=BasicStates.ONE) assert np.allclose(vec.psi, np.array([0, 1]), rtol=0, atol=1e-15) - # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # plus_i vec = Statevec(nqubit=1, data=BasicStates.PLUS_I) assert np.allclose(vec.psi, np.array([1, 1j] / np.sqrt(2))) - # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # minus_i vec = Statevec(nqubit=1, data=BasicStates.MINUS_I) assert np.allclose(vec.psi, np.array([1, -1j] / np.sqrt(2))) - # assert vec.Nqubit == 1 assert len(vec.dims()) == 1 # even more tests? @@ -62,14 +56,12 @@ def test_default_tensor_success(self): print(f"nqb is {nqb}") vec = Statevec(nqubit=nqb) assert np.allclose(vec.psi, np.ones(((2,) * nqb)) / (np.sqrt(2)) ** nqb) - # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb vec = Statevec(nqubit=nqb, data=BasicStates.MINUS_I) sv_list = [BasicStates.MINUS_I.get_statevector() for _ in range(nqb)] sv = functools.reduce(np.kron, sv_list) assert np.allclose(vec.psi, sv.reshape((2,) * nqb)) - # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb # tensor of same state @@ -80,7 +72,6 @@ def test_default_tensor_success(self): sv_list = [state.get_statevector() for _ in range(nqb)] sv = functools.reduce(np.kron, sv_list) assert np.allclose(vec.psi, sv.reshape((2,) * nqb)) - # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb # tensor of different states @@ -91,7 +82,6 @@ def test_default_tensor_success(self): sv_list = [state.get_statevector() for state in states] sv = functools.reduce(np.kron, sv_list) assert np.allclose(vec.psi, sv.reshape((2,) * nqb)) - # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb def test_data_success(self): @@ -101,7 +91,6 @@ def test_data_success(self): rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) vec = Statevec(data=rand_vec) assert np.allclose(vec.psi, rand_vec.reshape((2,) * nqb)) - # assert vec.Nqubit == nqb assert len(vec.dims()) == nqb # fail: incorrect len @@ -112,9 +101,7 @@ def test_data_dim_fail(self): with pytest.raises(ValueError): vec = Statevec(data=rand_vec) - # with less qubit than number of qubits inferred from a correct state vect - # returns a truncated statevec that is hence not normalized - # NOTE weird behaviour?? + # fail: with less qubit than number of qubits inferred from a correct state vect def test_data_dim_fail_mismatch(self): nqb = 3 rand_vec = self.rng.random(2**nqb) + 1j * self.rng.random(2**nqb) @@ -145,7 +132,6 @@ def test_copy_success(self): vec = Statevec(data=test_vec) assert np.allclose(vec.psi, test_vec.psi) - # assert vec.Nqubit == test_vec.Nqubit assert len(vec.dims()) == len(test_vec.dims()) # try calling with incorrect number of qubits compared to inferred one From 7d488c378b1bf360130ebe7ea02bbea599e69693 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Fri, 3 May 2024 15:04:18 +0200 Subject: [PATCH 18/48] Update CHANGELOG.md, isort, black --- CHANGELOG.md | 6 ++++++ examples/noisy_mbqc.py | 1 - graphix/linalg_validations.py | 2 -- graphix/random_objects.py | 3 ++- graphix/sim/density_matrix.py | 16 +++++++--------- graphix/sim/statevec.py | 13 ++++++------- graphix/sim/tensornet.py | 2 -- graphix/states.py | 5 ++++- graphix/transpiler.py | 4 +--- tests/test_density_matrix.py | 15 +++++---------- tests/test_kraus.py | 4 ---- tests/test_pattern.py | 19 +++++++++---------- tests/test_random_utilities.py | 8 -------- tests/test_statevec.py | 7 ++++--- tests/test_statevec_backend.py | 10 ++++------ tests/test_tnsim.py | 2 +- 16 files changed, 49 insertions(+), 68 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 4b33e8ec9..6e3c2212f 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -14,6 +14,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 be used in any subsequent gate. - Added `gflow.find_pauliflow`, `gflow.verify_pauliflow` and `pauliflow_from_pattern` methods (#117) +- Allow arbitrary states for initializing input nodes in state vector + and density matrix backends. + ### Fixed ### Changed @@ -30,6 +33,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - Patterns are now allowed to measure all their nodes, and have an empty output set. +- Basic states are now defined in `states.BasicStates` and no longer + in `ops.States`. + ## [0.2.11] - 2024-03-16 ### Added diff --git a/examples/noisy_mbqc.py b/examples/noisy_mbqc.py index 01fce8186..f6413bba3 100644 --- a/examples/noisy_mbqc.py +++ b/examples/noisy_mbqc.py @@ -43,7 +43,6 @@ class NoisyGraphState(NoiseModel): - def __init__(self, p_z=0.1): self.p_z = p_z diff --git a/graphix/linalg_validations.py b/graphix/linalg_validations.py index b3cd0a525..4e914ebcf 100644 --- a/graphix/linalg_validations.py +++ b/graphix/linalg_validations.py @@ -68,7 +68,6 @@ def check_data_normalization(data: Union[list, tuple, np.ndarray]) -> bool: def check_data_dims(data: Union[list, tuple, np.ndarray]) -> bool: - # convert to set to remove duplicates dims = set([i["operator"].shape for i in data]) @@ -83,7 +82,6 @@ def check_data_dims(data: Union[list, tuple, np.ndarray]) -> bool: def check_data_values_type(data: Union[list, tuple, np.ndarray]) -> bool: - if not all( isinstance(i, dict) for i in data ): # ni liste ni ensemble mais iterable (lazy) pas stocké, executé au besoin diff --git a/graphix/random_objects.py b/graphix/random_objects.py index 2ed2e61bc..b811437e7 100644 --- a/graphix/random_objects.py +++ b/graphix/random_objects.py @@ -1,5 +1,6 @@ -import numpy as np from typing import Optional + +import numpy as np import scipy.linalg from scipy.stats import unitary_group diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index 6abbb67ee..5f8b3fdf7 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -3,22 +3,21 @@ Simulate MBQC with density matrix representation. """ -from copy import deepcopy -import typing import functools import numbers +import typing +from copy import deepcopy import numpy as np -import typing -from graphix.linalg_validations import check_square, check_hermitian, check_unit_trace, check_psd -from graphix.channels import KrausChannel -from graphix.ops import Ops -from graphix.clifford import CLIFFORD -from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, meas_op, Statevec import graphix.sim.base_backend import graphix.states import graphix.types +from graphix.channels import KrausChannel +from graphix.clifford import CLIFFORD +from graphix.linalg_validations import check_hermitian, check_psd, check_square, check_unit_trace +from graphix.ops import Ops +from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, Statevec, meas_op Data = typing.Union[ graphix.states.State, @@ -38,7 +37,6 @@ def __init__( data: typing.Optional[Data] = graphix.states.BasicStates.PLUS, nqubit: typing.Optional[graphix.types.PositiveOrNullInt] = None, ): - """ rewrite! Parameters diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index cc42b06fd..dbcfa5c53 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -1,18 +1,18 @@ -from copy import deepcopy +import functools import numbers import typing +import warnings +from copy import deepcopy import numpy as np -import functools import pydantic -import warnings +import graphix.pauli import graphix.sim.base_backend -from graphix.clifford import CLIFFORD, CLIFFORD_CONJ, CLIFFORD_MUL -from graphix.ops import Ops import graphix.states -import graphix.pauli import graphix.types +from graphix.clifford import CLIFFORD, CLIFFORD_CONJ, CLIFFORD_MUL +from graphix.ops import Ops # Python >= 3.9 # from collections.abc import Iterable # or use Protocols? @@ -217,7 +217,6 @@ def __init__( data: typing.Optional[SV_Data] = graphix.states.BasicStates.PLUS, nqubit: typing.Optional[graphix.types.PositiveOrNullInt] = None, ): - """Initialize statevector Parameters diff --git a/graphix/sim/tensornet.py b/graphix/sim/tensornet.py index a40c1b263..985b0d187 100644 --- a/graphix/sim/tensornet.py +++ b/graphix/sim/tensornet.py @@ -8,8 +8,6 @@ from graphix.clifford import CLIFFORD, CLIFFORD_CONJ, CLIFFORD_MUL from graphix.ops import Ops from graphix.states import BasicStates -import string -from copy import deepcopy class TensorNetworkBackend: diff --git a/graphix/states.py b/graphix/states.py index 81ecd4877..c64a48510 100644 --- a/graphix/states.py +++ b/graphix/states.py @@ -2,10 +2,13 @@ quantum states and operators """ +import abc + import numpy as np import pydantic + import graphix.pauli -import abc + # generic class State for all States class State(abc.ABC): diff --git a/graphix/transpiler.py b/graphix/transpiler.py index 315fc857f..57da6bd16 100644 --- a/graphix/transpiler.py +++ b/graphix/transpiler.py @@ -16,10 +16,8 @@ import graphix.sim.base_backend from graphix.ops import Ops from graphix.pattern import Pattern -import graphix.sim.base_backend +from graphix.sim.density_matrix import Data, DensityMatrix from graphix.sim.statevec import Statevec, SV_Data -from graphix.sim.density_matrix import DensityMatrix, Data -import graphix.pauli @dataclasses.dataclass diff --git a/tests/test_density_matrix.py b/tests/test_density_matrix.py index 049117f80..1014affcf 100644 --- a/tests/test_density_matrix.py +++ b/tests/test_density_matrix.py @@ -1,21 +1,21 @@ +import functools import random import unittest -import pydantic from copy import deepcopy import numpy as np +import pydantic import pytest +import graphix.pauli import graphix.random_objects as randobj +import graphix.states +import tests.random_circuit from graphix import Circuit from graphix.channels import KrausChannel, dephasing_channel, depolarising_channel from graphix.ops import Ops from graphix.sim.density_matrix import DensityMatrix, DensityMatrixBackend from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, Statevec, StatevectorBackend -import graphix.states -import graphix.pauli -import tests.random_circuit -import functools class TestDensityMatrix(unittest.TestCase): @@ -116,7 +116,6 @@ def test_init_with_state_sucess(self): assert np.allclose(dm.rho, expected_dm) def test_init_with_state_fail(self): - nqb = 2 rand_angles = self.rng.random(nqb) * 2 * np.pi rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) @@ -885,7 +884,6 @@ def setUp(self): # test initialization only def test_init_success(self): - # plus state (default) backend = DensityMatrixBackend(self.hadamardpattern) dm = DensityMatrix(nqubit=1) @@ -914,7 +912,6 @@ def test_init_success(self): assert backend.Nqubit == self.nqb def test_init_fail(self): - rand_angles = self.rng.random(self.nqb + 1) * 2 * np.pi rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), self.nqb + 1) states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] @@ -927,7 +924,6 @@ def test_init_fail(self): with pytest.raises(TypeError): DensityMatrixBackend() - def test_init_success(self): circ = Circuit(1) circ.rx(0, np.pi / 2) @@ -939,7 +935,6 @@ def test_init_success(self): assert backend.Nqubit == 1 assert backend.max_qubit_num == 12 - def test_add_nodes(self): circ = Circuit(1) pattern = circ.transpile().pattern diff --git a/tests/test_kraus.py b/tests/test_kraus.py index 732c6dc7e..1d6e1a303 100644 --- a/tests/test_kraus.py +++ b/tests/test_kraus.py @@ -123,7 +123,6 @@ def test_init_with_data_fail(self): randobj.rand_channel_kraus(dim=2**2, rank=20) def test_dephasing_channel(self): - prob = np.random.rand() data = [ {"coef": np.sqrt(1 - prob), "operator": np.array([[1.0, 0.0], [0.0, 1.0]])}, @@ -140,7 +139,6 @@ def test_dephasing_channel(self): np.testing.assert_allclose(dephase_channel.kraus_ops[i]["operator"], data[i]["operator"]) def test_depolarising_channel(self): - prob = np.random.rand() data = [ {"coef": np.sqrt(1 - prob), "operator": np.eye(2)}, @@ -161,7 +159,6 @@ def test_depolarising_channel(self): np.testing.assert_allclose(depol_channel.kraus_ops[i]["operator"], data[i]["operator"]) def test_2_qubit_depolarising_channel(self): - prob = np.random.rand() data = [ {"coef": np.sqrt(1 - prob), "operator": np.kron(np.eye(2), np.eye(2))}, @@ -194,7 +191,6 @@ def test_2_qubit_depolarising_channel(self): np.testing.assert_allclose(depol_channel_2_qubit.kraus_ops[i]["operator"], data[i]["operator"]) def test_2_qubit_depolarising_tensor_channel(self): - prob = np.random.rand() data = [ {"coef": 1 - prob, "operator": np.kron(np.eye(2), np.eye(2))}, diff --git a/tests/test_pattern.py b/tests/test_pattern.py index 96a88da61..a2bea581f 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -4,16 +4,15 @@ import numpy as np from parameterized import parameterized +import graphix.ops +import graphix.states +import tests.random_circuit import tests.random_circuit as rc from graphix.pattern import CommandNode, Pattern -from graphix.transpiler import Circuit -import graphix.states -from graphix.simulator import PatternSimulator -from graphix.sim.statevec import Statevec from graphix.sim.density_matrix import DensityMatrix -import graphix.ops -import tests.random_circuit - +from graphix.sim.statevec import Statevec +from graphix.simulator import PatternSimulator +from graphix.transpiler import Circuit SEED = 42 rc.set_seed(SEED) @@ -519,17 +518,17 @@ def test_sv_sim(self): rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), self.nqb) states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] - out = self.randpattern.simulate_pattern(backend = "statevector", input_state = states) + out = self.randpattern.simulate_pattern(backend="statevector", input_state=states) - out_circ = self.rand_circ.simulate_statevector(input_state = states).statevec + out_circ = self.rand_circ.simulate_statevector(input_state=states).statevec # MBQC is up to a global phase! np.testing.assert_almost_equal(np.abs(np.dot(out.psi.flatten().conjugate(), out_circ.psi.flatten())), 1) # assert np.allclose(out.psi, out_circ.psi) - def test_dm_sim(self): pass + if __name__ == "__main__": unittest.main() diff --git a/tests/test_random_utilities.py b/tests/test_random_utilities.py index 089d05f4b..e9f6654d5 100644 --- a/tests/test_random_utilities.py +++ b/tests/test_random_utilities.py @@ -26,7 +26,6 @@ def test_rand_unit(self): np.testing.assert_allclose(tmp.conj().T @ tmp, np.eye(d), atol=1e-15) def test_random_channel_success(self): - nqb = np.random.randint(1, 5) dim = 2**nqb # np.random.randint(2, 8) @@ -51,7 +50,6 @@ def test_random_channel_success(self): assert channel.is_normalized def test_random_channel_fail(self): - # incorrect rank type with self.assertRaises(TypeError): mychannel = randobj.rand_channel_kraus(dim=2**2, rank=3.0) @@ -61,7 +59,6 @@ def test_random_channel_fail(self): mychannel = randobj.rand_channel_kraus(dim=2**2, rank=0) def test_rand_gauss_cpx(self): - nsample = int(1e4) dim = np.random.randint(2, 20) @@ -72,7 +69,6 @@ def test_rand_gauss_cpx(self): assert list(dimset)[0] == (dim, dim) def test_check_psd_success(self): - # Generate a random mixed state from state vectors with same probability # We know this is PSD @@ -92,7 +88,6 @@ def test_check_psd_success(self): assert check_psd(dm) def test_check_psd_fail(self): - # not hermitian # don't use dim = 2, too easy to have a PSD matrix. # NOTE useless test since eigvalsh treats the matrix as hermitian and takes only the L or U part @@ -128,7 +123,6 @@ def test_rand_dm_fail(self): dm = randobj.rand_dm(2 ** np.random.randint(2, 5) + 1) def test_rand_dm_rank(self): - rk = 3 dm = randobj.rand_dm(2 ** np.random.randint(2, 5), rank=rk) @@ -156,7 +150,6 @@ def test_pauli_tensor_ops(self): assert np.all(dims == (2**nqb, 2**nqb)) def test_pauli_tensor_ops_fail(self): - with self.assertRaises(TypeError): Pauli_tensor_ops = Ops.build_tensor_Pauli_ops(np.random.randint(2, 6) + 0.5) @@ -164,7 +157,6 @@ def test_pauli_tensor_ops_fail(self): Pauli_tensor_ops = Ops.build_tensor_Pauli_ops(0) def test_random_pauli_channel_success(self): - nqb = np.random.randint(2, 6) rk = np.random.randint(1, 2**nqb + 1) Pauli_channel = randobj.rand_Pauli_channel_kraus(dim=2**nqb, rank=rk) # default is full rank diff --git a/tests/test_statevec.py b/tests/test_statevec.py index e822f5efa..ce01f89df 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -1,12 +1,13 @@ +import functools import unittest -import pytest + import numpy as np +import pytest -from graphix.states import BasicStates, PlanarState import graphix.pauli import graphix.random_objects as randobj from graphix.sim.statevec import Statevec -import functools +from graphix.states import BasicStates, PlanarState class TestStatevec(unittest.TestCase): diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index ebdc3935c..ef9828487 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -1,13 +1,13 @@ import unittest -import pytest - from copy import deepcopy import numpy as np +import pytest + +import graphix.pauli from graphix import Circuit +from graphix.sim.statevec import Statevec, StatevectorBackend, meas_op from graphix.states import BasicStates, PlanarState -from graphix.sim.statevec import Statevec, meas_op, StatevectorBackend -import graphix.pauli class TestStatevecBackend(unittest.TestCase): @@ -56,7 +56,6 @@ def test_measurement_into_minus_state(self): class TestStatevecNew(unittest.TestCase): - # more tests not really needed since redundant with Statevec constructor tests def setUp(self): @@ -69,7 +68,6 @@ def setUp(self): # test initialization only def test_init_success(self): - # plus state (default) backend = StatevectorBackend(self.hadamardpattern) vec = Statevec(nqubit=1) diff --git a/tests/test_tnsim.py b/tests/test_tnsim.py index 7ed7bb0b6..1d156ee47 100644 --- a/tests/test_tnsim.py +++ b/tests/test_tnsim.py @@ -7,8 +7,8 @@ import tests.random_circuit as rc from graphix.clifford import CLIFFORD from graphix.ops import Ops -from graphix.states import BasicStates from graphix.sim.tensornet import MBQCTensorNet, gen_str +from graphix.states import BasicStates from graphix.transpiler import Circuit SEED = 42 From 209e77cd790464e5139443d790e5115787dc3528 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Fri, 3 May 2024 15:18:32 +0200 Subject: [PATCH 19/48] Use typing.Union instead of | to please Python 3.8/3.9 --- graphix/random_objects.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/graphix/random_objects.py b/graphix/random_objects.py index b811437e7..ee66f92ab 100644 --- a/graphix/random_objects.py +++ b/graphix/random_objects.py @@ -1,4 +1,4 @@ -from typing import Optional +from typing import Optional, Union import numpy as np import scipy.linalg @@ -30,7 +30,7 @@ def rand_unit(l: int): UNITS = np.array([1, 1j]) -def rand_dm(dim: int, rank: Optional[int] = None, dm_dtype=True) -> DensityMatrix | np.ndarray: +def rand_dm(dim: int, rank: Optional[int] = None, dm_dtype=True) -> Union[DensityMatrix, np.ndarray]: """Utility to generate random density matrices (positive semi-definite matrices with unit trace). Returns either a :class:`graphix.sim.density_matrix.DensityMatrix` or a :class:`np.ndarray` depending on the parameter `dm_dtype`. From 714dc123b96b38876687e3e339d4026bad797a0c Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Mon, 6 May 2024 11:40:54 +0200 Subject: [PATCH 20/48] add docs + cleanup --- graphix/sim/density_matrix.py | 40 +++++++++++++++++++++++++---------- graphix/sim/statevec.py | 37 ++++++++++++++++++-------------- graphix/states.py | 13 ------------ 3 files changed, 50 insertions(+), 40 deletions(-) diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index 5f8b3fdf7..7541354ae 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -37,14 +37,32 @@ def __init__( data: typing.Optional[Data] = graphix.states.BasicStates.PLUS, nqubit: typing.Optional[graphix.types.PositiveOrNullInt] = None, ): - """ - rewrite! - Parameters - ---------- - data : DensityMatrix, list, tuple, np.ndarray or None - Density matrix of shape (2**nqubits, 2**nqubits). - nqubit : int - Number of qubits. Default is 1. If both `data` and `nqubit` are specified, consistency is checked. + """Initialize density matrix objects. The behaviour builds on theo ne of `graphix.statevec.Statevec`. + `data` can be: + - a single :class:`graphix.states.State` (classical description of a quantum state) + - an iterable of :class:`graphix.states.State` objects + - an iterable of iterable of scalars (A 2**n x 2**n numerical density matrix) + - a `graphix.statevec.DensityMatrix` object + - a `graphix.statevec.Statevector` object + + If `nqubit` is not provided, the number of qubit is inferred from `data` and checked for consistency. + If only one :class:`graphix.states.State` is provided and nqubit is a valid integer, initialize the statevector + in the tensor product state. + If both `nqubit` and `data` are provided, consistency of the dimensions is checked. + If a `graphix.statevec.Statevec` or `graphix.statevec.DensityMatrix` is passed, returns a copy. + + + :param data: input data to prepare the state. Can be a classical description or a numerical input, defaults to graphix.states.BasicStates.PLUS + :type data: typing.Union[ + graphix.states.State, + "DensityMatrix", + Statevec, + typing.Iterable[graphix.states.State], + typing.Iterable[numbers.Number], + typing.Iterable[typing.Iterable[numbers.Number]], + ], optional + :param nqubit: number of qubits to prepare, defaults to None + :type nqubit: int, optional """ assert nqubit is None or isinstance(nqubit, numbers.Integral) and nqubit >= 0 @@ -63,8 +81,7 @@ def check_size_consistency(mat): if isinstance(data, typing.Iterable): input_list = list(data) if len(input_list) != 0: - # do try except else? - # needed since Object are iterable but not subscribable! + # needed since Object is iterable but not subscribable! try: if isinstance(input_list[0], typing.Iterable) and isinstance(input_list[0][0], numbers.Number): self.rho = np.array(input_list) @@ -77,7 +94,7 @@ def check_size_consistency(mat): except TypeError: pass statevec = Statevec(data, nqubit) - # NOTE sthis works since np.outer flattens the inputs! + # NOTE this works since np.outer flattens the inputs! self.rho = np.outer(statevec.psi, statevec.psi.conj()) self.Nqubit = len(statevec.dims()) @@ -326,6 +343,7 @@ def __init__(self, pattern, max_qubit_num=12, pr_calc=True, input_state: Data = pr_calc : bool whether or not to compute the probability distribution before choosing the measurement result. if False, measurements yield results 0/1 with 50% probabilities each. + input_state: same syntax as `graphix.statevec.DensityMatrix` constructor. """ # check that pattern has output nodes configured # assert len(pattern.output_nodes) > 0 diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index dbcfa5c53..a0fa00f50 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -38,8 +38,7 @@ def __init__( ----------- pattern : :class:`graphix.pattern.Pattern` object MBQC pattern to be simulated. - backend : str, 'statevector' - optional argument for simulation. + input_state: same syntax as `graphix.statevec.Statevec` constructor. max_qubit_num : int optional argument specifying the maximum number of qubits to be stored in the statevector at a time. @@ -217,22 +216,28 @@ def __init__( data: typing.Optional[SV_Data] = graphix.states.BasicStates.PLUS, nqubit: typing.Optional[graphix.types.PositiveOrNullInt] = None, ): - """Initialize statevector - Parameters - ---------- - data : is either - - a single state (:class:`graphix.states.State` object). THen prepares all nodes in that state (tensor product) - - a dictionary mapping the inputs to a :class:`graphix.states.State` object - - an arbitrary :class:`graphix.statevec.Statevec` object (arbitrary input) # TODO work on that since just copy? - nqubit : int, optional: ignored if iterable passed (State, direct data) - number of qubits. Defaults to 1. - # plus_states : bool, optional - whether or not to start all qubits in + state or 0 state. Defaults to + - - Defaults to |+> states and 1 qubit. - If nqubit > 1 and only one state : tensor all of them. Use the tensor method instead of hard code. + """Initialize statevector objects. The behaviour is as follows. `data` can be: + - a single :class:`graphix.states.State` (classical description of a quantum state) + - an iterable of :class:`graphix.states.State` objects + - an iterable of scalars (A 2**n numerical statevector) + - a `graphix.statevec.Statevec` object + + If `nqubit` is not provided, the number of qubit is inferred from `data` and checked for consistency. + If only one :class:`graphix.states.State` is provided and nqubit is a valid integer, initialize the statevector + in the tensor product state. + If both `nqubit` and `data` are provided, consistency of the dimensions is checked. + If a `graphix.statevec.Statevec` is passed, returns a copy. + + + :param data: input data to prepare the state. Can be a classical description or a numerical input, defaults to graphix.states.BasicStates.PLUS + :type data: typing.Union[ + graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[numbers.Number] + ], optional + :param nqubit: number of qubits to prepare, defaults to None + :type nqubit: int, optional """ + assert nqubit is None or isinstance(nqubit, numbers.Integral) and nqubit >= 0 if isinstance(data, Statevec): diff --git a/graphix/states.py b/graphix/states.py index c64a48510..446853d46 100644 --- a/graphix/states.py +++ b/graphix/states.py @@ -75,16 +75,3 @@ class BasicStates: # remove that in the end # need in TN backend VEC = [PLUS, MINUS, ZERO, ONE, PLUS_I, MINUS_I] - - -# Plane.cos.value Plane.cos is an Axis, Axis.value = 0,1,2 (enum) - -# Everywhere this is called. use StateVec(State)) -# class States: -# plus = np.array([1.0 / np.sqrt(2), 1.0 / np.sqrt(2)]) # plus -# minus = np.array([1.0 / np.sqrt(2), -1.0 / np.sqrt(2)]) # minus -# zero = np.array([1.0, 0.0]) # zero -# one = np.array([0.0, 1.0]) # one -# iplus = np.array([1.0 / np.sqrt(2), 1.0j / np.sqrt(2)]) # +1 eigenstate of Pauli Y -# iminus = np.array([1.0 / np.sqrt(2), -1.0j / np.sqrt(2)]) # -1 eigenstate of Pauli Y -# vec = [plus, minus, zero, one, iplus, iminus] From 4aafe65e663c41f3cfcdddafbce80b088f8f7b2f Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Mon, 6 May 2024 13:16:11 +0200 Subject: [PATCH 21/48] Update statevec.py --- graphix/sim/statevec.py | 1 - 1 file changed, 1 deletion(-) diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index a0fa00f50..fd84894a0 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -216,7 +216,6 @@ def __init__( data: typing.Optional[SV_Data] = graphix.states.BasicStates.PLUS, nqubit: typing.Optional[graphix.types.PositiveOrNullInt] = None, ): - """Initialize statevector objects. The behaviour is as follows. `data` can be: - a single :class:`graphix.states.State` (classical description of a quantum state) - an iterable of :class:`graphix.states.State` objects From 9d9302713a117f6f4e181281b5d00e50791903f5 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Tue, 7 May 2024 16:56:56 +0200 Subject: [PATCH 22/48] integrated review comments --- .vscode/settings.json | 16 ---------------- graphix/ops.py | 2 -- graphix/pattern.py | 7 ------- graphix/random_objects.py | 3 +-- graphix/sim/density_matrix.py | 2 -- graphix/sim/graphix.code-workspace | 7 ------- graphix/sim/statevec.py | 4 ---- graphix/states.py | 3 --- tests/test_density_matrix.py | 3 --- tests/test_statevec_backend.py | 4 +--- 10 files changed, 2 insertions(+), 49 deletions(-) delete mode 100644 .vscode/settings.json delete mode 100644 graphix/sim/graphix.code-workspace diff --git a/.vscode/settings.json b/.vscode/settings.json deleted file mode 100644 index 7a4707700..000000000 --- a/.vscode/settings.json +++ /dev/null @@ -1,16 +0,0 @@ -{ - "workbench.colorCustomizations": { - "activityBar.activeBackground": "#287d9e", - "activityBar.background": "#287d9e", - "activityBar.foreground": "#e7e7e7", - "activityBar.inactiveForeground": "#e7e7e799", - "activityBarBadge.background": "#e599d0", - "activityBarBadge.foreground": "#15202b" - }, - "peacock.color": "#1e5d75", - "python.testing.pytestArgs": [ - "tests" - ], - "python.testing.unittestEnabled": false, - "python.testing.pytestEnabled": true -} \ No newline at end of file diff --git a/graphix/ops.py b/graphix/ops.py index b74366333..9404ebe10 100644 --- a/graphix/ops.py +++ b/graphix/ops.py @@ -11,8 +11,6 @@ class Ops: """Basic single- and two-qubits operators""" - # class attributes. DOn' need to instantiate. Inherited by all class members - x = np.array([[0, 1], [1, 0]]) y = np.array([[0, -1j], [1j, 0]]) z = np.array([[1, 0], [0, -1]]) diff --git a/graphix/pattern.py b/graphix/pattern.py index fcbd80c5b..31706f21e 100644 --- a/graphix/pattern.py +++ b/graphix/pattern.py @@ -814,10 +814,6 @@ def connected_edges(self, node, edges): connected: set of tuple set of connected edges """ - # TODO modify that by using the graph nx.graph.edges(node)? and cached get_graph()? - # like in def get_measurement_order_from_flow(self): with self.get_graph() - # FIXME - # BUG connected = set() for edge in edges: @@ -1017,9 +1013,6 @@ def get_max_degree(self): max_degree = max([i for i in dict(degree).values()]) return max_degree - # TODO functools.cache() It is called in get measurement order from (g)flow - # - # It is called in get measurement order from (g)flow def get_graph(self): """returns the list of nodes and edges from the command sequence, extracted from 'N' and 'E' commands. diff --git a/graphix/random_objects.py b/graphix/random_objects.py index ee66f92ab..f99077323 100644 --- a/graphix/random_objects.py +++ b/graphix/random_objects.py @@ -49,9 +49,8 @@ def rand_dm(dim: int, rank: Optional[int] = None, dm_dtype=True) -> Union[Densit Note that setting `dm_dtype=False` allows to generate "density matrices" inconsistent with qubits i.e. with dimensions not being powers of 2. """ - # if not provided, use a random value. if rank is None: - rank = dim # np.random.randint(1, dim + 1) + rank = dim evals = np.random.rand(rank) diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index 7541354ae..8e5500533 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -194,8 +194,6 @@ def expectation_single(self, op, i): return np.trace(st1.rho) - # TODO - # @property def dims(self): return self.rho.shape diff --git a/graphix/sim/graphix.code-workspace b/graphix/sim/graphix.code-workspace deleted file mode 100644 index 9e68e72b8..000000000 --- a/graphix/sim/graphix.code-workspace +++ /dev/null @@ -1,7 +0,0 @@ -{ - "folders": [ - { - "path": ".." - } - ] -} \ No newline at end of file diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index fd84894a0..2fc9154b4 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -261,7 +261,6 @@ def __init__( if nqubit is not None and nqubit != 0: raise ValueError("nqubit is not null but input state is empty.") - # warnings.warn(f"Called Statevec with 0 qubits. Ignoring the state.") self.psi = np.array(1, dtype=np.complex128) # self.Nqubit = 0 else: @@ -287,14 +286,11 @@ def __init__( psi = np.array(input_list) if not np.allclose(np.sqrt(np.sum(np.abs(psi) ** 2)), 1): raise ValueError("Input state is not normalized") - # just reshape - # NOTE too many conversions to numpy arrays? self.psi = psi.reshape((2,) * nqubit) else: raise TypeError( f"First element of data has type {type(input_list[0])} whereas Number or State is expected" ) - # self.Nqubit = state.Nqubit def __repr__(self): return f"Statevec object with statevector {self.psi} and length {self.dims()}." diff --git a/graphix/states.py b/graphix/states.py index 446853d46..e5ec9299b 100644 --- a/graphix/states.py +++ b/graphix/states.py @@ -27,9 +27,6 @@ def get_densitymatrix(self) -> np.ndarray: return np.outer(self.get_statevector(), self.get_statevector().conj()) -# don't turn it into Statevec here -# Weird not to allow all states? -# Made it inherit from more generic State class. class PlanarState(pydantic.BaseModel, State): """Light object used to instantiate backends. doesn't cover all possible states but this is diff --git a/tests/test_density_matrix.py b/tests/test_density_matrix.py index 1014affcf..a44493b45 100644 --- a/tests/test_density_matrix.py +++ b/tests/test_density_matrix.py @@ -181,9 +181,6 @@ def test_init_with_densitymatrix_sucess(self): assert np.allclose(dm2.rho, expected_dm) assert np.allclose(dm2.rho, dm.rho) - # assert dm.Nqubit == n - # assert dm.rho.shape == (2**n, 2**n) - def test_evolve_single_fail(self): dm = DensityMatrix(nqubit=2) # generate random 4 x 4 unitary matrix diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index ef9828487..3a94b1d15 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -28,7 +28,6 @@ def test_remove_one_qubit(self): np.testing.assert_almost_equal(np.abs(sv.psi.flatten().dot(sv2.psi.flatten().conj())), 1) - # TODO This is a weird test! def test_measurement_into_each_XYZ_basis(self): n = 3 k = 0 @@ -36,9 +35,8 @@ def test_measurement_into_each_XYZ_basis(self): # NOTE weird choice (MINUS is orthogonal to PLUS so zero) for state in [BasicStates.PLUS, BasicStates.ZERO, BasicStates.ONE, BasicStates.PLUS_I, BasicStates.MINUS_I]: m_op = np.outer(state.get_statevector(), state.get_statevector().T.conjugate()) - # print(m_op) + sv = Statevec(nqubit=n) - # print(sv) sv.evolve(m_op, [k]) sv.remove_qubit(k) From ad2dc71881a72dd3b605857e744ce66e8fa2f8c7 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Tue, 7 May 2024 17:02:39 +0200 Subject: [PATCH 23/48] additional comment handling --- .pre-commit-config.yaml | 2 +- graphix/sim/statevec.py | 7 ------- 2 files changed, 1 insertion(+), 8 deletions(-) diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index ff5c14609..1b367f6e7 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -1,6 +1,6 @@ repos: - repo: https://github.com/psf/black - rev: 22.8.0 + rev: 24.4.0 hooks: - id: black args: [--line-length=120] diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 2fc9154b4..bdcb68a30 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -14,12 +14,6 @@ from graphix.clifford import CLIFFORD, CLIFFORD_CONJ, CLIFFORD_MUL from graphix.ops import Ops -# Python >= 3.9 -# from collections.abc import Iterable # or use Protocols? -# https://stackoverflow.com/questions/49427944/typehints-for-sized-iterable-in-python -# Python >= 3.8 -# typing.Iterable[T] - class StatevectorBackend(graphix.sim.base_backend.Backend): """MBQC simulator with statevector method.""" @@ -210,7 +204,6 @@ def meas_op(angle, vop=0, plane="XY", choice=0): class Statevec: """Statevector object""" - # TODO at this stage no need for indices just be careful of the ordering in add_nodes def __init__( self, data: typing.Optional[SV_Data] = graphix.states.BasicStates.PLUS, From 6181ef39114cfeb8aee88cba1ef0d0158c0d4016 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Tue, 7 May 2024 17:22:35 +0200 Subject: [PATCH 24/48] Delete .pre-commit-config.yaml --- .pre-commit-config.yaml | 7 ------- 1 file changed, 7 deletions(-) delete mode 100644 .pre-commit-config.yaml diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml deleted file mode 100644 index 1b367f6e7..000000000 --- a/.pre-commit-config.yaml +++ /dev/null @@ -1,7 +0,0 @@ -repos: -- repo: https://github.com/psf/black - rev: 24.4.0 - hooks: - - id: black - args: [--line-length=120] - files: ^(graphix|test)/ From f3102e4589ea8fa443209de0cccf36c27ac0d9be Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Mon, 13 May 2024 17:22:11 +0200 Subject: [PATCH 25/48] isort --- examples/noisy_mbqc.py | 3 ++- tests/conftest.py | 3 ++- tests/test_density_matrix.py | 2 +- 3 files changed, 5 insertions(+), 3 deletions(-) diff --git a/examples/noisy_mbqc.py b/examples/noisy_mbqc.py index e9108834c..835bf2680 100644 --- a/examples/noisy_mbqc.py +++ b/examples/noisy_mbqc.py @@ -9,9 +9,10 @@ First, let us import relevant modules and define a pattern """ +import matplotlib.pyplot as plt + # %% import numpy as np -import matplotlib.pyplot as plt from graphix import Circuit from graphix.channels import KrausChannel, dephasing_channel diff --git a/tests/conftest.py b/tests/conftest.py index d3d13b0c5..8b51963e3 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -1,7 +1,8 @@ import pytest from numpy.random import PCG64, Generator -import tests.random_circuit + import graphix.transpiler +import tests.random_circuit SEED = 25 diff --git a/tests/test_density_matrix.py b/tests/test_density_matrix.py index 257065034..ded866db8 100644 --- a/tests/test_density_matrix.py +++ b/tests/test_density_matrix.py @@ -1,6 +1,6 @@ from __future__ import annotations -import functools +import functools import random from copy import deepcopy from typing import TYPE_CHECKING From 36615d2df846c7c385bdfd76a1b42b464b10e47b Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Tue, 14 May 2024 14:46:47 +0200 Subject: [PATCH 26/48] updates --- graphix/sim/density_matrix.py | 2 -- graphix/sim/statevec.py | 2 +- graphix/simulator.py | 5 ----- tests/test_pattern.py | 3 --- tests/test_statevec.py | 1 - tests/test_statevec_backend.py | 10 +--------- 6 files changed, 2 insertions(+), 21 deletions(-) diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index 8e5500533..da8b6c826 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -343,8 +343,6 @@ def __init__(self, pattern, max_qubit_num=12, pr_calc=True, input_state: Data = if False, measurements yield results 0/1 with 50% probabilities each. input_state: same syntax as `graphix.statevec.DensityMatrix` constructor. """ - # check that pattern has output nodes configured - # assert len(pattern.output_nodes) > 0 self.pattern = pattern self.results = deepcopy(pattern.results) self.state = None diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index bdcb68a30..d195ea525 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -255,7 +255,7 @@ def __init__( raise ValueError("nqubit is not null but input state is empty.") self.psi = np.array(1, dtype=np.complex128) - # self.Nqubit = 0 + else: if isinstance(input_list[0], graphix.states.State): graphix.types.check_list_elements(input_list, graphix.states.State) diff --git a/graphix/simulator.py b/graphix/simulator.py index d606d9a60..65d3fdbdd 100644 --- a/graphix/simulator.py +++ b/graphix/simulator.py @@ -50,8 +50,6 @@ def __init__(self, pattern, backend="statevector", noise_model=None, **kwargs): ) if noise_model is not None: self.set_noise_model(noise_model) - # if noise: have to compute the probabilities - # NOTE : could remove, pr_calc defaults to True now. self.backend = DensityMatrixBackend(pattern, pr_calc=True, **kwargs) elif backend in {"tensornetwork", "mps"} and noise_model is None: self.noise_model = None @@ -84,9 +82,6 @@ def run(self): the output quantum state, in the representation depending on the backend used. """ - # use add_nodes or write a new method? - # self.backend.initialize_inputs(self.pattern.input_nodes, option, ...) - # self.backend.add_nodes(self.pattern.input_nodes, input_state=state) if self.noise_model is None: for cmd in self.pattern: if cmd[0] == "N": diff --git a/tests/test_pattern.py b/tests/test_pattern.py index be44581a2..753ce2f92 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -520,10 +520,7 @@ def test_sv_sim(self, fx_rng: Generator, nqb, rand_circ): out = randpattern.simulate_pattern(backend="statevector", input_state=states) out_circ = rand_circ.simulate_statevector(input_state=states).statevec - - # MBQC is up to a global phase! np.testing.assert_almost_equal(np.abs(np.dot(out.psi.flatten().conjugate(), out_circ.psi.flatten())), 1) - # assert np.allclose(out.psi, out_circ.psi) def test_dm_sim(self): pass diff --git a/tests/test_statevec.py b/tests/test_statevec.py index ce01f89df..a366897ca 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -17,7 +17,6 @@ def setUp(self): # set up the random numbers self.rng = np.random.default_rng() # seed=422 - # Errors: types, size, # test injitializing one qubit in plus state def test_default_success(self): diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index cf4728c90..fd466ad5f 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -31,11 +31,10 @@ def test_remove_one_qubit(self) -> None: sv2.normalize() assert np.abs(sv.psi.flatten().dot(sv2.psi.flatten().conj())) == pytest.approx(1) - @pytest.mark.parametrize( "state", [BasicStates.PLUS, BasicStates.ZERO, BasicStates.ONE, BasicStates.PLUS_I, BasicStates.MINUS_I] ) - def test_measurement_into_each_xyz_basis(self, state: BasicStates) -> None: + def test_measurement_into_each_XYZ_basis(self, state: BasicStates) -> None: n = 3 k = 0 # for measurement into |-> returns [[0, 0], ..., [0, 0]] (whose norm is zero) @@ -67,14 +66,12 @@ def test_init_success(self, hadamardpattern, fx_rng: Generator): backend = StatevectorBackend(hadamardpattern) vec = Statevec(nqubit=1) assert np.allclose(vec.psi, backend.state.psi) - # assert backend.state.Nqubit == 1 assert len(backend.state.dims()) == 1 # minus state backend = StatevectorBackend(hadamardpattern, input_state=BasicStates.MINUS) vec = Statevec(nqubit=1, data=BasicStates.MINUS) assert np.allclose(vec.psi, backend.state.psi) - # assert backend.state.Nqubit == 1 assert len(backend.state.dims()) == 1 # random planar state @@ -90,11 +87,6 @@ def test_init_success(self, hadamardpattern, fx_rng: Generator): # data input and Statevec input def test_init_fail(self, hadamardpattern, fx_rng: Generator): - # incorrect number of dimensions for State input - # only one input node, two states provided - # doesn't fail! just takes the first qubit! - # Discard second qubit so can be whatever - rand_angle = fx_rng.random(2) * 2 * np.pi rand_plane = fx_rng.choice(np.array([i for i in graphix.pauli.Plane]), 2) From 0a4fa528736a5ef8c4fcd3c926a9724db395508c Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Tue, 14 May 2024 14:50:38 +0200 Subject: [PATCH 27/48] black --- tests/test_statevec.py | 1 - 1 file changed, 1 deletion(-) diff --git a/tests/test_statevec.py b/tests/test_statevec.py index a366897ca..5a969afb0 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -17,7 +17,6 @@ def setUp(self): # set up the random numbers self.rng = np.random.default_rng() # seed=422 - # test injitializing one qubit in plus state def test_default_success(self): vec = Statevec(nqubit=1) From ff9769e01da52c3862bba08783717a95e00dfa03 Mon Sep 17 00:00:00 2001 From: mgarnier59 <50111289+mgarnier59@users.noreply.github.com> Date: Tue, 14 May 2024 14:53:35 +0200 Subject: [PATCH 28/48] black --- tests/test_statevec_backend.py | 1 + 1 file changed, 1 insertion(+) diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index fd466ad5f..8e0fc3254 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -31,6 +31,7 @@ def test_remove_one_qubit(self) -> None: sv2.normalize() assert np.abs(sv.psi.flatten().dot(sv2.psi.flatten().conj())) == pytest.approx(1) + @pytest.mark.parametrize( "state", [BasicStates.PLUS, BasicStates.ZERO, BasicStates.ONE, BasicStates.PLUS_I, BasicStates.MINUS_I] ) From d4d3de832604b69839daedfd7b8e300b1fb597fb Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Tue, 21 May 2024 14:53:24 +0200 Subject: [PATCH 29/48] Ruff-suggested fixes --- graphix/sim/density_matrix.py | 5 ++--- graphix/sim/statevec.py | 4 +--- 2 files changed, 3 insertions(+), 6 deletions(-) diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index da8b6c826..4ccd78f10 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -3,7 +3,6 @@ Simulate MBQC with density matrix representation. """ -import functools import numbers import typing from copy import deepcopy @@ -15,9 +14,9 @@ import graphix.types from graphix.channels import KrausChannel from graphix.clifford import CLIFFORD -from graphix.linalg_validations import check_hermitian, check_psd, check_square, check_unit_trace +from graphix.linalg_validations import check_psd, check_square, check_unit_trace from graphix.ops import Ops -from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, Statevec, meas_op +from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, Statevec Data = typing.Union[ graphix.states.State, diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index d195ea525..decb486aa 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -1,17 +1,15 @@ import functools import numbers import typing -import warnings from copy import deepcopy import numpy as np -import pydantic import graphix.pauli import graphix.sim.base_backend import graphix.states import graphix.types -from graphix.clifford import CLIFFORD, CLIFFORD_CONJ, CLIFFORD_MUL +from graphix.clifford import CLIFFORD, CLIFFORD_CONJ from graphix.ops import Ops From f377001769c2a3ff6ace2a84fdbe7826a6287889 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Tue, 21 May 2024 14:53:33 +0200 Subject: [PATCH 30/48] Check |+> input states for Pauli-preprocessed patterns --- graphix/sim/density_matrix.py | 2 ++ graphix/sim/statevec.py | 2 ++ 2 files changed, 4 insertions(+) diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index 4ccd78f10..2af3fcd46 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -343,6 +343,8 @@ def __init__(self, pattern, max_qubit_num=12, pr_calc=True, input_state: Data = input_state: same syntax as `graphix.statevec.DensityMatrix` constructor. """ self.pattern = pattern + if pattern._pauli_preprocessed and not input_state != graphix.states.BasicStates.PLUS: + raise ValueError("Pauli preprocessing is currently only available when inputs are initialized in |+> state") self.results = deepcopy(pattern.results) self.state = None self.node_index = [] diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index decb486aa..523e1be2c 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -41,6 +41,8 @@ def __init__( # check that pattern has output nodes configured # assert len(pattern.output_nodes) > 0 self.pattern = pattern + if pattern._pauli_preprocessed and input_state != graphix.states.BasicStates.PLUS: + raise ValueError("Pauli preprocessing is currently only available when inputs are initialized in |+> state") self.results = deepcopy(pattern.results) self.state = None self.node_index = [] From 8836a61e8374c7966bff31f9b3021cf3c8763a8d Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Tue, 21 May 2024 15:45:08 +0200 Subject: [PATCH 31/48] Fix check for density matrix Reported by Maxime Garnier. https://github.com/TeamGraphix/graphix/commit/f377001769c2a3ff6ace2a84fdbe7826a6287889#r142237724 The test for Pauli preprocessing in density_matrix.py should be ``` pattern._pauli_preprocessed and input_state != graphix.states.BasicStates.PLUS ``` without `not`. This bug wasn't caught by `pytest` because the density matrix backend was never tested with Pauli-preprocessed patterns. I think a natural place to add such a test is in `test_pattern.py`, by "parametrizing" `test_pauli_measurment` with `backend`. This was not enough because `test_pauli_measurment` was defined twice in the same file and the second definition was hiding the first one. Therefore, this commit: - rename the first `test_pauli_measurment` into `test_pauli_measurement_random_circuit` (fixing the typo by the way); - "parametrize" `test_pauli_measurement_random_circuit` with `backend`, testing `"statevector"` and `"densitymatrix"`; - a TODO is added for tensor network backend, since the default `graph_prep` mode (`auto`) selects the mode `parallel` despite this mode is incompatible with non-standardize patterns. I propose to leave that for another PR. --- graphix/sim/density_matrix.py | 2 +- tests/test_pattern.py | 24 +++++++++++++++++------- 2 files changed, 18 insertions(+), 8 deletions(-) diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index 2af3fcd46..be6e40345 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -343,7 +343,7 @@ def __init__(self, pattern, max_qubit_num=12, pr_calc=True, input_state: Data = input_state: same syntax as `graphix.statevec.DensityMatrix` constructor. """ self.pattern = pattern - if pattern._pauli_preprocessed and not input_state != graphix.states.BasicStates.PLUS: + if pattern._pauli_preprocessed and input_state != graphix.states.BasicStates.PLUS: raise ValueError("Pauli preprocessing is currently only available when inputs are initialized in |+> state") self.results = deepcopy(pattern.results) self.state = None diff --git a/tests/test_pattern.py b/tests/test_pattern.py index 753ce2f92..9744eaf91 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -12,6 +12,7 @@ from graphix.pattern import CommandNode, Pattern from graphix.simulator import PatternSimulator from graphix.transpiler import Circuit +import graphix.sim.base_backend if TYPE_CHECKING: from collections.abc import Sequence @@ -122,7 +123,11 @@ def test_shift_signals(self, fx_bg: PCG64, jumps: int) -> None: assert np.abs(np.dot(state_mbqc.flatten().conjugate(), state.flatten())) == pytest.approx(1) @pytest.mark.parametrize("jumps", range(1, 11)) - def test_pauli_measurment(self, fx_bg: PCG64, jumps: int, use_rustworkx: bool = True) -> None: + @pytest.mark.parametrize("backend", ["statevector", "densitymatrix"]) + # TODO: tensor network backend is excluded because "parallel preparation strategy does not support not-standardized pattern". + def test_pauli_measurement_random_circuit( + self, fx_bg: PCG64, jumps: int, backend: graphix.sim.base_backend.Backend, use_rustworkx: bool = True + ) -> None: rng = Generator(fx_bg.jumped(jumps)) nqubits = 3 depth = 3 @@ -133,11 +138,16 @@ def test_pauli_measurment(self, fx_bg: PCG64, jumps: int, use_rustworkx: bool = pattern.perform_pauli_measurements(use_rustworkx=use_rustworkx) pattern.minimize_space() state = circuit.simulate_statevector().statevec - state_mbqc = pattern.simulate_pattern() - assert np.abs(np.dot(state_mbqc.flatten().conjugate(), state.flatten())) == pytest.approx(1) + state_mbqc = pattern.simulate_pattern(backend) + if backend == "statevector": + assert np.abs(np.dot(state_mbqc.flatten().conjugate(), state.flatten())) == pytest.approx(1) + elif backend == "densitymatrix": + assert np.allclose(state_mbqc.rho, graphix.sim.density_matrix.DensityMatrix(state.flatten()).rho) @pytest.mark.parametrize("jumps", range(1, 11)) - def test_pauli_measurment_leave_input(self, fx_bg: PCG64, jumps: int, use_rustworkx: bool = True) -> None: + def test_pauli_measurement_leave_input_random_circuit( + self, fx_bg: PCG64, jumps: int, use_rustworkx: bool = True + ) -> None: rng = Generator(fx_bg.jumped(jumps)) nqubits = 3 depth = 3 @@ -152,7 +162,7 @@ def test_pauli_measurment_leave_input(self, fx_bg: PCG64, jumps: int, use_rustwo assert np.abs(np.dot(state_mbqc.flatten().conjugate(), state.flatten())) == pytest.approx(1) @pytest.mark.parametrize("jumps", range(1, 11)) - def test_pauli_measurment_opt_gate(self, fx_bg: PCG64, jumps: int, use_rustworkx: bool = True) -> None: + def test_pauli_measurement_opt_gate(self, fx_bg: PCG64, jumps: int, use_rustworkx: bool = True) -> None: rng = Generator(fx_bg.jumped(jumps)) nqubits = 3 depth = 3 @@ -167,7 +177,7 @@ def test_pauli_measurment_opt_gate(self, fx_bg: PCG64, jumps: int, use_rustworkx assert np.abs(np.dot(state_mbqc.flatten().conjugate(), state.flatten())) == pytest.approx(1) @pytest.mark.parametrize("jumps", range(1, 11)) - def test_pauli_measurment_opt_gate_transpiler(self, fx_bg: PCG64, jumps: int, use_rustworkx: bool = True) -> None: + def test_pauli_measurement_opt_gate_transpiler(self, fx_bg: PCG64, jumps: int, use_rustworkx: bool = True) -> None: rng = Generator(fx_bg.jumped(jumps)) nqubits = 3 depth = 3 @@ -182,7 +192,7 @@ def test_pauli_measurment_opt_gate_transpiler(self, fx_bg: PCG64, jumps: int, us assert np.abs(np.dot(state_mbqc.flatten().conjugate(), state.flatten())) == pytest.approx(1) @pytest.mark.parametrize("jumps", range(1, 11)) - def test_pauli_measurment_opt_gate_transpiler_without_signalshift( + def test_pauli_measurement_opt_gate_transpiler_without_signalshift( self, fx_bg: PCG64, jumps: int, From 6ebfcbe7afadbb65d62aaad8c863999c87a39899 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Tue, 21 May 2024 16:10:46 +0200 Subject: [PATCH 32/48] isort --- tests/test_pattern.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_pattern.py b/tests/test_pattern.py index 9744eaf91..e831acf88 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -7,12 +7,12 @@ import pytest import graphix.ops +import graphix.sim.base_backend import graphix.states import tests.random_circuit as rc from graphix.pattern import CommandNode, Pattern from graphix.simulator import PatternSimulator from graphix.transpiler import Circuit -import graphix.sim.base_backend if TYPE_CHECKING: from collections.abc import Sequence From 8bf21434d7a294c0c3a93ab8851e268e42066bbb Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Wed, 29 May 2024 17:22:48 +0200 Subject: [PATCH 33/48] WIP --- graphix/random_objects.py | 18 +++++++-------- graphix/sim/density_matrix.py | 42 +++++++++++++++++------------------ graphix/sim/statevec.py | 30 +++++++++++++------------ graphix/transpiler.py | 4 ++-- 4 files changed, 47 insertions(+), 47 deletions(-) diff --git a/graphix/random_objects.py b/graphix/random_objects.py index f99077323..e88353f4e 100644 --- a/graphix/random_objects.py +++ b/graphix/random_objects.py @@ -1,4 +1,4 @@ -from typing import Optional, Union +from __future__ import annotations import numpy as np import scipy.linalg @@ -9,28 +9,28 @@ from graphix.sim.density_matrix import DensityMatrix -def rand_herm(l: int): +def rand_herm(sz: int): """ - generate random hermitian matrix of size l*l + generate random hermitian matrix of size sz*sz """ - tmp = np.random.rand(l, l) + 1j * np.random.rand(l, l) + tmp = np.random.rand(sz, sz) + 1j * np.random.rand(sz, sz) return tmp + tmp.conj().T -def rand_unit(l: int): +def rand_unit(sz: int): """ - generate haar random unitary matrix of size l*l + generate haar random unitary matrix of size sz*sz """ - if l == 1: + if sz == 1: return np.array([np.exp(1j * np.random.rand(1) * 2 * np.pi)]) else: - return unitary_group.rvs(l) + return unitary_group.rvs(sz) UNITS = np.array([1, 1j]) -def rand_dm(dim: int, rank: Optional[int] = None, dm_dtype=True) -> Union[DensityMatrix, np.ndarray]: +def rand_dm(dim: int, rank: int | None = None, dm_dtype=True) -> DensityMatrix | np.ndarray: """Utility to generate random density matrices (positive semi-definite matrices with unit trace). Returns either a :class:`graphix.sim.density_matrix.DensityMatrix` or a :class:`np.ndarray` depending on the parameter `dm_dtype`. diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index be6e40345..f4b7719c4 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -3,8 +3,10 @@ Simulate MBQC with density matrix representation. """ +from __future__ import annotations + +import collections import numbers -import typing from copy import deepcopy import numpy as np @@ -18,23 +20,14 @@ from graphix.ops import Ops from graphix.sim.statevec import CNOT_TENSOR, CZ_TENSOR, SWAP_TENSOR, Statevec -Data = typing.Union[ - graphix.states.State, - "DensityMatrix", - Statevec, - typing.Iterable[graphix.states.State], - typing.Iterable[numbers.Number], - typing.Iterable[typing.Iterable[numbers.Number]], -] - class DensityMatrix: """DensityMatrix object.""" def __init__( self, - data: typing.Optional[Data] = graphix.states.BasicStates.PLUS, - nqubit: typing.Optional[graphix.types.PositiveOrNullInt] = None, + data: Data = graphix.states.BasicStates.PLUS, + nqubit: graphix.types.PositiveOrNullInt | None = None, ): """Initialize density matrix objects. The behaviour builds on theo ne of `graphix.statevec.Statevec`. `data` can be: @@ -52,14 +45,7 @@ def __init__( :param data: input data to prepare the state. Can be a classical description or a numerical input, defaults to graphix.states.BasicStates.PLUS - :type data: typing.Union[ - graphix.states.State, - "DensityMatrix", - Statevec, - typing.Iterable[graphix.states.State], - typing.Iterable[numbers.Number], - typing.Iterable[typing.Iterable[numbers.Number]], - ], optional + :type data: graphix.states.State | "DensityMatrix" | Statevec | collections.abc.Iterable[graphix.states.State] |collections.abc.Iterable[numbers.Number] | collections.abc.Iterable[collections.abc.Iterable[numbers.Number]], optional :param nqubit: number of qubits to prepare, defaults to None :type nqubit: int, optional """ @@ -77,12 +63,14 @@ def check_size_consistency(mat): self.rho = data.rho.copy() self.Nqubit = data.Nqubit return - if isinstance(data, typing.Iterable): + if isinstance(data, collections.abc.Iterable): input_list = list(data) if len(input_list) != 0: # needed since Object is iterable but not subscribable! try: - if isinstance(input_list[0], typing.Iterable) and isinstance(input_list[0][0], numbers.Number): + if isinstance(input_list[0], collections.abc.Iterable) and isinstance( + input_list[0][0], numbers.Number + ): self.rho = np.array(input_list) assert check_square(self.rho) check_size_consistency(self.rho) @@ -448,3 +436,13 @@ def finalize(self): """To be run at the end of pattern simulation.""" self.sort_qubits() self.state.normalize() + + +Data = ( + graphix.states.State + | DensityMatrix + | Statevec + | collections.abc.Iterable[graphix.states.State] + | collections.abc.Iterable[numbers.Number] + | collections.abc.Iterable[collections.abc.Iterable[numbers.Number]] +) diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 523e1be2c..4bfd23469 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -1,6 +1,8 @@ +from __future__ import annotations + +import collections import functools import numbers -import typing from copy import deepcopy import numpy as np @@ -19,9 +21,7 @@ class StatevectorBackend(graphix.sim.base_backend.Backend): def __init__( self, pattern, - input_state: typing.Union[ - graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[numbers.Number] - ] = graphix.states.BasicStates.PLUS, + input_state: SV_Data = graphix.states.BasicStates.PLUS, max_qubit_num=20, pr_calc=True, ): @@ -196,18 +196,14 @@ def meas_op(angle, vop=0, plane="XY", choice=0): dtype=np.complex128, ) -SV_Data = typing.Union[ - graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[numbers.Number] -] - class Statevec: """Statevector object""" def __init__( self, - data: typing.Optional[SV_Data] = graphix.states.BasicStates.PLUS, - nqubit: typing.Optional[graphix.types.PositiveOrNullInt] = None, + data: SV_Data = graphix.states.BasicStates.PLUS, + nqubit: graphix.types.PositiveOrNullInt | None = None, ): """Initialize statevector objects. The behaviour is as follows. `data` can be: - a single :class:`graphix.states.State` (classical description of a quantum state) @@ -223,9 +219,7 @@ def __init__( :param data: input data to prepare the state. Can be a classical description or a numerical input, defaults to graphix.states.BasicStates.PLUS - :type data: typing.Union[ - graphix.states.State, "Statevec", typing.Iterable[graphix.states.State], typing.Iterable[numbers.Number] - ], optional + :type data: SV_Data, optional :param nqubit: number of qubits to prepare, defaults to None :type nqubit: int, optional """ @@ -245,7 +239,7 @@ def __init__( if nqubit is None: nqubit = 1 input_list = [data] * nqubit - elif isinstance(data, typing.Iterable): + elif isinstance(data, collections.abc.Iterable): input_list = list(data) else: raise TypeError(f"Incorrect type for data: {type(data)}") @@ -502,3 +496,11 @@ def expectation_value(self, op, qargs): def _get_statevec_norm(psi): """returns norm of the state""" return np.sqrt(np.sum(psi.flatten().conj() * psi.flatten())) + + +SV_Data = ( + graphix.states.State + | Statevec + | collections.abc.Iterable[graphix.states.State] + | collections.abc.Iterable[numbers.Number] +) diff --git a/graphix/transpiler.py b/graphix/transpiler.py index 9d7ae6244..1b668fd7e 100644 --- a/graphix/transpiler.py +++ b/graphix/transpiler.py @@ -8,7 +8,7 @@ import dataclasses from copy import deepcopy -from typing import Optional, Sequence +from typing import Sequence import numpy as np @@ -1358,7 +1358,7 @@ def _sort_outputs(self, pattern: Pattern, output_nodes: Sequence[int]): elif cmd[1] in old_out: cmd[1] = output_nodes[old_out.index(cmd[1])] - def simulate_statevector(self, input_state: Optional[SV_Data] = None) -> SimulateResult: + def simulate_statevector(self, input_state: SV_Data | None = None) -> SimulateResult: """Run statevector simulation of the gate sequence, using graphix.Statevec Parameters From 54c03380702b8cd3e2f1e811191ca9f321c34af4 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Thu, 30 May 2024 12:08:30 +0200 Subject: [PATCH 34/48] Fix None result and rng fixtures --- tests/test_pattern.py | 4 +- tests/test_statevec.py | 76 ++++++++++++++++------------------ tests/test_statevec_backend.py | 4 +- 3 files changed, 39 insertions(+), 45 deletions(-) diff --git a/tests/test_pattern.py b/tests/test_pattern.py index 6682b45fa..8ebb14d81 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -530,7 +530,7 @@ def test_pauli_measurement_end_with_measure(self) -> None: # for testing with arbitrary inputs # SV and DM backend class TestPatternSim: - def test_sv_sim(self, fx_rng: Generator, nqb, rand_circ): + def test_sv_sim(self, fx_rng: Generator, nqb, rand_circ) -> None: rand_angles = fx_rng.random(nqb) * 2 * np.pi rand_planes = fx_rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] @@ -540,7 +540,7 @@ def test_sv_sim(self, fx_rng: Generator, nqb, rand_circ): out_circ = rand_circ.simulate_statevector(input_state=states).statevec np.testing.assert_almost_equal(np.abs(np.dot(out.psi.flatten().conjugate(), out_circ.psi.flatten())), 1) - def test_dm_sim(self): + def test_dm_sim(self) -> None: pass diff --git a/tests/test_statevec.py b/tests/test_statevec.py index 5a969afb0..92575b098 100644 --- a/tests/test_statevec.py +++ b/tests/test_statevec.py @@ -1,29 +1,23 @@ import functools -import unittest import numpy as np import pytest import graphix.pauli -import graphix.random_objects as randobj from graphix.sim.statevec import Statevec from graphix.states import BasicStates, PlanarState -class TestStatevec(unittest.TestCase): +class TestStatevec: """Test for Statevec class. Particularly new constructor.""" - def setUp(self): - # set up the random numbers - self.rng = np.random.default_rng() # seed=422 - # test injitializing one qubit in plus state - def test_default_success(self): + def test_default_success(self) -> None: vec = Statevec(nqubit=1) assert np.allclose(vec.psi, np.array([1, 1] / np.sqrt(2))) assert len(vec.dims()) == 1 - def test_basicstates_success(self): + def test_basicstates_success(self) -> None: # minus vec = Statevec(nqubit=1, data=BasicStates.MINUS) assert np.allclose(vec.psi, np.array([1, -1] / np.sqrt(2))) @@ -50,8 +44,8 @@ def test_basicstates_success(self): assert len(vec.dims()) == 1 # even more tests? - def test_default_tensor_success(self): - nqb = self.rng.integers(2, 5) + def test_default_tensor_success(self, fx_rng: np.random.Generator) -> None: + nqb = fx_rng.integers(2, 5) print(f"nqb is {nqb}") vec = Statevec(nqubit=nqb) assert np.allclose(vec.psi, np.ones(((2,) * nqb)) / (np.sqrt(2)) ** nqb) @@ -64,8 +58,8 @@ def test_default_tensor_success(self): assert len(vec.dims()) == nqb # tensor of same state - rand_angle = self.rng.random() * 2 * np.pi - rand_plane = self.rng.choice(np.array([i for i in graphix.pauli.Plane])) + rand_angle = fx_rng.random() * 2 * np.pi + rand_plane = fx_rng.choice(np.array([i for i in graphix.pauli.Plane])) state = PlanarState(plane=rand_plane, angle=rand_angle) vec = Statevec(nqubit=nqb, data=state) sv_list = [state.get_statevector() for _ in range(nqb)] @@ -74,8 +68,8 @@ def test_default_tensor_success(self): assert len(vec.dims()) == nqb # tensor of different states - rand_angles = self.rng.random(nqb) * 2 * np.pi - rand_planes = self.rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) + rand_angles = fx_rng.random(nqb) * 2 * np.pi + rand_planes = fx_rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) states = [PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] vec = Statevec(nqubit=nqb, data=states) sv_list = [state.get_statevector() for state in states] @@ -83,48 +77,48 @@ def test_default_tensor_success(self): assert np.allclose(vec.psi, sv.reshape((2,) * nqb)) assert len(vec.dims()) == nqb - def test_data_success(self): - nqb = self.rng.integers(2, 5) - l = 2**nqb - rand_vec = self.rng.random(l) + 1j * self.rng.random(l) + def test_data_success(self, fx_rng: np.random.Generator) -> None: + nqb = fx_rng.integers(2, 5) + length = 2**nqb + rand_vec = fx_rng.random(length) + 1j * fx_rng.random(length) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) vec = Statevec(data=rand_vec) assert np.allclose(vec.psi, rand_vec.reshape((2,) * nqb)) assert len(vec.dims()) == nqb # fail: incorrect len - def test_data_dim_fail(self): - l = 5 - rand_vec = self.rng.random(l) + 1j * self.rng.random(l) + def test_data_dim_fail(self, fx_rng: np.random.Generator) -> None: + length = 5 + rand_vec = fx_rng.random(length) + 1j * fx_rng.random(length) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) with pytest.raises(ValueError): - vec = Statevec(data=rand_vec) + _vec = Statevec(data=rand_vec) # fail: with less qubit than number of qubits inferred from a correct state vect - def test_data_dim_fail_mismatch(self): + def test_data_dim_fail_mismatch(self, fx_rng: np.random.Generator) -> None: nqb = 3 - rand_vec = self.rng.random(2**nqb) + 1j * self.rng.random(2**nqb) + rand_vec = fx_rng.random(2**nqb) + 1j * fx_rng.random(2**nqb) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) with pytest.raises(ValueError): - vec = Statevec(nqubit=2, data=rand_vec) + _vec = Statevec(nqubit=2, data=rand_vec) # fail: not normalized - def test_data_norm_fail(self): - nqb = self.rng.integers(2, 5) - l = 2**nqb - rand_vec = self.rng.random(l) + 1j * self.rng.random(l) + def test_data_norm_fail(self, fx_rng: np.random.Generator) -> None: + nqb = fx_rng.integers(2, 5) + length = 2**nqb + rand_vec = fx_rng.random(length) + 1j * fx_rng.random(length) with pytest.raises(ValueError): - vec = Statevec(data=rand_vec) + _vec = Statevec(data=rand_vec) - def test_defaults_to_one(self): + def test_defaults_to_one(self) -> None: vec = Statevec() assert len(vec.dims()) == 1 # try copying Statevec input - def test_copy_success(self): - nqb = self.rng.integers(2, 5) - l = 2**nqb - rand_vec = self.rng.random(l) + 1j * self.rng.random(l) + def test_copy_success(self, fx_rng: np.random.Generator) -> None: + nqb = fx_rng.integers(2, 5) + length = 2**nqb + rand_vec = fx_rng.random(length) + 1j * fx_rng.random(length) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) test_vec = Statevec(data=rand_vec) # try to copy it @@ -134,12 +128,12 @@ def test_copy_success(self): assert len(vec.dims()) == len(test_vec.dims()) # try calling with incorrect number of qubits compared to inferred one - def test_copy_fail(self): - nqb = self.rng.integers(2, 5) - l = 2**nqb - rand_vec = self.rng.random(l) + 1j * self.rng.random(l) + def test_copy_fail(self, fx_rng: np.random.Generator) -> None: + nqb = fx_rng.integers(2, 5) + length = 2**nqb + rand_vec = fx_rng.random(length) + 1j * fx_rng.random(length) rand_vec /= np.sqrt(np.sum(np.abs(rand_vec) ** 2)) test_vec = Statevec(data=rand_vec) with pytest.raises(ValueError): - vec = Statevec(nqubit=l - 1, data=test_vec) + _vec = Statevec(nqubit=length - 1, data=test_vec) diff --git a/tests/test_statevec_backend.py b/tests/test_statevec_backend.py index 8e0fc3254..f49fbc44d 100644 --- a/tests/test_statevec_backend.py +++ b/tests/test_statevec_backend.py @@ -62,7 +62,7 @@ class TestStatevecNew: # more tests not really needed since redundant with Statevec constructor tests # test initialization only - def test_init_success(self, hadamardpattern, fx_rng: Generator): + def test_init_success(self, hadamardpattern, fx_rng: Generator) -> None: # plus state (default) backend = StatevectorBackend(hadamardpattern) vec = Statevec(nqubit=1) @@ -87,7 +87,7 @@ def test_init_success(self, hadamardpattern, fx_rng: Generator): # data input and Statevec input - def test_init_fail(self, hadamardpattern, fx_rng: Generator): + def test_init_fail(self, hadamardpattern, fx_rng: Generator) -> None: rand_angle = fx_rng.random(2) * 2 * np.pi rand_plane = fx_rng.choice(np.array([i for i in graphix.pauli.Plane]), 2) From 7e56152b76ad275a8bfed2498e7ec7901b30afb4 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Fri, 31 May 2024 18:13:30 +0200 Subject: [PATCH 35/48] Use NDArray in types and remove useless coercion to np.ndarray --- graphix/states.py | 9 +++++---- tests/test_density_matrix.py | 2 +- 2 files changed, 6 insertions(+), 5 deletions(-) diff --git a/graphix/states.py b/graphix/states.py index e5ec9299b..0b0b8b232 100644 --- a/graphix/states.py +++ b/graphix/states.py @@ -5,6 +5,7 @@ import abc import numpy as np +import numpy.typing as npt import pydantic import graphix.pauli @@ -19,10 +20,10 @@ class State(abc.ABC): """ @abc.abstractmethod - def get_statevector(self) -> np.ndarray: + def get_statevector(self) -> npt.NDArray: pass - def get_densitymatrix(self) -> np.ndarray: + def get_densitymatrix(self) -> npt.NDArray: # return DM in 2**n x 2**n dim (2x2 here) return np.outer(self.get_statevector(), self.get_statevector().conj()) @@ -45,10 +46,10 @@ class PlanarState(pydantic.BaseModel, State): plane: graphix.pauli.Plane angle: float - def __repr__(self): + def __repr__(self) -> str: return f"PlanarState object defined in plane {self.plane} with angle {self.angle}." - def get_statevector(self) -> np.ndarray: + def get_statevector(self) -> npt.NDArray: if self.plane == graphix.pauli.Plane.XY: return np.array([1, np.exp(1j * self.angle)]) / np.sqrt(2) diff --git a/tests/test_density_matrix.py b/tests/test_density_matrix.py index ded866db8..1ed06a112 100644 --- a/tests/test_density_matrix.py +++ b/tests/test_density_matrix.py @@ -103,7 +103,7 @@ def test_init_with_state_sucess(self, fx_rng: Generator) -> None: nqb = fx_rng.integers(2, 5) print(f"nqb is {nqb}") rand_angles = fx_rng.random(nqb) * 2 * np.pi - rand_planes = fx_rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) + rand_planes = fx_rng.choice([i for i in graphix.pauli.Plane], nqb) states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] vec = Statevec(data=states) # flattens input! From 28f6aaae246e024e2c6d651bab1efca320833fd3 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Fri, 31 May 2024 18:20:47 +0200 Subject: [PATCH 36/48] Use `pytest_configure` instead of a fixture for constant depth --- tests/conftest.py | 13 ++++++------- 1 file changed, 6 insertions(+), 7 deletions(-) diff --git a/tests/conftest.py b/tests/conftest.py index 8b51963e3..38d1e25c2 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -7,6 +7,10 @@ SEED = 25 +def pytest_configure(): + pytest.depth = 1 + + @pytest.fixture() def fx_rng() -> Generator: return Generator(PCG64(SEED)) @@ -30,13 +34,8 @@ def nqb(fx_rng: Generator): @pytest.fixture -def depth(): - return 1 - - -@pytest.fixture -def rand_circ(nqb, depth, fx_rng: Generator): - return tests.random_circuit.get_rand_circuit(nqb, depth, fx_rng) +def rand_circ(nqb, fx_rng: Generator): + return tests.random_circuit.get_rand_circuit(nqb, pytest.depth, fx_rng) @pytest.fixture From 1adc0370821dcf1ccfc6ce25b93db058e50d73b8 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Fri, 31 May 2024 18:27:30 +0200 Subject: [PATCH 37/48] Add missing type annotations --- tests/conftest.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/tests/conftest.py b/tests/conftest.py index 38d1e25c2..5f0817267 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -7,7 +7,7 @@ SEED = 25 -def pytest_configure(): +def pytest_configure() -> None: pytest.depth = 1 @@ -22,22 +22,22 @@ def fx_bg() -> PCG64: @pytest.fixture -def hadamardpattern(): +def hadamardpattern() -> graphix.pattern.Pattern: circ = graphix.transpiler.Circuit(1) circ.h(0) return circ.transpile().pattern @pytest.fixture -def nqb(fx_rng: Generator): +def nqb(fx_rng: Generator) -> int: return fx_rng.integers(2, 5) @pytest.fixture -def rand_circ(nqb, fx_rng: Generator): +def rand_circ(nqb, fx_rng: Generator) -> graphix.transpiler.Circuit: return tests.random_circuit.get_rand_circuit(nqb, pytest.depth, fx_rng) @pytest.fixture -def randpattern(rand_circ): +def randpattern(rand_circ) -> graphix.pattern.Pattern: return rand_circ.transpile().pattern From 51335dcf152ddd168a2b429d78fd6070c8dcfb4a Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Fri, 31 May 2024 18:38:03 +0200 Subject: [PATCH 38/48] Truncate matrix in error messages --- graphix/linalg_validations.py | 9 ++++++++- 1 file changed, 8 insertions(+), 1 deletion(-) diff --git a/graphix/linalg_validations.py b/graphix/linalg_validations.py index 4e914ebcf..a9d3fb61e 100644 --- a/graphix/linalg_validations.py +++ b/graphix/linalg_validations.py @@ -18,6 +18,13 @@ def check_square(matrix: np.ndarray) -> bool: return True +def truncate(s: str, max_length: int = 80, ellipsis: str = "...") -> str: + "Auxilliary function to truncate a long string for formatting error messages." + if len(s) <= max_length: + return s + return s[: max_length - len(ellipsis)] + ellipsis + + def check_psd(matrix: np.ndarray, tol: float = 1e-15) -> bool: """ check if a density matrix is positive semidefinite by diagonalizing. @@ -33,7 +40,7 @@ def check_psd(matrix: np.ndarray, tol: float = 1e-15) -> bool: evals = np.linalg.eigvalsh(matrix) if not all(evals >= -tol): - raise ValueError("The matrix {matrix} is not positive semi-definite.") + raise ValueError("The matrix {truncate(str(matrix))} is not positive semi-definite.") return True From cabd019db292da2032f440be7c3e492851a1a997 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Fri, 31 May 2024 18:41:12 +0200 Subject: [PATCH 39/48] Rename SV_Data into Data --- graphix/sim/statevec.py | 6 +++--- graphix/transpiler.py | 10 +++++----- 2 files changed, 8 insertions(+), 8 deletions(-) diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 4bfd23469..1f8e496f6 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -21,7 +21,7 @@ class StatevectorBackend(graphix.sim.base_backend.Backend): def __init__( self, pattern, - input_state: SV_Data = graphix.states.BasicStates.PLUS, + input_state: Data = graphix.states.BasicStates.PLUS, max_qubit_num=20, pr_calc=True, ): @@ -202,7 +202,7 @@ class Statevec: def __init__( self, - data: SV_Data = graphix.states.BasicStates.PLUS, + data: Data = graphix.states.BasicStates.PLUS, nqubit: graphix.types.PositiveOrNullInt | None = None, ): """Initialize statevector objects. The behaviour is as follows. `data` can be: @@ -498,7 +498,7 @@ def _get_statevec_norm(psi): return np.sqrt(np.sum(psi.flatten().conj() * psi.flatten())) -SV_Data = ( +Data = ( graphix.states.State | Statevec | collections.abc.Iterable[graphix.states.State] diff --git a/graphix/transpiler.py b/graphix/transpiler.py index 1b668fd7e..961d8c565 100644 --- a/graphix/transpiler.py +++ b/graphix/transpiler.py @@ -16,7 +16,7 @@ import graphix.sim.base_backend from graphix.ops import Ops from graphix.pattern import Pattern -from graphix.sim.statevec import Statevec, SV_Data +import graphix.sim.statevec @dataclasses.dataclass @@ -41,7 +41,7 @@ class SimulateResult: classical_measures : tuple[int,...], classical measures """ - statevec: Statevec + statevec: graphix.sim.statevec.Statevec classical_measures: tuple[int, ...] @@ -1358,7 +1358,7 @@ def _sort_outputs(self, pattern: Pattern, output_nodes: Sequence[int]): elif cmd[1] in old_out: cmd[1] = output_nodes[old_out.index(cmd[1])] - def simulate_statevector(self, input_state: SV_Data | None = None) -> SimulateResult: + def simulate_statevector(self, input_state: graphix.sim.statevec.Data | None = None) -> SimulateResult: """Run statevector simulation of the gate sequence, using graphix.Statevec Parameters @@ -1372,9 +1372,9 @@ def simulate_statevector(self, input_state: SV_Data | None = None) -> SimulateRe """ if input_state is None: - state = Statevec(nqubit=self.width) + state = graphix.sim.statevec.Statevec(nqubit=self.width) else: - state = Statevec(nqubit=self.width, data=input_state) + state = graphix.sim.statevec.Statevec(nqubit=self.width, data=input_state) classical_measures = [] From 55331b5179c05408ad41febea75d03c21cf16d00 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Mon, 3 Jun 2024 09:32:43 +0200 Subject: [PATCH 40/48] Use DEPTH constant instead of pytest.depth --- tests/conftest.py | 7 ++----- 1 file changed, 2 insertions(+), 5 deletions(-) diff --git a/tests/conftest.py b/tests/conftest.py index 5f0817267..055509c1e 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -5,10 +5,7 @@ import tests.random_circuit SEED = 25 - - -def pytest_configure() -> None: - pytest.depth = 1 +DEPTH = 1 @pytest.fixture() @@ -35,7 +32,7 @@ def nqb(fx_rng: Generator) -> int: @pytest.fixture def rand_circ(nqb, fx_rng: Generator) -> graphix.transpiler.Circuit: - return tests.random_circuit.get_rand_circuit(nqb, pytest.depth, fx_rng) + return tests.random_circuit.get_rand_circuit(nqb, DEPTH, fx_rng) @pytest.fixture From c406422747d9f25400ba1b941c407bb76887a37f Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Mon, 17 Jun 2024 14:24:07 +0200 Subject: [PATCH 41/48] Black and isort --- examples/MBQCvqe.py | 9 ++++++++- graphix/transpiler.py | 2 +- 2 files changed, 9 insertions(+), 2 deletions(-) diff --git a/examples/MBQCvqe.py b/examples/MBQCvqe.py index 05eefdbc2..37c0a5e7b 100644 --- a/examples/MBQCvqe.py +++ b/examples/MBQCvqe.py @@ -22,9 +22,11 @@ import networkx as nx import numpy as np from scipy.optimize import minimize + from graphix import Circuit from graphix.simulator import PatternSimulator + # %% # Define the Hamiltonian for the VQE problem (Example: H = Z0Z1 + X0 + X1) def create_hamiltonian(): @@ -33,6 +35,7 @@ def create_hamiltonian(): H = np.kron(Z, Z) + np.kron(X, np.eye(2)) + np.kron(np.eye(2), X) return H + # %% # Function to build the VQE circuit def build_vqe_circuit(n_qubits, params): @@ -45,6 +48,7 @@ def build_vqe_circuit(n_qubits, params): circuit.cnot(i, i + 1) return circuit + # %% class MBQCVQE: def __init__(self, n_qubits, hamiltonian): @@ -85,6 +89,7 @@ def compute_energy(self, params): energy = tn.expectation_value(self.hamiltonian, qubit_indices=range(self.n_qubits)) return energy + # %% # Set parameters for VQE n_qubits = 2 @@ -94,18 +99,20 @@ def compute_energy(self, params): # Instantiate the MBQCVQE class mbqc_vqe = MBQCVQE(n_qubits, hamiltonian) + # %% # Define the cost function def cost_function(params): return mbqc_vqe.compute_energy(params) + # %% # Random initial parameters initial_params = np.random.rand(n_qubits * 3) # %% # Perform the optimization using COBYLA -result = minimize(cost_function, initial_params, method='COBYLA', options={'maxiter': 100}) +result = minimize(cost_function, initial_params, method="COBYLA", options={"maxiter": 100}) print(f"Optimized parameters: {result.x}") print(f"Optimized energy: {result.fun}") diff --git a/graphix/transpiler.py b/graphix/transpiler.py index 961d8c565..d7c31858e 100644 --- a/graphix/transpiler.py +++ b/graphix/transpiler.py @@ -14,9 +14,9 @@ import graphix.pauli import graphix.sim.base_backend +import graphix.sim.statevec from graphix.ops import Ops from graphix.pattern import Pattern -import graphix.sim.statevec @dataclasses.dataclass From 858bf0ab2f66970542d05fd72e557e637403900f Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Mon, 17 Jun 2024 14:40:09 +0200 Subject: [PATCH 42/48] Fix for Python <3.10 --- graphix/sim/density_matrix.py | 29 ++++++++++++++++++++--------- graphix/sim/statevec.py | 23 ++++++++++++++++------- 2 files changed, 36 insertions(+), 16 deletions(-) diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index f4b7719c4..4f6ac6ca3 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -7,6 +7,7 @@ import collections import numbers +import typing from copy import deepcopy import numpy as np @@ -437,12 +438,22 @@ def finalize(self): self.sort_qubits() self.state.normalize() - -Data = ( - graphix.states.State - | DensityMatrix - | Statevec - | collections.abc.Iterable[graphix.states.State] - | collections.abc.Iterable[numbers.Number] - | collections.abc.Iterable[collections.abc.Iterable[numbers.Number]] -) +## Python <3.10: +## TypeError: unsupported operand type(s) for |: 'ABCMeta' and 'type' +## TypeError: 'ABCMeta' object is not subscriptable +#Data = ( +# graphix.states.State +# | DensityMatrix +# | Statevec +# | collections.abc.Iterable[graphix.states.State] +# | collections.abc.Iterable[numbers.Number] +# | collections.abc.Iterable[collections.abc.Iterable[numbers.Number]] +#) +Data = typing.Union[ + graphix.states.State, + DensityMatrix, + Statevec, + typing.Iterable[graphix.states.State], + typing.Iterable[numbers.Number], + typing.Iterable[typing.Iterable[numbers.Number]], +] diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 1f8e496f6..37e996066 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -3,6 +3,7 @@ import collections import functools import numbers +import typing from copy import deepcopy import numpy as np @@ -497,10 +498,18 @@ def _get_statevec_norm(psi): """returns norm of the state""" return np.sqrt(np.sum(psi.flatten().conj() * psi.flatten())) - -Data = ( - graphix.states.State - | Statevec - | collections.abc.Iterable[graphix.states.State] - | collections.abc.Iterable[numbers.Number] -) +## Python <3.10: +## TypeError: unsupported operand type(s) for |: 'ABCMeta' and 'type' +## TypeError: 'ABCMeta' object is not subscriptable +#Data = ( +# graphix.states.State +# | Statevec +# | collections.abc.Iterable[graphix.states.State] +# | collections.abc.Iterable[numbers.Number] +#) +Data = typing.Union[ + graphix.states.State, + Statevec, + typing.Iterable[graphix.states.State], + typing.Iterable[numbers.Number], +] From dcf529e21af12298c89cc0a244e83cad0b78541b Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Mon, 17 Jun 2024 14:45:31 +0200 Subject: [PATCH 43/48] Black --- graphix/sim/density_matrix.py | 5 +++-- graphix/sim/statevec.py | 5 +++-- 2 files changed, 6 insertions(+), 4 deletions(-) diff --git a/graphix/sim/density_matrix.py b/graphix/sim/density_matrix.py index 4f6ac6ca3..fa4d21ac2 100644 --- a/graphix/sim/density_matrix.py +++ b/graphix/sim/density_matrix.py @@ -438,17 +438,18 @@ def finalize(self): self.sort_qubits() self.state.normalize() + ## Python <3.10: ## TypeError: unsupported operand type(s) for |: 'ABCMeta' and 'type' ## TypeError: 'ABCMeta' object is not subscriptable -#Data = ( +# Data = ( # graphix.states.State # | DensityMatrix # | Statevec # | collections.abc.Iterable[graphix.states.State] # | collections.abc.Iterable[numbers.Number] # | collections.abc.Iterable[collections.abc.Iterable[numbers.Number]] -#) +# ) Data = typing.Union[ graphix.states.State, DensityMatrix, diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 37e996066..0ef2ce9ec 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -498,15 +498,16 @@ def _get_statevec_norm(psi): """returns norm of the state""" return np.sqrt(np.sum(psi.flatten().conj() * psi.flatten())) + ## Python <3.10: ## TypeError: unsupported operand type(s) for |: 'ABCMeta' and 'type' ## TypeError: 'ABCMeta' object is not subscriptable -#Data = ( +# Data = ( # graphix.states.State # | Statevec # | collections.abc.Iterable[graphix.states.State] # | collections.abc.Iterable[numbers.Number] -#) +# ) Data = typing.Union[ graphix.states.State, Statevec, From 20d02379bb9167c10958bb497fb0c0f2144b2acd Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Tue, 18 Jun 2024 13:26:09 +0200 Subject: [PATCH 44/48] Remove reference to `SV_Data` in doc-comment --- graphix/sim/statevec.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 0ef2ce9ec..073f62b1a 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -220,7 +220,7 @@ def __init__( :param data: input data to prepare the state. Can be a classical description or a numerical input, defaults to graphix.states.BasicStates.PLUS - :type data: SV_Data, optional + :type data: Data, optional :param nqubit: number of qubits to prepare, defaults to None :type nqubit: int, optional """ From 26408450d87f731b32ce0561fbada996b68e9c43 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Tue, 18 Jun 2024 14:01:50 +0200 Subject: [PATCH 45/48] Add tests for density matrix and tensor network backends --- tests/test_pattern.py | 37 +++++++++++++++++++++++-------------- 1 file changed, 23 insertions(+), 14 deletions(-) diff --git a/tests/test_pattern.py b/tests/test_pattern.py index 8ebb14d81..70e6f0df9 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -11,6 +11,8 @@ import graphix.states import tests.random_circuit as rc from graphix.pattern import CommandNode, Pattern +from graphix.sim.density_matrix import DensityMatrix +from graphix.sim.statevec import Statevec from graphix.simulator import PatternSimulator from graphix.transpiler import Circuit @@ -20,6 +22,15 @@ from numpy.random import PCG64, Generator +def compare_backend_result_with_statevec(backend: str, backend_state, statevec: Statevec) -> float: + if backend == "statevector": + return np.abs(np.dot(backend_state.flatten().conjugate(), statevec.flatten())) + elif backend == "densitymatrix": + return np.abs(np.dot(backend_state.rho.flatten().conjugate(), DensityMatrix(statevec).rho.flatten())) + else: + raise NotImplementedError(backend) + + class TestPattern: # this fails without behaviour modification def test_manual_generation(self) -> None: @@ -140,10 +151,7 @@ def test_pauli_measurement_random_circuit( pattern.minimize_space() state = circuit.simulate_statevector().statevec state_mbqc = pattern.simulate_pattern(backend) - if backend == "statevector": - assert np.abs(np.dot(state_mbqc.flatten().conjugate(), state.flatten())) == pytest.approx(1) - elif backend == "densitymatrix": - assert np.allclose(state_mbqc.rho, graphix.sim.density_matrix.DensityMatrix(state.flatten()).rho) + assert compare_backend_result_with_statevec(backend, state_mbqc, state) == pytest.approx(1) @pytest.mark.parametrize("jumps", range(1, 11)) def test_pauli_measurement_leave_input_random_circuit( @@ -526,22 +534,23 @@ def test_pauli_measurement_end_with_measure(self) -> None: p.add(["M", 1, "XY", 0, [], []]) p.perform_pauli_measurements() - -# for testing with arbitrary inputs -# SV and DM backend -class TestPatternSim: - def test_sv_sim(self, fx_rng: Generator, nqb, rand_circ) -> None: + @pytest.mark.parametrize("backend", ["statevector", "densitymatrix"]) + def test_arbitrary_inputs(self, fx_rng: Generator, nqb: int, rand_circ: Circuit, backend: str) -> None: rand_angles = fx_rng.random(nqb) * 2 * np.pi rand_planes = fx_rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] randpattern = rand_circ.transpile().pattern - out = randpattern.simulate_pattern(backend="statevector", input_state=states) - + out = randpattern.simulate_pattern(backend=backend, input_state=states) out_circ = rand_circ.simulate_statevector(input_state=states).statevec - np.testing.assert_almost_equal(np.abs(np.dot(out.psi.flatten().conjugate(), out_circ.psi.flatten())), 1) + assert compare_backend_result_with_statevec(backend, out, out_circ) == pytest.approx(1) - def test_dm_sim(self) -> None: - pass + def test_arbitrary_inputs_tn(self, fx_rng: Generator, nqb: int, rand_circ: Circuit) -> None: + rand_angles = fx_rng.random(nqb) * 2 * np.pi + rand_planes = fx_rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) + states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] + randpattern = rand_circ.transpile().pattern + with pytest.raises(TypeError): # __init__() got an unexpected keyword argument 'input_state' + randpattern.simulate_pattern(backend="tensornetwork", graph_prep="sequential", input_state=states) def assert_equal_edge(edge: Sequence[int], ref: Sequence[int]) -> bool: From 9e363b15e91d485ffedbb790d3f588abe998402f Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Tue, 18 Jun 2024 14:19:10 +0200 Subject: [PATCH 46/48] Increase number of shots in test_transpiler to reduce failures --- tests/test_transpiler.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_transpiler.py b/tests/test_transpiler.py index ff2c0a878..d6647cfe2 100644 --- a/tests/test_transpiler.py +++ b/tests/test_transpiler.py @@ -165,6 +165,6 @@ def simulate_and_measure() -> int: assert circuit_simulate.classical_measures[0] == (circuit_simulate.statevec.psi[0][1].imag > 0) return circuit_simulate.classical_measures[0] - nb_shots = 1000 + nb_shots = 10000 count = sum(1 for _ in range(nb_shots) if simulate_and_measure()) assert abs(count - nb_shots / 2) < nb_shots / 20 From 05aa7c949cad0daa4bc4a0f066603bbf46abb374 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Tue, 18 Jun 2024 22:58:37 +0200 Subject: [PATCH 47/48] NotImplementedError for Pauli preprocessing and TN --- graphix/sim/statevec.py | 4 +++- graphix/sim/tensornet.py | 7 ++++++- tests/test_pattern.py | 2 +- 3 files changed, 10 insertions(+), 3 deletions(-) diff --git a/graphix/sim/statevec.py b/graphix/sim/statevec.py index 073f62b1a..755565cc1 100644 --- a/graphix/sim/statevec.py +++ b/graphix/sim/statevec.py @@ -43,7 +43,9 @@ def __init__( # assert len(pattern.output_nodes) > 0 self.pattern = pattern if pattern._pauli_preprocessed and input_state != graphix.states.BasicStates.PLUS: - raise ValueError("Pauli preprocessing is currently only available when inputs are initialized in |+> state") + raise NotImplementedError( + "Pauli preprocessing is currently only available when inputs are initialized in |+> state (see https://github.com/TeamGraphix/graphix/issues/168 )." + ) self.results = deepcopy(pattern.results) self.state = None self.node_index = [] diff --git a/graphix/sim/tensornet.py b/graphix/sim/tensornet.py index 985b0d187..fc9928cb8 100644 --- a/graphix/sim/tensornet.py +++ b/graphix/sim/tensornet.py @@ -16,7 +16,7 @@ class TensorNetworkBackend: Executes the measurement pattern using TN expression of graph states. """ - def __init__(self, pattern, graph_prep="auto", **kwargs): + def __init__(self, pattern, graph_prep="auto", input_state=BasicStates.PLUS, **kwargs): """ Parameters @@ -33,8 +33,13 @@ def __init__(self, pattern, graph_prep="auto", **kwargs): In this strategy, All N and E commands executed sequentially. 'auto'(default) : Automatically select a preparation strategy based on the max degree of a graph + input_state : preparation for input states (only BasicStates.PLUS is supported for tensor networks yet), **kwargs : Additional keyword args to be passed to quimb.tensor.TensorNetwork. """ + if input_state != BasicStates.PLUS: + raise NotImplementedError( + "TensorNetworkBackend currently only supports |+> input state (see https://github.com/TeamGraphix/graphix/issues/167 )." + ) self.pattern = pattern self.output_nodes = pattern.output_nodes self.results = deepcopy(pattern.results) diff --git a/tests/test_pattern.py b/tests/test_pattern.py index 70e6f0df9..3608acd02 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -549,7 +549,7 @@ def test_arbitrary_inputs_tn(self, fx_rng: Generator, nqb: int, rand_circ: Circu rand_planes = fx_rng.choice(np.array([i for i in graphix.pauli.Plane]), nqb) states = [graphix.states.PlanarState(plane=i, angle=j) for i, j in zip(rand_planes, rand_angles)] randpattern = rand_circ.transpile().pattern - with pytest.raises(TypeError): # __init__() got an unexpected keyword argument 'input_state' + with pytest.raises(NotImplementedError): randpattern.simulate_pattern(backend="tensornetwork", graph_prep="sequential", input_state=states) From 19dfd8a272997970840dd516031e31bb68c28803 Mon Sep 17 00:00:00 2001 From: Thierry Martinez Date: Wed, 19 Jun 2024 11:10:29 +0200 Subject: [PATCH 48/48] Lower requirement for test_empty_output_nodes --- tests/test_pattern.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_pattern.py b/tests/test_pattern.py index 3608acd02..bf44b3695 100644 --- a/tests/test_pattern.py +++ b/tests/test_pattern.py @@ -96,7 +96,7 @@ def simulate_and_measure(): nb_shots = 1000 nb_ones = sum(1 for _ in range(nb_shots) if simulate_and_measure()) - assert abs(nb_ones - nb_shots / 2) < nb_shots / 20 + assert abs(nb_ones - nb_shots / 2) < nb_shots / 10 def test_minimize_space_graph_maxspace_with_flow(self, fx_rng: Generator) -> None: max_qubits = 20