linear operators

lib_q_computer_math.py

@@ -113,6 +113,18 @@ class Matrix(object):
         return str(self.m)
 
 
+class HorizontalVector(Matrix):
+    """Horizontal vector is basically a list"""
+
+    def __init__(self, m: ListOrNdarray = None, *args, **kwargs):
+        super().__init__(m, *args, **kwargs)
+        if not self._is_vector():
+            raise TypeError("Not a vector")
+
+    def _is_vector(self):
+        return len(self.m.shape) == 1
+
+
 class Vector(Matrix):
     def __init__(self, m: ListOrNdarray = None, *args, **kwargs):
         super().__init__(m, *args, **kwargs)
@@ -337,38 +349,48 @@ def humanize(m):
 
 
 class Operator(object):
-    def __init__(self, func=None, *args, **kwargs):
+    def __init__(self, func: sp.Lambda, *args, **kwargs):
         """An Operator turns one function into another"""
         self.func = func
 
-    def on(self, *args, **kwargs):
+    def on(self, *args):
-        return self(*args, **kwargs)
+        return self(*args)
 
-    def __call__(self, *args, **kwargs):
+    def __call__(self, *args):
-        return self.func(*args, **kwargs)
+        return self.func(*args)
 
 
 class LinearOperator(Operator):
-    def __init__(self, func=None, *args, **kwargs):
+    def __init__(self, func: sp.Lambda, *args, **kwargs):
         """Linear operators satisfy f(x+y) = f(x) + f(y) and a*f(x) = f(a*x)"""
         super().__init__(func, *args, **kwargs)
         if not self._is_linear():
             raise TypeError("Not a linear operator")
 
     def _is_linear(self):
-        # TODO: How to verify if the func is linear?
+        # A function is (jointly) linear in a given set of variables
-        #       in case of Unitary Operator, self.func is a lambda that takes a Matrix (assumes has .m component)
+        # if all second-order derivatives are identically zero (including mixed ones).
+        # https://stackoverflow.com/questions/36283548/check-if-an-equation-is-linear-for-a-specific-set-of-variables
+        expr, vars_ = self.func.expr, self.func.variables
+        for x in vars_:
+            for y in vars_:
+                try:
+                    if not sp.Eq(sp.diff(expr, x, y), 0):
+                        return False
+                except TypeError:
+                    return False
         return True
-        # a, b = sympy.symbols('a, b')
+
-        # expr, vars_ = a+b, [a, b]
+
-        # for x in vars_:
+def test_linear_operator():
-        #     for y in vars_:
+    # Operators turn one vector into another
-        #         try:
+    # the times 2 operator should return the times two multiplication
-        #             if not sympy.Eq(sympy.diff(expr, x, y), 0):
+    _x = sp.Symbol('x')
-        #                 return False
+    _times_2 = LinearOperator(sp.Lambda(_x, 2 * _x))
-        #         except TypeError:
+    assert _times_2.on(5) == 10
-        #             return False
+    assert _times_2(5) == 10
-        # return True
+
+    assert_raises(TypeError, "Not a linear operator", LinearOperator, sp.Lambda(_x, _x ** 2))
 
 
 class SquareMatrix(Matrix):
@@ -381,6 +403,39 @@ class SquareMatrix(Matrix):
         return self.m.shape[0] == self.m.shape[1]
 
 
+class LinearTransformation(LinearOperator, Matrix):
+    def __init__(self, m: ListOrNdarray, func=None, name: str = '', *args, **kwargs):
+        """LinearTransformation (or linear map) inherits from both LinearOperator and a Matrix
+        It is used to act on a State vector by defining the operator to be the dot product"""
+        self.name = name
+        Matrix.__init__(self, m=m, *args, **kwargs)
+        self.func = func or self.operator_func
+        LinearOperator.__init__(self, func=func, *args, **kwargs)
+
+    def _is_linear(self):
+        # Every matrix transformation is a linear transformation
+        # https://www.mathbootcamps.com/proof-every-matrix-transformation-is-a-linear-transformation/
+        return True
+
+    def operator_func(self, other):
+        return Vector(np.dot(self.m, other.m))
+
+    def get_eigens(self):
+        """ Returns (eigenvalue, eigenvector)
+        M|v> = λ|v>  ->
+          M   :is the operator (self)
+          |v> :is eigenstate of M
+          λ   :is the corresponding eigenvalue
+        """
+        eigenvalues, eigenvectors = np.linalg.eig(self.m)
+        rv = []
+        for i in range(0, len(eigenvectors)):
+            eigenvalue = eigenvalues[i]
+            eigenvector = HorizontalVector(eigenvectors[:, i])
+            rv.append((eigenvalue, eigenvector))
+        return rv
+
+
 class UnitaryMatrix(SquareMatrix):
     def __init__(self, m: ListOrNdarray, *args, **kwargs):
         """Represents a Unitary matrix that satisfies UU+ = I"""
@@ -407,13 +462,13 @@ class HermitianMatrix(SquareMatrix):
         return np.isclose(self.m, self._conjugate_transpose()).all()
 
 
-class UnitaryOperator(LinearOperator, UnitaryMatrix):
+class UnitaryOperator(LinearTransformation, UnitaryMatrix):
     def __init__(self, m: ListOrNdarray, name: str = '', *args, **kwargs):
-        """UnitaryOperator inherits from both LinearOperator and a Unitary matrix
+        """UnitaryOperator inherits from both LinearTransformation and a Unitary matrix
         It is used to act on a State vector by defining the operator to be the dot product"""
         self.name = name
         UnitaryMatrix.__init__(self, m=m, *args, **kwargs)
-        LinearOperator.__init__(self, func=self.operator_func, *args, **kwargs)
+        LinearTransformation.__init__(self, m=m, func=self.operator_func, *args, **kwargs)
 
     def operator_func(self, other):
         return State(np.dot(self.m, other.m))
@@ -424,9 +479,9 @@ class UnitaryOperator(LinearOperator, UnitaryMatrix):
         return str(self.m)
 
 
-class HermitianOperator(LinearOperator, HermitianMatrix):
+class HermitianOperator(LinearTransformation, HermitianMatrix):
     def __init__(self, m: ListOrNdarray, name: str = '', *args, **kwargs):
-        """HermitianMatrix inherits from both LinearOperator and a Hermitian matrix
+        """HermitianMatrix inherits from both LinearTransformation and a Hermitian matrix
         It is used to act on a State vector by defining the operator to be the dot product"""
         self.name = name
         HermitianMatrix.__init__(self, m=m, *args, **kwargs)
@@ -596,6 +651,7 @@ def Rz(theta):
                             [0, np.power(np.e, 1j * theta / 2)]],
                            name="Rz")
 
+
 H = UnitaryOperator([[1 / np.sqrt(2), 1 / np.sqrt(2)],
                      [1 / np.sqrt(2), -1 / np.sqrt(2)], ],
                     name="H")
@@ -604,7 +660,7 @@ CNOT = TwoQubitOperator([[1, 0, 0, 0],
                          [0, 1, 0, 0],
                          [0, 0, 0, 1],
                          [0, 0, 1, 0], ],
-                         C_partial, x_partial, I, X)
+                        C_partial, x_partial, I, X)
 
 
 def assert_raises(exception, msg, callable, *args, **kwargs):
@@ -667,7 +723,11 @@ class MeasurementOpeartor(HermitianOperator):
         return cls(matrix.conjugate_transpose().x(matrix).m)
 
     def get_prob(self, state: State):
+        """Returns result of <ψ|M†_m M_m|ψ>
+        """
+        # <ψ|
         state_ct = state.conjugate_transpose()
+        #  This is: <ψ|   .   M†_m M_m . |ψ>
         return state_ct.m.dot(self.m.dot(state.m)).item()
 
 
@@ -707,6 +767,11 @@ def testRotPauli():
     assert np.allclose(abs_squared(Rz(np.pi).m), abs_squared(Z.m))
 
 
+def test_eigenstuff():
+    assert LinearTransformation(m=[[1, 0], [0, 0]]).get_eigens() == \
+           [(1.0, HorizontalVector([1., 0.])), (0., HorizontalVector([0., 1.]))]
+
+
 def test():
     # Test properties of Hilbert vector space
     # The four postulates of Quantum Mechanics
@@ -734,11 +799,9 @@ def test():
 
     # II: Dynamics | The evolution of a closed system is described by a unitary transformation
     #
-    # Operators turn one vector into another
+    test_linear_operator()
-    # the times 2 operator should return the times two multiplication
+
-    _times_2 = Operator(lambda x: 2 * x)
+    test_eigenstuff()
-    assert _times_2.on(5) == 10
-    assert _times_2(5) == 10
 
     # Understanding the difference between unitary and hermitians
     test_unitary_hermitian()