light experiment works

06_krisis.py

@@ -2,7 +2,7 @@ import cirq
 import numpy as np
 from collections import defaultdict
 
-from lib_q_computer_math import State, QuantumCircuit, QuantumProcessor, C, H, x, _, _0, _1
+from lib import State, QuantumCircuit, QuantumProcessor, C, H, x, _, _0, _1
 
 
 def from_angles_1(theta, phi):

compiler.py

@@ -1,4 +1,4 @@
-from lib_q_computer_math import s
+from lib import s
 
 
 def int_to_state(i):

lib.py

@@ -2,6 +2,7 @@ import random
 from collections import defaultdict
 from functools import reduce
 from numbers import Complex
+from pprint import pprint
 from typing import Union
 
 import numpy as np
@@ -46,6 +47,9 @@ class Matrix(object):
     def __add__(self, other):
         return Matrix(self.m + other.m)
 
+    def __sub__(self, other):
+        return Matrix(self.m - other.m)
+
     def __eq__(self, other):
         if isinstance(other, Complex):
             return bool(np.allclose(self.m, other))
@@ -86,8 +90,30 @@ class Matrix(object):
     def __len__(self):
         return len(self.m)
 
+    def inner(self, other):
+        """Define inner product - a.k.a dot product, scalar product - <0|0>"""
+        return Matrix(np.dot(self._conjugate_transpose(), other.m))
+
+    def dot(self, other):
+        """Alias to inner"""
+        return self.inner(other)
+
+    def scalar(self, other):
+        """Alias to inner"""
+        return self.inner(other)
+
     def outer(self, other):
-        """Define outer product |0><0|"""
+        """Define outer product |0><0|
+        https://en.wikipedia.org/wiki/Outer_product
+
+        Given two vectors u and v, their outer product u ⊗ v is defined
+        as the m × n matrix A obtained by multiplying each element of u
+        by each element of v
+
+        The outer product u ⊗ v is equivalent to a matrix multiplication uvT,
+        provided that u is represented as a m × 1 column vector and v as a
+        n × 1 column vector (which makes vT a row vector).
+        """
         return Matrix(np.outer(self.m, other.m))
 
     def x(self, other):
@@ -136,10 +162,16 @@ class Vector(Matrix):
         return self.m.shape[1] == 1
 
 
+VectorOrList = Union[Vector, ListOrNdarray]
+
+
 class State(Vector):
-    def __init__(self, m: ListOrNdarray = None, name: str = '', *args, **kwargs):
+    def __init__(self, m: VectorOrList = None, name: str = '', *args, **kwargs):
         """State vector representing quantum state"""
-        super().__init__(m, *args, **kwargs)
+        if type(m) is Vector:
+            super().__init__(m.m, *args, **kwargs)
+        else:
+            super().__init__(m, *args, **kwargs)
         self.name = name
         self.measurement_result = None
         # TODO: SHOULD WE NORMALIZE?
@@ -190,6 +222,9 @@ class State(Vector):
         assert 0 <= phi <= 2 * np.pi
         return theta, phi
 
+    def to_ket(self):
+        return self.conjugate_transpose()
+
     def rotate_x(self, theta):
         return Rx(theta).on(self)
 
@@ -276,7 +311,7 @@ class State(Vector):
         Measures with a measurement operator mo
         TODO: Can't define `mo: MeasurementOperator` because in python you can't declare classes before defining them
         """
-        m = mo.on(self) / np.sqrt(mo.get_prob(self))
+        m = mo.on(self).m / np.sqrt(mo.get_prob(self))
         return State(m)
 
     def measure_partial(self, qubit_n):
@@ -342,7 +377,7 @@ class State(Vector):
 
 
 def test_measure_partial():
-    state = b_00
+    state = b_phi_p
     state.measure_partial(1)
 
 
@@ -635,31 +670,31 @@ _m = State([[1 / np.sqrt(2)],
            name='-')
 
 # 4 Bell states
-b_00 = State([[1 / np.sqrt(2)],
+b_phi_p = State([[1 / np.sqrt(2)],
-              [0],
+                 [0],
-              [0],
+                 [0],
-              [1 / np.sqrt(2)]],
+                 [1 / np.sqrt(2)]],
-             name='B00')
+                name="{}+".format(REPR_GREEK_PHI))
 
-b_01 = State([[0],
+b_psi_p = State([[0],
-              [1 / np.sqrt(2)],
+                 [1 / np.sqrt(2)],
-              [1 / np.sqrt(2)],
+                 [1 / np.sqrt(2)],
-              [0]],
+                 [0]],
-             name='B01')
+                name="{}+".format(REPR_GREEK_PSI))
 
-b_10 = State([[1 / np.sqrt(2)],
+b_phi_m = State([[1 / np.sqrt(2)],
-              [0],
+                 [0],
-              [0],
+                 [0],
-              [-1 / np.sqrt(2)]],
+                 [-1 / np.sqrt(2)]],
-             name='B10')
+                name="{}-".format(REPR_GREEK_PHI))
 
-b_11 = State([[0],
+b_psi_m = State([[0],
-              [1 / np.sqrt(2)],
+                 [1 / np.sqrt(2)],
-              [-1 / np.sqrt(2)],
+                 [-1 / np.sqrt(2)],
-              [0]],
+                 [0]],
-             name='B11')
+                name="{}-".format(REPR_GREEK_PSI))
 
-well_known_states = [_p, _m, b_00, b_01, b_10, b_11]
+well_known_states = [_p, _m, b_phi_p, b_psi_p, b_phi_m, b_psi_m]
 
 _ = I = UnitaryOperator([[1, 0],
                          [0, 1]],
@@ -847,17 +882,17 @@ def test():
     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.inner(_0) == 1  # parallel have 1 product
-    assert _0 | _1 == 0  # orthogonal have 0 product
+    assert _0.inner(_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 _1.inner(8 * _0) == 8 * _1.inner(_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 _p.inner(_1 + _0) == _p.inner(_1) + _p.inner(_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
+    assert _0.inner(_1) == _1.inner(_0).complex_conjugate()  # non-commutative inner product
 
     test_to_from_angles()
 
@@ -871,14 +906,14 @@ def test():
     test_unitary_hermitian()
 
     # Pauli X gate flips the |0> to |1> and the |1> to |0>
-    assert X | _1 == _0
+    assert X.on(_1) == _0
-    assert X | _0 == _1
+    assert X.on(_0) == _1
 
     # Test the Y Pauli operator with complex number literal notation
-    assert Y | _0 == State([[0],
+    assert Y.on(_0) == State([[0],
-                            [1j]])
+                              [1j]])
-    assert Y | _1 == State([[-1j],
+    assert Y.on(_1) == State([[-1j],
-                            [0]])
+                              [0]])
 
     # Test Pauli rotation gates
     testRotPauli()
@@ -928,11 +963,11 @@ def test():
     # 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
+    superpos = H.on(_0)
     assert superpos == _p
 
     # Then CNOT the superposition with a |0> qubit
-    bell = CNOT | (superpos * _0)
+    bell = CNOT.on(superpos * _0)
     assert bell == State([[1 / np.sqrt(2)],
                           [0.],
                           [0.],
@@ -944,7 +979,6 @@ def test():
     assert np.isclose(bell.get_prob(0b11), 0.5)  # Probability for bell in 11 is 0.5
 
     ################################
-    # TODO: Don't know where outer product fits - something about density operator?
     assert _0.x(_0) == Matrix([[1, 0],
                                [0, 0]])
     assert _0.x(_1) == Matrix([[0, 1],
@@ -1196,33 +1230,44 @@ def test_quantum_processor():
 
 def test_light():
     # http://alienryderflex.com/polarizer/
-    # TODO: Are these measurement operators the correct way to represent hor/ver/diag filter?
-    #       No, becuase they are not unitaries
     hor_filter = MeasurementOperator.create_from_prob(Matrix(_0.m), name='h')
     diag_filter = MeasurementOperator.create_from_prob(Matrix(_p.m), name='d')
     ver_filter = MeasurementOperator.create_from_prob(Matrix(_1.m), name='v')
 
-    # create random light/photon
+    def random_light():
-    random_pol = Vector([[np.random.uniform(0, 1)], [np.random.uniform(0, 1)]])
+        random_pol = Vector([[np.random.uniform(0, 1)], [np.random.uniform(0, 1)]])
-    random_light = normalize_state(random_pol)
+        return normalize_state(random_pol)
 
-    # TODO: Doesn't work...
+    def experiment(id, random_ls, filters):
-    qc = QuantumCircuit(1, initial_steps=[[random_light]])
+        total_light = defaultdict(int)
-    qc.add_row([hor_filter, ver_filter])
+        human_state = "Light"
-    qc.print()
+        for filter in filters:
-    qp = QuantumProcessor(qc)
+            human_state += "->{}".format(filter.name)
-    qp.print_sample(qp.get_sample(100))
 
-    # add three filters as operators
+        for r in random_ls:
-    qc = QuantumCircuit(1, initial_steps=[[random_light]])
+            st = filters[0].on(r)
-    qc.add_row([hor_filter, diag_filter, ver_filter])
+            for filter in filters:
-    qc.print()
+                st = filter.on(st)
-    qp = QuantumProcessor(qc)
+            st = State(st)
-    qp.print_sample(qp.get_sample(100))
+            total_light[st.measure()] += 1
+        print("{}. {}:\n    {}".format(id, human_state, dict(total_light)))
+
+    its = 100
+    random_lights = [random_light() for _ in range(its)]
+
+    # Just horizonal - should result in some light
+    experiment(1, random_lights, [hor_filter, ])
+
+    # Vertical after horizonal - should result in 0
+    experiment(2, random_lights, [ver_filter, hor_filter])
+
+    # TODO: Something is wrong here...
+    # Vertical after diagonal after horizontal - should result in ~ 50% compared to only Horizontal
+    experiment(3, random_lights, [ver_filter, diag_filter, hor_filter])
 
 
 if __name__ == "__main__":
     test()
-    # test_quantum_processor()
+    test_quantum_processor()
-    # test_light()
+    test_light()
     # test_measure_partial()

play.py

@@ -1,20 +0,0 @@
-import numpy as np
-from lib_q_computer_math import Vector, Matrix, MeasurementOperator, _0, _1, _p, normalize_state, State
-
-
-def main():
-    random_pol = Vector([[np.random.uniform(0, 1)], [np.random.uniform(0, 1)]])
-    random_light = normalize_state(random_pol)
-    hor_filter = MeasurementOperator.create_from_prob(Matrix(_0.m), name='h')
-    diag_filter = MeasurementOperator.create_from_prob(Matrix(_p.m), name='d')
-    ver_filter = MeasurementOperator.create_from_prob(Matrix(_1.m), name='v')
-    State(hor_filter.on(random_light).m).measure()
-    print(hor_filter.on(_0))
-    print(ver_filter.on(_0))
-    print(ver_filter.on(hor_filter.on(_0)))
-    print(ver_filter.on(diag_filter.on(hor_filter.on(_0))))
-    print()
-
-
-if __name__ == "__main__":
-    main()