From b9c14103b821d16204d446795fd621d6e7ff5e8b Mon Sep 17 00:00:00 2001 From: Synge Todo Date: Thu, 16 Nov 2023 23:56:35 +0900 Subject: [PATCH] remove mps_t --- src/qailo/__init__.py | 3 +- src/qailo/mps_p/__init__.py | 4 +- src/qailo/mps_p/mps_p.py | 4 +- src/qailo/mps_p/projector.py | 50 +---------- src/qailo/mps_t/__init__.py | 9 -- src/qailo/mps_t/mps_t.py | 138 ------------------------------- src/qailo/mps_t/product_state.py | 19 ----- src/qailo/mps_t/tpool.py | 119 -------------------------- test/mps/test_apply.py | 56 +++++-------- test/mps/test_canonical.py | 4 +- test/mps/test_move.py | 4 +- test/mps/test_mps.py | 3 - test/mps_p/test_projector.py | 114 +------------------------ test/mps_t/test_mps_t.py | 40 --------- 14 files changed, 34 insertions(+), 533 deletions(-) delete mode 100644 src/qailo/mps_t/__init__.py delete mode 100644 src/qailo/mps_t/mps_t.py delete mode 100644 src/qailo/mps_t/product_state.py delete mode 100644 src/qailo/mps_t/tpool.py delete mode 100644 test/mps_t/test_mps_t.py diff --git a/src/qailo/__init__.py b/src/qailo/__init__.py index 4273535..6475731 100644 --- a/src/qailo/__init__.py +++ b/src/qailo/__init__.py @@ -1,4 +1,4 @@ -from . import alg, mps, mps_p, mps_t, operator, state_vector, util +from . import alg, mps, mps_p, operator, state_vector, util from . import operator as op from . import state_vector as sv from ._version import version @@ -8,7 +8,6 @@ alg, mps, mps_p, - mps_t, operator, state_vector, util, diff --git a/src/qailo/mps_p/__init__.py b/src/qailo/mps_p/__init__.py index 6b4a888..1152c5e 100644 --- a/src/qailo/mps_p/__init__.py +++ b/src/qailo/mps_p/__init__.py @@ -1,10 +1,10 @@ from .mps_p import MPS_P from .product_state import one, product_state, zero -from .projector import compact_projector +from .projector import projector __all__ = [ MPS_P, - compact_projector, + projector, one, product_state, zero, diff --git a/src/qailo/mps_p/mps_p.py b/src/qailo/mps_p/mps_p.py index ac9136c..d6cd39e 100644 --- a/src/qailo/mps_p/mps_p.py +++ b/src/qailo/mps_p/mps_p.py @@ -5,7 +5,7 @@ from ..mps.svd import tensor_svd from ..mps.type import MPS from ..operator import type as op -from .projector import compact_projector +from .projector import projector class MPS_P(MPS): @@ -126,6 +126,6 @@ def _apply_two(self, p, s, maxdim=None, reverse=False): t1 = np.einsum(t1, [0, 4, 3], p0, [2, 4, 1]) tt0 = np.einsum(self.env[s], [0, 4], t0, [4, 1, 2, 3]) tt1 = np.einsum(t1, [0, 1, 2, 4], self.env[s + 2], [4, 3]) - _, WLh, WR = compact_projector(tt0, [0, 1, 4, 5], tt1, [5, 4, 2, 3], maxdim) + _, WLh, WR = projector(tt0, [0, 1, 4, 5], tt1, [5, 4, 2, 3], maxdim) self.tensors[s] = np.einsum(t0, [0, 1, 3, 4], WR, [3, 4, 2]) self.tensors[s + 1] = np.einsum(WLh, [3, 4, 0], t1, [4, 3, 1, 2]) diff --git a/src/qailo/mps_p/projector.py b/src/qailo/mps_p/projector.py index 3b6b2f4..67f6301 100644 --- a/src/qailo/mps_p/projector.py +++ b/src/qailo/mps_p/projector.py @@ -52,7 +52,7 @@ def collect_legs(ss0, ss1): return legs0L, legs0R, legs1L, legs1R -def compact_projector(T0, ss0_in, T1, ss1_in, nkeep=None, tol=1e-12): +def projector(T0, ss0_in, T1, ss1_in, nkeep=None, tol=1e-12): ss0, ss1 = normalize_ss(ss0_in, ss1_in) legs0L, legs0R, legs1L, legs1R = collect_legs(ss0, ss1) dim0L = np.prod([T0.shape[i] for i in legs0L]) @@ -70,51 +70,3 @@ def compact_projector(T0, ss0_in, T1, ss1_in, nkeep=None, tol=1e-12): WL = np.einsum(TT0.conj(), [0] + ss_sum, U, [0, max(ss_sum) + 1]) WR = np.einsum(TT1, [0] + ss_sum, V, [0, max(ss_sum) + 1]) return S, WL.conj(), WR - - -def _biorth(VL, VR): - assert len(VL.shape) == 2 and VL.shape == VR.shape and VL.shape[0] >= VL.shape[1] - r = VL.shape[1] - for i in range(r): - for j in range(i): - VL[:, i] = VL[:, i] - (VR[:, j].conj().T @ VL[:, i]) * VL[:, j] - VR[:, i] = VR[:, i] - (VL[:, j].conj().T @ VR[:, i]) * VR[:, j] - norm = VL[:, i].conj().T @ VR[:, i] - VL[:, i] = VL[:, i] / (norm.conj() / np.sqrt(np.abs(norm))) - VR[:, i] = VR[:, i] / np.sqrt(np.abs(norm)) - return VL, VR - - -def _full_projector(T0, ss0_in, T1, ss1_in, tol=1e-12): - """ - Return: - S: singular values - d: number of non-zero singular values - WL, WR: full projector WL* @ WR = I - """ - S, WLdh, WRd = compact_projector(T0, ss0_in, T1, ss1_in, tol=tol) - WLd = WLdh.conj() - d = S.shape[0] - assert d == WLd.shape[-1] and d == WRd.shape[-1] - dimsWL0 = WLd.shape[:-1] - dimsWR0 = WRd.shape[:-1] - m = np.prod(dimsWL0) - assert m == np.prod(dimsWR0) - if d < m: - WLd = WLd.reshape((m, d)) - WRd = WRd.reshape((m, d)) - P = np.identity(m) - WRd @ WLd.conj().T - VL = np.random.rand(m, m - d) - VR = np.random.rand(m, m - d) - VL = P.conj().T @ VL - VR = P @ VR - VL, VR = _biorth(VL, VR) - S = np.block([S, np.zeros(m - d)]) - WL = np.block([WLd, VL]) - WR = np.block([WRd, VR]) - assert S.shape == (m,) and WL.shape == (m, m) and WR.shape == (m, m) - WL = WL.reshape(dimsWL0 + (m,)) - WR = WR.reshape(dimsWR0 + (m,)) - return S, d, WL.conj(), WR - else: - return S, d, WLd.conj(), WRd diff --git a/src/qailo/mps_t/__init__.py b/src/qailo/mps_t/__init__.py deleted file mode 100644 index 17bff05..0000000 --- a/src/qailo/mps_t/__init__.py +++ /dev/null @@ -1,9 +0,0 @@ -from .mps_t import MPS_T -from .product_state import one, product_state, zero - -__all__ = [ - MPS_T, - one, - product_state, - zero, -] diff --git a/src/qailo/mps_t/mps_t.py b/src/qailo/mps_t/mps_t.py deleted file mode 100644 index 553dda9..0000000 --- a/src/qailo/mps_t/mps_t.py +++ /dev/null @@ -1,138 +0,0 @@ -from copy import deepcopy - -import numpy as np - -from ..mps.svd import tensor_svd -from ..mps.type import MPS -from ..operator import type as op -from .tpool import tpool - - -class MPS_T(MPS): - """ - MPS representation of quantum pure state - - shape of tensors: [du, dp, dl] - du: dimension of upper leg (1 for top tensor) - dp: dimension of physical leg (typically 2) - dl: dimension of lower leg (1 for bottom tensor) - - canonical position: cp in range(n) - 0 <= cp(0) <= cp(1) < n - tensors [0...cp(0)-1]: top canonical - tensors [cp(1)+1...n-1]: bottom canonical - """ - - def __init__(self, tensors): - assert isinstance(tensors, list) - n = len(tensors) - self.tp = tpool() - self.tcurrent = [] - for t in tensors: - self.tcurrent.append(self.tp._register(deepcopy(t))) - self.q2t = list(range(n)) - self.t2q = list(range(n)) - self.cp = [0, n - 1] - # canonicalization matrices - # put sentinels (1x1 identities) at t = 0 and t = n - self.env = [np.identity(1)] + [None] * (n - 1) + [np.identity(1)] - - def _num_qubits(self): - return len(self.tcurrent) - - def _state_vector(self): - n = self._num_qubits() - v = self._tensor(0) - for t in range(1, n): - ss0 = list(range(t + 1)) + [t + 3] - ss1 = [t + 3, t + 1, t + 2] - v = np.einsum(v, ss0, self._tensor(t), ss1) - v = v.reshape((2,) * n) - return np.einsum(v, self.t2q).reshape((2,) * n + (1,)) - - def _tensor(self, t): - return self.tp.tpool[self.tcurrent[t]][0] - - def _canonicalize(self, p0, p1=None): - p1 = p0 if p1 is None else p1 - assert 0 <= p0 and p0 <= p1 and p1 < self._num_qubits() - if self.cp[0] < p0: - for t in range(self.cp[0], p0): - A = np.einsum(self.env[t], [0, 3], self._tensor(t), [3, 1, 2]) - _, self.env[t + 1] = tensor_svd(A, [[0, 1], [2]], "left") - self.cp[0] = p0 - self.cp[1] = max(p0, self.cp[1]) - if self.cp[1] > p1: - for t in range(self.cp[1], p1, -1): - A = np.einsum(self._tensor(t), [0, 1, 3], self.env[t + 1], [3, 2]) - self.env[t], _ = tensor_svd(A, [[0], [1, 2]], "right") - self.cp[1] = p1 - - def _is_canonical(self): - # tensor shape - n = len(self.tcurrent) - dims = [] - assert self._tensor(0).shape[0] == 1 - dims.append(self._tensor(0).shape[0]) - for t in range(1, n - 1): - dims.append(self._tensor(t).shape[0]) - assert self._tensor(t).shape[0] == self._tensor(t - 1).shape[2] - assert self._tensor(t).shape[2] == self._tensor(t + 1).shape[0] - assert self._tensor(n - 1).shape[2] == 1 - dims.append(self._tensor(n - 1).shape[0]) - dims.append(self._tensor(n - 1).shape[2]) - - # qubit <-> tensor mapping - for q in range(n): - assert self.t2q[self.q2t[q]] == q - for t in range(n): - assert self.q2t[self.t2q[t]] == t - - # canonicality - assert self.cp[0] in range(n) - assert self.cp[1] in range(n) - A = np.identity(1) - for t in range(0, self.cp[0]): - A = np.einsum(A, [0, 3], self._tensor(t), [3, 1, 2]) - A = np.einsum(A, [2, 3, 1], self._tensor(t).conj(), [2, 3, 0]) - B = np.einsum(self.env[t + 1], [2, 1], self.env[t + 1].conj(), [2, 0]) - assert A.shape == B.shape - assert np.allclose(A, B) - A = np.identity(1) - for t in range(n - 1, self.cp[1], -1): - A = np.einsum(self._tensor(t), [0, 1, 3], A, [3, 2]) - A = np.einsum(self._tensor(t).conj(), [1, 2, 3], A, [0, 2, 3]) - B = np.einsum(self.env[t], [0, 2], self.env[t].conj(), [1, 2]) - assert np.allclose(A, B) - return True - - def _apply_one(self, p, s): - assert op.num_qubits(p) == 1 - pid = self.tp._register(p) - self.tcurrent[s] = self.tp._contract(self.tcurrent[s], [0, 3, 2], pid, [1, 3]) - self.cp[0] = min(self.cp[0], s) - self.cp[1] = max(self.cp[1], s) - - def _apply_two(self, p, s, maxdim=None, reverse=False): - """ - apply 2-qubit operator on neighboring tensors, s and s+1 - """ - self._canonicalize(s, s + 1) - tid0 = self.tcurrent[s] - tid1 = self.tcurrent[s + 1] - p0, p1 = tensor_svd(p, [[0, 2], [1, 3]]) - pid0 = self.tp._register(p0) - pid1 = self.tp._register(p1) - if not reverse: - tid0 = self.tp._contract(tid0, [0, 4, 3], pid0, [1, 4, 2]) - tid1 = self.tp._contract(tid1, [0, 4, 3], pid1, [1, 2, 4]) - else: - tid0 = self.tp._contract(tid0, [0, 4, 3], pid1, [2, 1, 4]) - tid1 = self.tp._contract(tid1, [0, 4, 3], pid0, [2, 4, 1]) - tt0 = np.einsum(self.env[s], [0, 4], self.tp.tpool[tid0][0], [4, 1, 2, 3]) - tt1 = np.einsum(self.tp.tpool[tid1][0], [0, 1, 2, 4], self.env[s + 2], [4, 3]) - _, WLhid, WRid = self.tp._projector( - tt0, [0, 1, 4, 5], tt1, [5, 4, 2, 3], maxdim - ) - self.tcurrent[s] = self.tp._contract(tid0, [0, 1, 3, 4], WRid, [3, 4, 2]) - self.tcurrent[s + 1] = self.tp._contract(WLhid, [3, 4, 0], tid1, [4, 3, 1, 2]) diff --git a/src/qailo/mps_t/product_state.py b/src/qailo/mps_t/product_state.py deleted file mode 100644 index adfa10d..0000000 --- a/src/qailo/mps_t/product_state.py +++ /dev/null @@ -1,19 +0,0 @@ -from ..mps.product_state import tensor_decomposition -from ..state_vector.state_vector import one as sv_one -from ..state_vector.state_vector import zero as sv_zero -from .mps_t import MPS_T - - -def product_state(states): - tensors = [] - for s in states: - tensors = tensors + tensor_decomposition(s) - return MPS_T(tensors) - - -def zero(n=1): - return product_state([sv_zero()] * n) - - -def one(n=1): - return product_state([sv_one()] * n) diff --git a/src/qailo/mps_t/tpool.py b/src/qailo/mps_t/tpool.py deleted file mode 100644 index 0598e3c..0000000 --- a/src/qailo/mps_t/tpool.py +++ /dev/null @@ -1,119 +0,0 @@ -import numpy as np - -from ..mps_p.projector import _full_projector - - -class tpool: - def __init__(self): - self.tpool = [] - self.gpool = [] - - def _register(self, tensor): - id = len(self.tpool) - self.tpool.append([tensor, "initial", None]) - return id - - def _contract(self, tid0, ss0, tid1, ss1, ss2=None): - id = len(self.tpool) - if ss2 is None: - t = np.einsum(self.tpool[tid0][0], ss0, self.tpool[tid1][0], ss1) - else: - t = np.einsum(self.tpool[tid0][0], ss0, self.tpool[tid1][0], ss1, ss2) - self.tpool.append([t, "product", tid0, ss0, tid1, ss1, ss2]) - return id - - def _trim(self, T, d): - """trim last dimension of tensor""" - assert d <= T.shape[-1] - shape = T.shape - T = T.reshape((np.prod(shape[:-1]), shape[-1])) - return T[:, :d].reshape(shape[:-1] + (d,)) - - def _projector(self, t0, ss0, t1, ss1, maxdim=None): - S, d, WLh, WR = _full_projector(t0, ss0, t1, ss1) - d = d if maxdim is None else min(d, maxdim) - gid = len(self.gpool) - self.gpool.append([S, d, WLh, WR]) - lid = len(self.tpool) - self.tpool.append([self._trim(WLh, d), "squeezer L", gid]) - rid = len(self.tpool) - self.tpool.append([self._trim(WR, d), "squeezer R", gid]) - return gid, lid, rid - - def _dump(self, prefix): - import json - - dic = {"prefix": prefix} - tlist = [] - for id, tp in enumerate(self.tpool): - m = {} - m["id"] = id - m["type"] = tp[1] - m["shape"] = tp[0].shape - if tp[1] == "product": - m["from 0"] = tp[2] - m["subscripts 0"] = tp[3] - m["from 1"] = tp[4] - m["subscripts 1"] = tp[5] - if tp[6] is not None: - m["subscripts"] = tp[6] - elif tp[1] == "squeezer L" or tp[1] == "squeezer R": - m["from"] = tp[2] - tlist.append(m) - if tp[1] == "initial": - np.save(f"{prefix}-tensor-{id}", tp[0]) - dic["tensor"] = tlist - glist = [] - for id, gp in enumerate(self.gpool): - m = {} - m["id"] = id - m["dim from"] = gp[0].shape[0] - m["dim to"] = gp[1] - m["shape L"] = gp[2].shape - m["shape R"] = gp[3].shape - glist.append(m) - np.savez(f"{prefix}-generator-{id}", s=gp[0], wl=gp[1], wr=gp[2]) - dic["generator"] = glist - - with open(prefix + "-graph.json", mode="w") as f: - json.dump(dic, f, indent=2) - - def _load(self, prefix): - import json - - with open(prefix + "-graph.json", mode="r") as f: - dic = json.load(f) - print(dic) - assert "prefix" in dic and dic["prefix"] == prefix - - self.tpool = [] - if "tensor" in dic: - for t in dic["tensor"]: - id = t["id"] - if t["type"] == "initial": - self.tpool.append( - [np.load(f"{prefix}-tensor-{id}.npy"), "initial", None] - ) - assert list(self.tpool[-1][0].shape) == t["shape"] - elif t["type"] == "product": - if "subscripts" not in dic: - t["subscripts"] = None - self.tpool.append( - [ - None, - "product", - t["from 0"], - t["subscripts 0"], - t["from 1"], - t["subscripts 1"], - t["subscripts"], - ] - ) - elif "squeezer L" in dic or "squeezer R" in dic: - self.tpool.append([None, t["type"], t["from"]]) - self.gpool = [] - if "generator" in dic: - for g in dic["generator"]: - id = g["id"] - data = np.load(f"{prefix}-generator-{id}.npz") - self.gpool.append([data["s"], g["dim to"], data["wl"], data["wr"]]) diff --git a/test/mps/test_apply.py b/test/mps/test_apply.py index 057d9f1..2b86dd4 100644 --- a/test/mps/test_apply.py +++ b/test/mps/test_apply.py @@ -4,16 +4,14 @@ from pytest import approx -def apply(op, m0, m1, m2, m3, m4, m5, v, seq, pos, maxdim=None): +def apply(op, m0, m1, m2, m3, v, seq, pos, maxdim=None): m0 = q.mps.apply(m0, op, pos) m1 = q.mps.apply(m1, op, pos) m2 = q.mps.apply(m2, op, pos) m3 = q.mps.apply(m3, op, pos, maxdim) - m4 = q.mps.apply(m4, op, pos, maxdim) - m5 = q.mps.apply(m5, op, pos, maxdim) v = q.sv.apply(v, op, pos) seq.append([op, pos]) - return m0, m1, m2, m3, m4, m5, v, seq + return m0, m1, m2, m3, v, seq def test_apply(): @@ -23,19 +21,15 @@ def test_apply(): m0 = q.mps.zero(n, q.mps.MPS_C) m1 = q.mps.zero(n, q.mps_p.MPS_P) - m2 = q.mps.zero(n, q.mps_t.MPS_T) - m3 = q.mps.zero(n, q.mps.MPS_C) - m4 = q.mps.zero(n, q.mps_p.MPS_P) - m5 = q.mps.zero(n, q.mps_t.MPS_T) + m2 = q.mps.zero(n, q.mps.MPS_C) + m3 = q.mps.zero(n, q.mps_p.MPS_P) v = q.sv.zero(n) seq = [] i = 4 j = 0 print("apply cz on {} and {}".format(i, j)) - m0, m1, m2, m3, m4, m5, v, seq = apply( - q.op.cz(), m0, m1, m2, m3, m4, m5, v, seq, [i, j], maxdim - ) + m0, m1, m2, m3, v, seq = apply(q.op.cz(), m0, m1, m2, m3, v, seq, [i, j], maxdim) for _ in range(p): i = random.randrange(n) @@ -44,44 +38,40 @@ def test_apply(): t = random.randrange(3) if t == 0: print("apply h on {}".format(i)) - m0, m1, m2, m3, m4, m5, v, seq = apply( - q.op.h(), m0, m1, m2, m3, m4, m5, v, seq, [i], maxdim + m0, m1, m2, m3, v, seq = apply( + q.op.h(), m0, m1, m2, m3, v, seq, [i], maxdim ) elif t == 1: print("apply x on {}".format(i)) - m0, m1, m2, m3, m4, m5, v, seq = apply( - q.op.s(), m0, m1, m2, m3, m4, m5, v, seq, [i], maxdim + m0, m1, m2, m3, v, seq = apply( + q.op.s(), m0, m1, m2, m3, v, seq, [i], maxdim ) elif t == 2: print("apply s on {}".format(i)) - m0, m1, m2, m3, m4, m5, v, seq = apply( - q.op.t(), m0, m1, m2, m3, m4, m5, v, seq, [i], maxdim + m0, m1, m2, m3, v, seq = apply( + q.op.t(), m0, m1, m2, m3, v, seq, [i], maxdim ) else: t = random.randrange(2) if t == 0: print("apply cx on {} and {}".format(i, j)) - m0, m1, m2, m3, m4, m5, v, seq = apply( - q.op.cx(), m0, m1, m2, m3, m4, m5, v, seq, [i, j], maxdim + m0, m1, m2, m3, v, seq = apply( + q.op.cx(), m0, m1, m2, m3, v, seq, [i, j], maxdim ) elif t == 1: print("apply cz on {} and {}".format(i, j)) - m0, m1, m2, m3, m4, m5, v, seq = apply( - q.op.cz(), m0, m1, m2, m3, m4, m5, v, seq, [i, j], maxdim + m0, m1, m2, m3, v, seq = apply( + q.op.cz(), m0, m1, m2, m3, v, seq, [i, j], maxdim ) m0._is_canonical() m1._is_canonical() m2._is_canonical() m3._is_canonical() - m4._is_canonical() - m5._is_canonical() f0 = q.sv.fidelity(q.mps.state_vector(m0), v) f1 = q.sv.fidelity(q.mps.state_vector(m1), v) f2 = q.sv.fidelity(q.mps.state_vector(m2), v) f3 = q.sv.fidelity(q.mps.state_vector(m3), v) - f4 = q.sv.fidelity(q.mps.state_vector(m4), v) - f5 = q.sv.fidelity(q.mps.state_vector(m5), v) - print("fidelity = {} {} {} {} {} {}".format(f0, f1, f2, f3, f4, f5)) + print("fidelity = {} {} {} {}".format(f0, f1, f2, f3)) assert f0 == approx(1) assert f1 == approx(1) assert f2 == approx(1) @@ -90,21 +80,15 @@ def test_apply(): f1 = q.sv.fidelity(q.mps.state_vector(m1), v) f2 = q.sv.fidelity(q.mps.state_vector(m2), v) f3 = q.sv.fidelity(q.mps.state_vector(m3), v) - f4 = q.sv.fidelity(q.mps.state_vector(m4), v) - f5 = q.sv.fidelity(q.mps.state_vector(m5), v) - print("final fidelity = {} {} {} {} {} {}".format(f0, f1, f2, f3, f4, f5)) + print("final fidelity = {} {} {} {}".format(f0, f1, f2, f3)) assert f0 == approx(1) assert f1 == approx(1) assert f2 == approx(1) # assert f2 == approx(f3) - m6 = q.mps.zero(n, q.mps_p.MPS_P) - f6 = q.sv.fidelity(q.mps.state_vector(q.mps.apply_seq(m6, seq)), v) - assert f6 == approx(1) - - m7 = q.mps.zero(n, q.mps_t.MPS_T) - f7 = q.sv.fidelity(q.mps.state_vector(q.mps.apply_seq(m7, seq)), v) - assert f7 == approx(1) + m4 = q.mps.zero(n, q.mps_p.MPS_P) + f4 = q.sv.fidelity(q.mps.state_vector(q.mps.apply_seq(m4, seq)), v) + assert f4 == approx(1) if __name__ == "__main__": diff --git a/test/mps/test_canonical.py b/test/mps/test_canonical.py index d8053d1..cf39395 100644 --- a/test/mps/test_canonical.py +++ b/test/mps/test_canonical.py @@ -14,7 +14,7 @@ def test_canonical(): tensors.append(np.random.random((d, 2, dn))) d = dn tensors.append(np.random.random((d, 2, 1))) - for mps in [q.mps.MPS_C, q.mps_p.MPS_P, q.mps_t.MPS_T]: + for mps in [q.mps.MPS_C, q.mps_p.MPS_P]: m = mps(tensors) norm = q.mps.norm(m) assert q.mps.norm(m) == approx(norm) @@ -36,7 +36,7 @@ def test_canonical(): v = np.random.random(2**n).reshape((2,) * n + (1,)) v /= np.linalg.norm(v) tensors = q.mps.tensor_decomposition(v, maxdim) - for mps in [q.mps.MPS_C, q.mps_p.MPS_P, q.mps_t.MPS_T]: + for mps in [q.mps.MPS_C, q.mps_p.MPS_P]: m = mps(tensors) norm = q.mps.norm(m) diff --git a/test/mps/test_move.py b/test/mps/test_move.py index 9d8a3bd..1c1f25a 100644 --- a/test/mps/test_move.py +++ b/test/mps/test_move.py @@ -16,7 +16,7 @@ def test_swap(): tensors.append(np.random.random((d, 2, dn))) d = dn tensors.append(np.random.random((d, 2, 1))) - for mps in [q.mps.MPS_C, q.mps_p.MPS_P, q.mps_t.MPS_T]: + for mps in [q.mps.MPS_C, q.mps_p.MPS_P]: m = mps(tensors) q.mps.is_canonical(m) norm = q.mps.norm(m) @@ -48,7 +48,7 @@ def test_move(): tensors.append(np.random.random((d, 2, dn))) d = dn tensors.append(np.random.random((d, 2, 1))) - for mps in [q.mps.MPS_C, q.mps_p.MPS_P, q.mps_t.MPS_T]: + for mps in [q.mps.MPS_C, q.mps_p.MPS_P]: m = mps(tensors) q.mps.is_canonical(m) norm = q.mps.norm(m) diff --git a/test/mps/test_mps.py b/test/mps/test_mps.py index 70888f6..22a7e8f 100644 --- a/test/mps/test_mps.py +++ b/test/mps/test_mps.py @@ -7,13 +7,10 @@ def test_mps(): v = q.sv.product_state(states) m0 = q.mps.product_state(states) m1 = q.mps_p.product_state(states) - m2 = q.mps_t.product_state(states) v0 = q.mps.state_vector(m0) v1 = q.mps.state_vector(m1) - v2 = q.mps.state_vector(m2) assert np.allclose(v0, v) assert np.allclose(v1, v) - assert np.allclose(v2, v) if __name__ == "__main__": diff --git a/test/mps_p/test_projector.py b/test/mps_p/test_projector.py index 2df2121..7f54568 100644 --- a/test/mps_p/test_projector.py +++ b/test/mps_p/test_projector.py @@ -1,7 +1,5 @@ import numpy as np import qailo as q -from pytest import approx -from qailo.mps_p.projector import _full_projector def test_projector(): @@ -11,7 +9,7 @@ def test_projector(): n0, n1, n2, d = np.random.randint(2, maxn, size=(4,)) T0 = np.random.random((n0, n1)) T1 = np.random.random((n1, n2)) - _, WLh, WR = q.mps_p.compact_projector(T0, [0, 2], T1, [2, 1], nkeep=d) + _, WLh, WR = q.mps_p.projector(T0, [0, 2], T1, [2, 1], nkeep=d) assert WLh.shape[1] <= d and WR.shape[1] <= d assert WLh.shape[0] == n1 and WR.shape[0] == n1 A = np.einsum(T0, [0, 2], WR, [2, 3], WLh.T, [3, 4], T1, [4, 1]) @@ -24,9 +22,7 @@ def test_projector(): n0, n1, n2, n3, n4, n5, d = np.random.randint(2, maxn, size=(7,)) T0 = np.random.random((n0, n1, n2, n3)) T1 = np.random.random((n3, n2, n4, n5)) - _, WLh, WR = q.mps_p.compact_projector( - T0, [0, 1, 4, 5], T1, [5, 4, 2, 3], nkeep=d - ) + _, WLh, WR = q.mps_p.projector(T0, [0, 1, 4, 5], T1, [5, 4, 2, 3], nkeep=d) assert WLh.shape[2] <= d and WR.shape[2] <= d assert WLh.shape[0] == n2 and WLh.shape[1] == n3 assert WR.shape[0] == n2 and WR.shape[1] == n3 @@ -52,7 +48,7 @@ def test_projector(): assert np.allclose(np.einsum(t0, [0, 1, 4, 5], t1, [5, 4, 2, 3]), B) tt0 = np.einsum(w0, [0, 4], t0, [4, 1, 2, 3]) tt1 = np.einsum(t1, [0, 1, 2, 4], w2, [4, 3]) - _, WLh, WR = q.mps_p.compact_projector(tt0, [0, 1, 4, 5], tt1, [5, 4, 2, 3]) + _, WLh, WR = q.mps_p.projector(tt0, [0, 1, 4, 5], tt1, [5, 4, 2, 3]) assert np.allclose( np.einsum(WLh, [2, 3, 0], WR, [2, 3, 1]), np.identity(WLh.shape[2]) ) @@ -76,7 +72,7 @@ def test_projector(): assert np.allclose(np.einsum(t0, [0, 1, 4, 5], t1, [4, 5, 2, 3]), B) tt0 = np.einsum(w0, [0, 4], t0, [4, 1, 2, 3]) tt1 = np.einsum(t1, [0, 1, 2, 4], w2, [4, 3]) - _, WLh, WR = q.mps_p.compact_projector(tt0, [0, 1, 4, 5], tt1, [4, 5, 2, 3]) + _, WLh, WR = q.mps_p.projector(tt0, [0, 1, 4, 5], tt1, [4, 5, 2, 3]) assert np.allclose( np.einsum(WLh, [2, 3, 0], WR, [2, 3, 1]), np.identity(WLh.shape[2]) ) @@ -87,107 +83,5 @@ def test_projector(): assert np.allclose(A, B) -def test_full_projector(): - maxn = 4 - nt = 16 - for _ in range(nt): - n0, n1, n2, d = np.random.randint(2, maxn, size=(4,)) - d = min(d, n0, n2) - T0 = np.random.random((n0, n1)) - T1 = np.random.random((n1, n2)) - S, _, WLh, WR = _full_projector(T0, [0, 2], T1, [2, 1]) - assert WLh.shape[0] == n1 and WR.shape[0] == n1 - A = np.einsum(T0, [0, 2], WR[:, :d], [2, 3], WLh.T[:d, :], [3, 4], T1, [4, 1]) - - S, U, V = q.mps.compact_svd(np.einsum(T0, [0, 2], T1, [2, 1]), nkeep=d) - B = np.einsum("k,ik,jk->ij", S, U, V) - assert np.allclose(A, B) - - for _ in range(nt): - n0, n1, n2, n3, n4, n5, d = np.random.randint(2, maxn, size=(7,)) - d = min(d, n0 * n1, n4 * n5) - T0 = np.random.random((n0, n1, n2, n3)) - T1 = np.random.random((n3, n2, n4, n5)) - S, WLh, WR = q.mps_p.compact_projector( - T0, [0, 1, 4, 5], T1, [5, 4, 2, 3], nkeep=d - ) - assert WLh.shape[0] == n2 and WLh.shape[1] == n3 - assert WR.shape[0] == n2 and WR.shape[1] == n3 - A = np.einsum( - T0, - [0, 1, 4, 5], - WR[:, :, :d], - [4, 5, 6], - WLh[:, :, :d], - [7, 8, 6], - T1, - [8, 7, 2, 3], - ) - - L, R = q.mps.tensor_svd( - np.einsum(T0, [0, 1, 4, 5], T1, [5, 4, 2, 3]), [[0, 1], [2, 3]], nkeep=d - ) - B = np.einsum(L, [0, 1, 4], R, [4, 2, 3]) - assert np.allclose(A, B) - - for _ in range(nt): - w0 = np.random.random((1, 1)) + 1.0j * np.random.random((1, 1)) - t0 = np.random.random((1, 2, 3)) + 1.0j * np.random.random((1, 2, 3)) - t1 = np.random.random((3, 2, 2)) + 1.0j * np.random.random((3, 2, 2)) - w2 = np.random.random((2, 2)) + 1.0j * np.random.random((2, 2)) - p = np.random.random((2, 2, 2, 2)) + 1.0j * np.random.random((2, 2, 2, 2)) - B = np.einsum(t0, [0, 4, 6], t1, [6, 5, 3], p, [1, 2, 4, 5]) - p0, p1 = q.mps.svd.tensor_svd(p, [[0, 2], [1, 3]]) - assert np.allclose(np.einsum(p0, [0, 2, 4], p1, [4, 1, 3]), p) - t0 = np.einsum(t0, [0, 4, 3], p0, [1, 4, 2]) - t1 = np.einsum(t1, [0, 4, 3], p1, [1, 2, 4]) - assert np.allclose(np.einsum(t0, [0, 1, 4, 5], t1, [5, 4, 2, 3]), B) - tt0 = np.einsum(w0, [0, 4], t0, [4, 1, 2, 3]) - tt1 = np.einsum(t1, [0, 1, 2, 4], w2, [4, 3]) - _, _, WLh, WR = _full_projector(tt0, [0, 1, 4, 5], tt1, [5, 4, 2, 3]) - tt0 = np.einsum(t0, [0, 1, 3, 4], WR, [3, 4, 2]) - tt1 = np.einsum(WLh, [3, 4, 0], t1, [4, 3, 1, 2]) - A = np.einsum(tt0, [0, 1, 4], tt1, [4, 2, 3]) - print(np.linalg.norm(A - B)) - assert np.allclose(A, B) - - assert np.allclose( - np.einsum(WLh, [2, 3, 0], WR, [2, 3, 1]), np.identity(WLh.shape[2]) - ) - n = WR.shape[0] * WR.shape[1] - assert np.linalg.norm( - np.einsum(WR, [0, 1, 4], WLh, [2, 3, 4]).reshape((n, n)) - np.identity(n) - ) == approx(0, abs=1e-8) - - for _ in range(nt): - w0 = np.random.random((1, 1)) + 1.0j * np.random.random((1, 1)) - t0 = np.random.random((1, 2, 3)) + 1.0j * np.random.random((1, 2, 3)) - t1 = np.random.random((3, 2, 2)) + 1.0j * np.random.random((3, 2, 2)) - w2 = np.random.random((2, 2)) + 1.0j * np.random.random((2, 2)) - p = np.random.random((2, 2, 2, 2)) + 1.0j * np.random.random((2, 2, 2, 2)) - B = np.einsum(t0, [0, 4, 6], t1, [6, 5, 3], p, [1, 2, 4, 5]) - p0, p1 = q.mps.svd.tensor_svd(p, [[0, 2], [1, 3]]) - assert np.allclose(np.einsum(p0, [0, 2, 4], p1, [4, 1, 3]), p) - t0 = np.einsum(t0, [0, 4, 3], p0, [1, 4, 2]) - t1 = np.einsum(t1, [1, 4, 3], p1, [0, 2, 4]) - assert np.allclose(np.einsum(t0, [0, 1, 4, 5], t1, [4, 5, 2, 3]), B) - tt0 = np.einsum(w0, [0, 4], t0, [4, 1, 2, 3]) - tt1 = np.einsum(t1, [0, 1, 2, 4], w2, [4, 3]) - _, _, WLh, WR = _full_projector(tt0, [0, 1, 4, 5], tt1, [4, 5, 2, 3]) - tt0 = np.einsum(t0, [0, 1, 3, 4], WR, [3, 4, 2]) - tt1 = np.einsum(WLh, [3, 4, 0], t1, [3, 4, 1, 2]) - A = np.einsum(tt0, [0, 1, 4], tt1, [4, 2, 3]) - print(np.linalg.norm(A - B)) - assert np.allclose(A, B) - - assert np.allclose( - np.einsum(WLh, [2, 3, 0], WR, [2, 3, 1]), np.identity(WLh.shape[2]) - ) - n = WR.shape[0] * WR.shape[1] - assert np.linalg.norm( - np.einsum(WR, [0, 1, 4], WLh, [2, 3, 4]).reshape((n, n)) - np.identity(n) - ) == approx(0, abs=1e-8) - - if __name__ == "__main__": test_projector() diff --git a/test/mps_t/test_mps_t.py b/test/mps_t/test_mps_t.py deleted file mode 100644 index 112bebd..0000000 --- a/test/mps_t/test_mps_t.py +++ /dev/null @@ -1,40 +0,0 @@ -import qailo as q -from qailo.mps_t.tpool import tpool - - -def test_mps_t(): - v = q.mps_t.zero(3) - - print("input:") - print("state vector:", q.vector(v)) - print("probabitily:", q.probability(v)) - - v = q.apply(v, q.op.h(), [0]) - v = q.apply(v, q.op.h(), [2]) - v = q.apply(v, q.op.cx(), [0, 1]) - v = q.apply(v, q.op.cz(), [1, 2]) - v = q.apply(v, q.op.h(), [2]) - - print("output:") - print("state vector:", q.vector(v)) - print("probabitily:", q.probability(v)) - - print("# tensor pool") - for id, tp in enumerate(v.tp.tpool): - print(f"{id} {tp[0].shape} {tp[1]} {tp[2]}") - assert len(v.tp.tpool) == 25 - - print("# generator pool") - for id, gp in enumerate(v.tp.gpool): - print(f"{id} {gp[0].shape} {gp[1]} {gp[2].shape} {gp[3].shape}") - assert len(v.tp.gpool) == 2 - - prefix = "test_mps_t" - v.tp._dump(prefix) - - tp = tpool() - tp._load(prefix) - - -if __name__ == "__main__": - test_mps_t()