From c8e03fc8b42f98712b1b308a18536ded540bff50 Mon Sep 17 00:00:00 2001 From: Daniel Tsvetkov Date: Tue, 17 Dec 2019 18:05:14 +0100 Subject: [PATCH] deep math understanding and following qm postulates --- lib_q_computer.py | 11 ++ lib_q_computer_math.py | 293 +++++++++++++++++++++++++++++++++++++++++ load_test.py | 2 +- 3 files changed, 305 insertions(+), 1 deletion(-) create mode 100644 lib_q_computer_math.py diff --git a/lib_q_computer.py b/lib_q_computer.py index a32fb33..1c11780 100644 --- a/lib_q_computer.py +++ b/lib_q_computer.py @@ -111,6 +111,17 @@ class State(object): choices = [format_str.format(i) for i in range(len(weights))] return random.choices(choices, weights)[0] + def partial_measure(self): + """ + Say we have the state a|00> + b|01> + c|10> + d|11>. + Measuring \0> on q1: + |0>(a|0> + b|1>) + |1>(c|0> + d|1>) -> + + We need to normalize for the other qubit: (a|0> + b|1>) / math.sqrt(|a|^2 + \b|^2) + :return: + """ + ... + @jit() def kron(_list): diff --git a/lib_q_computer_math.py b/lib_q_computer_math.py new file mode 100644 index 0000000..d6efe4f --- /dev/null +++ b/lib_q_computer_math.py @@ -0,0 +1,293 @@ +from __future__ import annotations + +from typing import Union + +import numpy as np +from numbers import Complex + +from load_test import sizeof_fmt + +ListOrNdarray = Union[list, np.ndarray] + + +class Matrix(object): + """Represents a quantum state as a vector in Hilbert space""" + + def __init__(self, m: ListOrNdarray = None, *args, **kwargs): + """ + Can be initialized with a matrix, e.g. for |0> this is [[0],[1]] + :param m: a matrix representing the quantum state + :param name: the name of the state (optional) + """ + if m is None: + self.m = np.array([]) + elif type(m) is list: + self.m = np.array(m) + elif type(m) is np.ndarray: + self.m = m + else: + raise TypeError("m needs to be a list or ndarray type") + + def __add__(self, other): + return other.__class__(self.m + other.m) + + def __eq__(self, other): + if isinstance(other, Complex): + return bool(np.allclose(self.m, other)) + return np.allclose(self.m, other.m) + + def __mul_or_kron__(self, other): + """ + Multiplication with a number is Linear op; + with another state is a composition via the Kronecker/Tensor product + """ + if isinstance(other, Complex): + return self.__class__(self.m * other) + elif isinstance(other, Matrix): + return other.__class__(np.kron(self.m, other.m)) + raise NotImplementedError + + def __rmul__(self, other): + return self.__mul_or_kron__(other) + + def __mul__(self, other): + return self.__mul_or_kron__(other) + + def __or__(self, other): + """Define inner product: + """ + return other.__class__(np.dot(self._conjugate_transpose(), other.m)) + + def __len__(self): + return len(self.m) + + def _conjugate_transpose(self): + return self.m.transpose().conjugate() + + def _complex_conjugate(self): + return self.m.conjugate() + + def conjugate_transpose(self): + return State(self._conjugate_transpose()) + + def complex_conjugate(self): + return State(self._complex_conjugate()) + + +class State(Matrix): + def __init__(self, m: ListOrNdarray = None, name: str = '', *args, **kwargs): + """An Operator turns one function into another""" + super().__init__(m) + self.name = name + + def __repr__(self): + if self.name: + return '|{}>'.format(self.name) + return str(self.m) + + def norm(self): + """Norm/Length of the vector = sqrt()""" + return self.length() + + def length(self): + """Norm/Length of the vector = sqrt()""" + return np.sqrt((self | self).m).item(0) + + def is_orthogonal(self, other): + """If the inner (dot) product is zero, this vector is orthogonal to other""" + return self | other == 0 + + def measure(self, j): + """pr(j) = ||^2""" + # j_th basis vector + e_j = State([[1] if i == int(j) else [0] for i in range(len(self.m))]) + return np.absolute((e_j | self).m.item(0)) ** 2 + + +class Operator(object): + def __init__(self, func=None, *args, **kwargs): + """An Operator turns one function into another""" + self.func = func + + def on(self, *args, **kwargs): + return self(*args, **kwargs) + + def __call__(self, *args, **kwargs): + return self.func(*args, **kwargs) + + +class UnitaryMatrix(Matrix): + def __init__(self, m: ListOrNdarray, *args, **kwargs): + """Represents a Unitary matrix that satisfies UU+ = I""" + super().__init__(m, *args, **kwargs) + if not self._is_unitary(): + raise TypeError("Not a Unitary matrix") + + def _is_unitary(self): + """Checks if the matrix is + 1. square and + 2. the product of itself with conjugate transpose is Identity UU+ = I""" + is_square = self.m.shape[0] == self.m.shape[1] + UU_ = np.dot(self._conjugate_transpose(), self.m) + I = np.eye(self.m.shape[0]) + return is_square and np.isclose(UU_, I).all() + + +class UnitaryOperator(Operator, UnitaryMatrix): + def __init__(self, m: ListOrNdarray, name: str = '', *args, **kwargs): + """UnitaryOperator inherits from both Operator and a Unitary matrix + It is used to act on a State vector by defining the operator to be the dot product""" + UnitaryMatrix.__init__(self, m=m, name=name, *args, **kwargs) + Operator.__init__(self, func=lambda other: State(np.dot(self.m, other.m)), *args, **kwargs) + + +def test(): + _0 = State([[1], [0]], name='0') + _1 = State([[0], [1]], name='1') + _p = State([[1 / np.sqrt(2)], [1 / np.sqrt(2)]], name='+') + m = State([[1 / np.sqrt(2)], [-1 / np.sqrt(2)]], name='-') + + # Test properties of Hilbert vector space + # The four postulates of Quantum Mechanics + # I: States | Associated to any physical system is a complex vector space + # known as the state space of the system. If the system is closed + # then the system is described completely by its state vector + # which is a unit vector in the space. + # Mathematically, this vector space is also a function space + assert _0 + _1 == _1 + _0 # commutativity of vector addition + assert _0 + (_1 + _p) == (_0 + _1) + _p # associativity of vector addition + assert 8 * (_0 + _1) == 8 * _0 + 8 * _1 # Linear when multiplying by constants + assert _0 | _0 == 1 # parallel have 1 product + assert _0 | _1 == 0 # orthogonal have 0 product + assert _0.is_orthogonal(_1) + assert _1 | (8 * _0) == 8 * (_1 | _0) # Inner product is linear multiplied by constants + assert _p | (_1 + _0) == (_p | _1) + (_p | _0) # Inner product is linear in superpos of vectors + + assert np.isclose(_1.length(), 1.0) # all of the vector lengths are normalized + assert np.isclose(_0.length(), 1.0) + assert np.isclose(_p.length(), 1.0) + + assert _0 | _1 == (_1 | _0).complex_conjugate() # non-commutative inner product + + # II: Dynamics | The evolution of a closed system is described by a unitary transformation + # + # Operators turn one vector into another + # Pauli X gate flips the |0> to |1> and the |1> to |0> + + # the times 2 operator should return the times two multiplication + _times_2 = Operator(lambda x: 2 * x) + assert _times_2.on(5) == 10 + assert _times_2(5) == 10 + + X = UnitaryOperator([[0, 1], + [1, 0]]) + assert X | _1 == _0 + assert X | _0 == _1 + + # Test the Y Pauli operator with complex number literal notation + Y = UnitaryOperator([[0, -1j], + [1j, 0], ]) + + assert Y | _0 == State([[0], + [1j]]) + assert Y | _1 == State([[-1j], + [0]]) + + # III: Measurement | A quantum measurement is described by an orthonormal basis |e_j> + # for state space. If the initial state of the system is |psi> + # then we get outcome j with probability pr(j) = ||^2 + + assert _0.measure(0) == 1 # Probability for |0> in 0 is 1 + assert _0.measure(1) == 0 # Probability for |0> in 1 is 0 + + assert _1.measure(0) == 0 # Probability for |1> in 0 is 0 + assert _1.measure(1) == 1 # Probability for |1> in 1 is 1 + + assert np.isclose(_p.measure(0), 0.5) # Probability for |+> in 0 is 0.5 + assert np.isclose(_p.measure(1), 0.5) # Probability for |+> in 1 is 0.5 + + # IV: Compositing | tensor/kronecker product when composing + assert _0 * _0 == State([[1], [0], [0], [0]]) + assert _0 * _1 == State([[0], [1], [0], [0]]) + assert _1 * _0 == State([[0], [0], [1], [0]]) + assert _1 * _1 == State([[0], [0], [0], [1]]) + + # CNOT applies a control qubit on a target. + # If the control is a |0>, target remains unchanged. + # If the control is a |1>, target is flipped. + CNOT = UnitaryOperator([[1, 0, 0, 0], + [0, 1, 0, 0], + [0, 0, 0, 1], + [0, 0, 1, 0], ]) + assert CNOT.on(_0 * _0) == _0 * _0 + assert CNOT.on(_0 * _1) == _0 * _1 + assert CNOT.on(_1 * _0) == _1 * _1 + assert CNOT.on(_1 * _1) == _1 * _0 + + # ALL FOUR NOW - Create a Bell state (I) that has H|0> (II), measures (III) and composition in CNOT (IV) + # Bell state + # First - create a superposition + H = UnitaryOperator([[1 / np.sqrt(2), 1 / np.sqrt(2)], + [1 / np.sqrt(2), -1 / np.sqrt(2)], ]) + superpos = H | _0 + assert superpos == _p + + # Then CNOT the superposition with a |0> qubit + bell = CNOT | (superpos * _0) + assert bell == State([[1 / np.sqrt(2)], + [0.], + [0.], + [1 / np.sqrt(2)], ]) + + assert np.isclose(bell.measure(0b00), 0.5) # Probability for bell in 00 is 0.5 + assert np.isclose(bell.measure(0b01), 0) # Probability for bell in 01 is 0. + assert np.isclose(bell.measure(0b10), 0) # Probability for bell in 10 is 0. + assert np.isclose(bell.measure(0b11), 0.5) # Probability for bell in 11 is 0.5 + + +def naive_load_test(N): + import os + import psutil + import gc + from time import time + from sys import getsizeof + + print("{:>10} {:>10} {:>10} {:>10} {:>10} {:>10} {:>10} {:>10} {:>10}".format( + "qbits", + "kron_len", + "mem_used", + "mem_per_q", + "getsizeof", + "getsiz/len", + "nbytes", + "nbytes/len", + "time")) + + _0 = State([[1], [0]], name='0') + + process = psutil.Process(os.getpid()) + mem_init = process.memory_info().rss + for i in range(2, N + 1): + start = time() + m = _0 + for _ in range(i): + m = m * _0 + len_m = len(m) + elapsed = time() - start + mem_b = process.memory_info().rss - mem_init + print("{:>10} {:>10} {:>10} {:>10} {:>10} {:>10} {:>10} {:>10} {:>10}".format( + i, + len(m), + sizeof_fmt(mem_b), + sizeof_fmt(mem_b / len_m), + sizeof_fmt(getsizeof(m)), + sizeof_fmt(getsizeof(m) / len_m), + sizeof_fmt(m.m.nbytes), + sizeof_fmt(m.m.nbytes / len_m), + np.round(elapsed, 2))) + + gc.collect() + + +if __name__ == "__main__": + naive_load_test(23) diff --git a/load_test.py b/load_test.py index 2b7e235..4d3140c 100644 --- a/load_test.py +++ b/load_test.py @@ -25,7 +25,7 @@ def run_load_test(N=20): (2**30 * 8 ) / 1024**3 = 8.0 - Example run output with N=27 on 64GB machine **(without jit)** + Example run output with N=27 on 64GB machine **(without jit, without cython optimizations)** qbits kron_len mem_used mem_per_q getsizeof getsiz/len nbytes nbytes/len time 2 4 0.0B 0.0B 144.0B 36.0B 32.0B 8.0B 0.0