Using Row Operations

The row operations that I implemented earlier can now be put to use - namely for making a particular element in a matrix zero 'MakeElementZero' takes a matrix M, a row A, a row B, a position P in the row, and does a row operation to make M[A][P] equal to zero. 'MakeColumnZero' takes a matrix M, a specified column, and calls MakeElementZero to insert zeros appropriately (into upper triangular form). The next step here is to recursively use MakeColumnZero to make a matrix upper triangular (good for finding the determinant and solving linear equations/inverting).

-- Abishek (AR-MA210)

Code/LinearAlgebra.idr

@@ -1,24 +1,35 @@
 module LinearAlgebra
 
-%access public export
-
 import BaseN
-import MultiSolver
-import MatrixNumeric
+import Linear.MultiSolver
+import Data.Matrix.Numeric
 import Data.Vect
 import Data.Fin
-import Matrix
+import Merge
+import Data.Matrix
 import ZZ
-import Rationals_2
+import Rationals
+
+%access public export
 
 -- Some auxillary functions for elementary operations
 
-SimpleCast:  (v1:(Vect (k+ (S Z)) elem)) ->( Vect (S k) elem)  -- necessary for type matching
+FST: ZZPair -> ZZ --some issue with fst
+FST (a, b) = a
+
+Pred: Nat -> Nat -- needed to get (n-1) for number of operations
+Pred Z = Z
+Pred (S k) = k
+
+SimpleCast:  (v1:(Vect (k+ (S Z)) elem)) ->( Vect (S k) elem)  -- necessary for type matching, taken from Shafil's code
 SimpleCast {k} v1 = (Vecttake (S k) v1)
 
 RowN: Matrix n n ZZPair -> (k: Nat) -> Vect n ZZPair -- returns the nth row (indexing from 0)
 RowN {n} x k = Data.Vect.index (tofinNat k n) x
 
+ij: Matrix n n ZZPair -> (i: Nat) -> (j: Nat) -> ZZPair
+ij {n} x i j = Data.Vect.index (tofinNat j n) (RowN x i)
+
 ColumnN: Matrix n n ZZPair -> (k: Nat) -> Vect n ZZPair --returns the nth column (indexing from 0)
 ColumnN {n} x k = Data.Vect.index (tofinNat k n) (transpose x)
 
@@ -53,7 +64,43 @@ First_k_rows_from_RowP n x p (S k) = SimpleCast ((First_k_rows_from_RowP n x p k
 SwapRows: (n: Nat) -> Matrix n n ZZPair -> (a: Nat) -> (b: Nat) -> Matrix n n ZZPair
 SwapRows n x a b = replaceAt (tofinNat b n) (RowN x a) (replaceAt (tofinNat a n) (RowN x b) (x))
 
--- Performs an elementary row operation, subtracts from ROW A the ROW B scaled by C (A -> A - C*B)
+-- Performs an elementary row operation, subtracts from row a the row b scaled by c
 
 RowOperation: (n: Nat) -> Matrix n n ZZPair -> (a: Nat) -> (b: Nat) -> (c: ZZPair) -> Matrix n n ZZPair
 RowOperation n x a b c = replaceAt (tofinNat a n) (axpy n (RowN x a) (RowN x b) c) (x)
+
+-- Makes the kth number in Row A zero by subtracting a scaling of Row B (for upper triangularizing)
+-- As usual, indexing starts from 0
+
+MakeElementZero: (n: Nat) -> (x: Matrix n n ZZPair) -> (row_a : Nat) -> (row_b : Nat) -> (pos : Nat) -> Matrix n n ZZPair
+MakeElementZero n x row_a row_b pos = case (FST(ij x row_a pos)) of
+                                      (Pos Z) => x
+                                      (Pos (S k)) => case ((fst (ij x row_b pos))) of
+                                                                          (Pos Z) => (SwapRows n x row_a row_b)
+                                                                          (Pos (S k)) => (RowOperation n x row_a row_b (divZZ (ij x row_a pos) (ij x row_b pos)) )
+                                                                          (NegS k) => (RowOperation n x row_a row_b (divZZ (ij x row_a pos) (ij x row_b pos)) )
+                                      (NegS k) => case ((fst (ij x row_b pos))) of
+                                                                          (Pos Z) => (SwapRows n x row_a row_b)
+                                                                          (Pos (S k)) => (RowOperation n x row_a row_b (divZZ (ij x row_a pos) (ij x row_b pos)) )
+                                                                          (NegS k) => (RowOperation n x row_a row_b (divZZ (ij x row_a pos) (ij x row_b pos)) )
+
+-- Explanation: We want to transform the element x[row_a][pos] into zero by a row operation. There are a few different cases
+--              1. If x[row_a][pos] is already zero, nothing needs to be done.
+--              2. If not, we scale row_b appropriately and do a row operation (A -> A - cB)
+--              3. If x[row_b][pos] is zero, however, no row operation will work. We simply swap the
+--                 row_a and row_b here; then, x[row_a][pos] will become zero.
+
+-- This algorithm is important to make a matrix upper-triangular.
+
+MakeColumnZero: (n: Nat) -> (x: Matrix n n ZZPair) -> (col : Nat) -> (iter : Nat) -> Matrix n n ZZPair
+MakeColumnZero n x col Z = x
+MakeColumnZero n x col (S k) = case (isLTE (S k) col) of
+                                   (Yes prf) => x
+                                   (No contra) => (MakeColumnZero n (MakeElementZero n x (S k) col col) col k)
+
+-- This function turns a column into what it should be in upper-triangular form by adding in the necessary zeros.
+-- iter is a variable to induct on (trick courtesy Sriram)
+-- When using this, make sure to set "iter" as n-1 (1 less than the number of rows)
+
+-- The next step here is to use the above function to make a matrix upper triangular. This can be
+-- done by inducting on the number of columns.