diff --git a/README b/README
index c883e80..da90d01 100644
--- a/README
+++ b/README
@@ -1,6 +1,50 @@
+Linear Algebra in Common Lisp
+
+Copyright (c) 2011-2014, Odonata Research LLC
+
+Permission is hereby granted, free  of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction,  including without limitation the rights
+to use, copy, modify,  merge,  publish,  distribute,  sublicense, and/or sell
+copies of the  Software,  and  to  permit  persons  to  whom  the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and  this  permission  notice shall be included in
+all copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED  "AS IS",  WITHOUT  WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT  NOT  LIMITED  TO  THE  WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE  AND  NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT  HOLDERS  BE  LIABLE  FOR  ANY  CLAIM,  DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.
+
+Copyright (c) 2021-2023, Ten Factor Growth, LLC
+Permission is hereby granted, free  of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction,  including without limitation the rights
+to use, copy, modify,  merge,  publish,  distribute,  sublicense, and/or sell
+copies of the  Software,  and  to  permit  persons  to  whom  the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and  this  permission  notice shall be included in
+all copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED  "AS IS",  WITHOUT  WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT  NOT  LIMITED  TO  THE  WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE  AND  NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT  HOLDERS  BE  LIABLE  FOR  ANY  CLAIM,  DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.
+
+
 linear-algebra
 
 linear-algebra is a Common Lisp library of numeric linear algebra routines.
 It approximates the functionality present in the Basic Linear Algebra
 Subroutines (BLAS) Fortran library. The objective of this project is to
 generate linear algebra routines in native ANSI Common Lisp.
+
+Use within quicklisp with (ql:quickload :linear-algebra)
diff --git a/lisp/data-vector.lisp b/lisp/data-vector.lisp
index de80077..3d39391 100644
--- a/lisp/data-vector.lisp
+++ b/lisp/data-vector.lisp
@@ -36,6 +36,10 @@
   (:documentation
    "A data vector."))
 
+(defmethod print-object ((vector data-vector) stream)
+  (print-unreadable-object (vector stream :type t :identity t)
+    (format stream "~a" (contents vector))))
+
 (defclass row-vector (data-vector)
   ()
   (:documentation
@@ -200,6 +204,9 @@ applying the function to each element of the vectors."
          function
          vectors))
 
+(defmethod reduce-vector ((function function) (vector data-vector) initial-value)
+  (reduce function (contents vector) :initial-value initial-value))
+
 ;;; Data vector transformations
 
 (defmethod apply-rotation :before
@@ -260,6 +267,10 @@ applying the function to each element of the vectors."
   "Return the p-norm of the vector."
   (norm-vector (contents vector) measure))
 
+(defmethod min-vector ((vector data-vector))
+  (loop for element across (contents vector) minimize (abs element)))
+
+
 (defmethod transpose ((vector column-vector))
   "Return a row vector."
   (make-instance
@@ -434,3 +445,101 @@ applying the function to each element of the vectors."
     ((vector1 row-vector) (vector2 column-vector) &optional scalar)
   "Return the dot product of vector1 and vector2."
   (inner-product-vector (contents vector1) (contents vector2) scalar))
+
+(defmethod product :before
+  ((vector1 data-vector) (vector2 data-vector) &optional scalar)
+  "Verify that the dimensions are equal."
+  (declare (ignore scalar))
+  (unless (= (vector-length vector1) (vector-length vector2))
+    (error "VECTOR1 and VECTOR2 are not of equal length.")))
+
+(defmethod product
+    ((vector1 data-vector) (vector2 data-vector) &optional scalar)
+  "Return the dot product of vector1 and vector2."
+  (inner-product-vector (contents vector1) (contents vector2) scalar))
+
+(defmethod vec-equal :before
+  ((vector1 data-vector) (vector2 data-vector))
+  "Verify that the dimensions are equal."
+  (declare (ignore scalar))
+  (unless (= (vector-length vector1) (vector-length vector2))
+    (error "VECTOR1 and VECTOR2 are not of equal length.")))
+
+(defmethod vec-equal
+    ((vector1 data-vector) (vector2 data-vector))
+  (equalp (contents vector1) (contents vector2)))
+
+
+(defmethod elem-divide :before
+  ((vector1 data-vector) (vector2 data-vector))
+  "Verify that the dimensions are equal."
+  (declare (ignore scalar))
+  (unless (= (vector-length vector1) (vector-length vector2))
+    (error "VECTOR1 and VECTOR2 are not of equal length.")))
+
+(defmethod elem-divide
+    ((vector1 data-vector) (vector2 data-vector))
+  (make-instance
+   (common-class-of vector1 vector2)
+   :contents
+   (element-divide-vector
+    (contents vector1) (contents vector2))
+   ))
+
+(defmethod elem-multiply :before
+  ((vector1 data-vector) (vector2 data-vector))
+  "Verify that the dimensions are equal."
+  (declare (ignore scalar))
+  (unless (= (vector-length vector1) (vector-length vector2))
+    (error "VECTOR1 and VECTOR2 are not of equal length.")))
+
+(defmethod elem-multiply
+    ((vector1 data-vector) (vector2 data-vector))
+  (make-instance
+   (common-class-of vector1 vector2)
+   :contents
+   (element-multiply-vector
+    (contents vector1) (contents vector2))
+   ))
+
+(defmethod elem-greater
+    ((vector1 data-vector) (vector2 data-vector))
+  (element-greater-vector
+   (contents vector1) (contents vector2)))
+
+(defmethod vec-every ((vector data-vector) predicate)
+  (every predicate (contents vector)))
+    
+
+(defmethod distance :before
+  ((vector1 data-vector) (vector2 data-vector) &optional (measure 2))
+  "Verify that the dimensions are equal."
+  (declare (ignore scalar))
+  (unless (= (vector-length vector1) (vector-length vector2))
+    (error "VECTOR1 and VECTOR2 are not of equal length.")))
+
+(defmethod distance
+    ((vector1 data-vector) (vector2 data-vector) &optional (measure 2))
+  (norm (subtract vector1 vector2) measure))
+
+(defmethod outer-product ((vector1 data-vector)
+                          (vector2 data-vector))
+  (let* ((len1 (vector-length vector1))
+         (len2 (vector-length vector2))
+         (mat-type (cond
+                     ((equal vector1 vector2) 'symmetric-matrix)
+                     ((= len1 len2) 'square-matrix)
+                     (t 'dense-matrix)))
+         (mat (make-instance mat-type
+                             :dimensions (list len1 len2)
+                             :initial-element 0
+                             :element-type (vector-element-type vector1)))
+         )
+    (dotimes (i len1)
+      (dotimes (j len2)
+        (setf
+         (mref mat i j)
+         (* (vref vector1 i) (vref vector2 j)))
+        ))
+    mat
+  ))
diff --git a/lisp/dense-matrix.lisp b/lisp/dense-matrix.lisp
index f06b28b..73e088d 100644
--- a/lisp/dense-matrix.lisp
+++ b/lisp/dense-matrix.lisp
@@ -34,6 +34,10 @@
   (:documentation
    "Dense matrix object."))
 
+(defmethod print-object ((matrix dense-matrix) stream)
+  (print-unreadable-object (matrix stream :type t :identity t)
+    (format stream "~a" (contents matrix))))
+
 (defmethod initialize-instance :after
   ((self dense-matrix) &rest initargs
    &key dimensions element-type initial-element initial-contents)
@@ -110,6 +114,10 @@
   "Return true if object is a dense matrix."
   (typep object 'dense-matrix))
 
+(defmethod mat-equal
+    ((matrix1 dense-matrix) (matrix2 dense-matrix))
+  (equalp (contents matrix1) (contents matrix2)))
+
 (defmethod matrix-in-bounds-p
     ((matrix dense-matrix) (row integer) (column integer))
   "Return true if row and column do not exceed the dimensions of matrix."
@@ -490,3 +498,44 @@ matrix with a column vector."
        :contents
        (gauss-invert (contents matrix)))
       (error "The number of rows does not equal columns.")))
+
+(defmethod add-diagonal ((value number) (matrix dense-matrix))
+  (let ((result (copy-matrix matrix))
+	(num    (apply #'min (matrix-dimensions matrix)))
+	)
+    (dotimes (index num result)
+      (setf
+       (mref result index index)
+       (+ (mref matrix index index) value)))
+    result))
+
+(defmethod add-diagonal :before
+    ((vec data-vector) (matrix dense-matrix))
+  (unless (compatible-dimensions-p :solve matrix vector)
+    (error "Matrix~A is incompatible with column vector(~D)."
+           (matrix-dimensions matrix) (vector-length vector))))
+
+(defmethod add-diagonal ((vec data-vector) (matrix dense-matrix))
+  (let ((result (copy-matrix matrix))
+	(num    (min (matrix-dimensions matrix)))
+	)
+    (dotimes (index num result)
+      (setf
+       (mref result index index)
+       (+ (mref matrix index index)
+	  (vref vec index))))
+    result))
+
+(defmethod matrix-diagonal ((matrix dense-matrix))
+  (let* ((num (apply #'min (matrix-dimensions matrix)))
+         (vec (make-instance 'data-vector :size num
+                                          :initial-element 0
+                                          :element-type (matrix-element-type matrix)))
+         )
+    (dotimes (index num )
+      (setf
+       (vref vec    index)
+       (mref matrix index index))
+      )
+    vec
+  ))
diff --git a/lisp/interface/fundamental-ops.lisp b/lisp/interface/fundamental-ops.lisp
index 3d142d1..44fabe6 100644
--- a/lisp/interface/fundamental-ops.lisp
+++ b/lisp/interface/fundamental-ops.lisp
@@ -32,6 +32,14 @@
   (:documentation
    "Return the norm according to measure."))
 
+(defgeneric min-vector (vector)
+  (:documentation
+   "Return the minimum absolute value"))
+
+(defgeneric distance (vector1 vector2 &optional measure)
+  (:documentation
+   "Return the distance between vector1 and vector2 using the measure"))
+
 (defgeneric transpose (vector-or-matrix)
   (:documentation
    "Transpose the vector or matrix."))
@@ -92,3 +100,8 @@
 (defgeneric ninvert (matrix)
   (:documentation
    "Return the invert of the matrix with in-place decomposition."))
+
+(defgeneric add-diagonal (scalar-or-vector matrix)
+  (:documentation
+   "Add a scalar or vector to every diagonal of the matrix"))
+
diff --git a/lisp/interface/matrix.lisp b/lisp/interface/matrix.lisp
index a8f84fa..1eb72d3 100644
--- a/lisp/interface/matrix.lisp
+++ b/lisp/interface/matrix.lisp
@@ -99,6 +99,10 @@ MATRIX."))
   (:documentation
    "Set the submatrix of the matrix."))
 
+(defgeneric matrix-trace (matrix)
+  (:documentation
+   "Return a column vector that is the trace of matrix"))
+
 (defgeneric replace-matrix
     (matrix1 matrix2 &key
      start-row1 end-row1
@@ -120,3 +124,11 @@ MATRIX."))
         (values start-row start-column end-row end-column)
         (error "The matrix range (~D:~D,~D:~D) is invalid."
                start-row start-column end-row end-column))))
+
+(defgeneric mat-equal (matrix1 matrix2)
+  (:documentation
+   "Return if the two arrays elments are the same"))
+
+(defgeneric matrix-diagonal (matrix)
+  (:documentation
+   "Return the diagonal of the matrix as a vector"))
diff --git a/lisp/interface/vector.lisp b/lisp/interface/vector.lisp
index 67a5d31..e802e29 100644
--- a/lisp/interface/vector.lisp
+++ b/lisp/interface/vector.lisp
@@ -89,6 +89,10 @@
    "Destructively modifies the result vector with the result of
 applying the function to each element of the vectors."))
 
+(defgeneric reduce-vector (function vector initial-value)
+  (:documentation
+   "Reduce a vector into a scalar using a function"))
+
 (defmacro dovector ((element vector &optional result) &body body)
   "Iterate over vector returning result."
   (let ((pos (gensym "POS-"))
@@ -108,3 +112,28 @@ applying the function to each element of the vectors."))
 (defgeneric napply-rotation (vector1 vector2 cc ss)
   (:documentation
    "Return the plane rotations of vector1 and vector2 by cc and ss."))
+
+(defgeneric vec-equal (vector1 vector2)
+  (:documentation
+   "Test if the contents of vector1 and vector2 are the same"))
+
+(defgeneric elem-divide (vector1 vector2)
+  (:documentation
+   "Return the element by element division of two vectors"))
+
+(defgeneric elem-multiply (vector1 vector2)
+  (:documentation
+   "Return the element by element multiplation of two vectors"))
+
+(defgeneric elem-greater (vector1 vector2)
+  (:documentation
+   "Return if vector1 is greater than vector 2 on every element"))
+
+(defgeneric vec-every (vector predicate)
+  (:documentation
+   "Return true if every element in the vector satisfies
+  predicate. Return nil on the first element that does not"))
+
+(defgeneric outer-product (vec1 vec2)
+  (:documentation
+   "Outer product of vec1 and vec2 returns an nxm matrix"))
diff --git a/lisp/kernel/binary-operations.lisp b/lisp/kernel/binary-operations.lisp
index 5cfa579..2dec6d4 100644
--- a/lisp/kernel/binary-operations.lisp
+++ b/lisp/kernel/binary-operations.lisp
@@ -75,6 +75,16 @@ the operation."))
   "Return the scaled operation."
   (lambda (n1 n2) (- (* scalar1 n1) (* scalar2 n2))))
 
+(defmethod binary-op
+    ((op (eql #'/)) )
+  "Return the scaled operation."
+  (lambda (n1 n2) (/ n1 n2)))
+
+(defmethod binary-op
+    ((op (eql #'*)) )
+  "Return the scaled operation."
+  (lambda (n1 n2) (* n1 n2)))
+
 ;;; Binary vector operations
 
 (defun %vector<-vector1-op-vector2 (operation vector1 vector2)
@@ -124,6 +134,27 @@ addition."
    (scaled-binary-op #'- scalar1 scalar2)
    vector1 vector2))
 
+(defun element-divide-vector (vector1 vector2)
+  "Vector pointwise division."
+  (%vector<-vector1-op-vector2
+   (binary-op #'/)
+   vector1 vector2))
+
+(defun element-multiply-vector (vector1 vector2)
+  "Vector pointwise multiplication."
+  (%vector<-vector1-op-vector2
+   (binary-op #'*)
+   vector1 vector2))
+
+(defun element-greater-vector (vector1 vector2)
+  (loop
+   for element1 across vector1
+   and element2 across vector2
+   when (<= element1 element2)
+     return nil
+   finally (return t)))
+
+
 (defun nsubtract-vector (vector1 vector2 scalar1 scalar2)
   "Destructive vector binary subtraction."
   (%vector1<-vector1-op-vector2
@@ -137,7 +168,7 @@ addition."
    and element2 across vector2
    sum (* element1 element2) into result
    finally
-   (return (if scalar (* scalar result) result))))
+      (return (if scalar (* scalar result) result))))
 
 ;;; Binary array/vector operations
 
diff --git a/lisp/kernel/linear-algebra-kernel.lisp b/lisp/kernel/linear-algebra-kernel.lisp
index 54b1bd0..5000233 100644
--- a/lisp/kernel/linear-algebra-kernel.lisp
+++ b/lisp/kernel/linear-algebra-kernel.lisp
@@ -52,6 +52,9 @@
            :subtract-vector :nsubtract-vector
            :add-array :nadd-array
            :subtract-array :nsubtract-array
+           :element-multiply-vector
+           :element-divide-vector
+           :element-greater-vector
            :inner-product-vector
            :product-vector-array
            :product-array-vector
diff --git a/lisp/kernel/unary-operations.lisp b/lisp/kernel/unary-operations.lisp
index 3421c20..2eedc29 100644
--- a/lisp/kernel/unary-operations.lisp
+++ b/lisp/kernel/unary-operations.lisp
@@ -193,6 +193,7 @@ array column."
     (dotimes (index (length vector) result)
       (setf (aref result index) (abs (aref vector index))))))
 
+
 (defmethod norm-vector ((data vector) (measure (eql 1)))
   "Return the Taxicab norm of the list."
   (loop for element across data sum (abs element)))
@@ -213,6 +214,7 @@ array column."
   "Return the infinity, or maximum, norm of vector."
   (loop for element across data maximize (abs element)))
 
+
 (defmethod norm-array ((data array) (measure (eql 1)))
   "Return the 1 norm of the array."
   (let ((m-rows (array-dimension data 0))
diff --git a/lisp/linear-algebra.lisp b/lisp/linear-algebra.lisp
index 1586a88..b31f16b 100644
--- a/lisp/linear-algebra.lisp
+++ b/lisp/linear-algebra.lisp
@@ -31,6 +31,7 @@
   (:use :floating-point :linear-algebra-kernel)
   ;; Fundamental operations
   (:export :norm
+	   :min-vector
            :transpose :ntranspose
            :permute
            :scale :nscale
@@ -51,8 +52,17 @@
            :replace-vector
            :map-vector
            :map-into-vector
+	   :reduce-vector
            :dovector
-           :apply-rotation :napply-rotation)
+           :vec-equal
+   	   :elem-divide
+           :elem-multiply
+           :vec-every
+           :outer-product
+           :elem-greater
+   	   :distance
+           :apply-rotation :napply-rotation
+   )
   ;; Matrix interface
   (:export :matrix-object
            :initialize-matrix
@@ -65,12 +75,18 @@
            :matrix-column-dimension
            :mref
            :copy-matrix
+           :mat-equal
+           :matrix-diagonal
            :submatrix
            :replace-matrix
-           :matrix-validated-range)
+           :matrix-validated-range
+	   :add-diagonal
+	   )
   ;; Identity matrix
-  (:export :identity-matrix
-           :identity-matrix-p)
+  (:export
+   :identity-matrix
+   :identity-matrix-p
+   )
   ;; Permutation matrix
   (:export :permutation-matrix
            :permutation-matrix-p)
diff --git a/lisp/sequence/vector.lisp b/lisp/sequence/vector.lisp
index dace026..75619de 100644
--- a/lisp/sequence/vector.lisp
+++ b/lisp/sequence/vector.lisp
@@ -110,3 +110,4 @@
       (inner-product-vector vector1 vector2 scalar)
       (error "VECTOR1(~D) and VECTOR2(~D) are not of equal length."
              (length vector1) (length vector2))))
+
diff --git a/test/data-vector.lisp b/test/data-vector.lisp
index e0b3589..d548025 100644
--- a/test/data-vector.lisp
+++ b/test/data-vector.lisp
@@ -1207,18 +1207,30 @@
    (linear-algebra:product
     (linear-algebra:row-vector 1 2 3)
     (linear-algebra:column-vector 1 2 3 4)))
-  (assert-error
-   'error
+  (assert-equal
+   14
    (linear-algebra:product
     (linear-algebra:column-vector 1 2 3)
     (linear-algebra:column-vector 1 2 3)))
-  (assert-error
-   'error
+  (assert-equal
+   14
    (linear-algebra:product
     (linear-algebra:row-vector 1 2 3)
     (linear-algebra:row-vector 1 2 3)))
-  (assert-error
-   'error
+  (assert-equal
+   14
    (linear-algebra:product
     (linear-algebra:column-vector 1 2 3)
     (linear-algebra:row-vector 1 2 3))))
+
+;; test the inner product
+(define-test inner-product-vector
+  (:tag :data-vector :product)
+  ;; Real vectors
+  (assert-rational-equal
+  14
+   (linear-algebra-kernel:inner-product-vector
+    (linear-algebra:row-vector 1 2 3)
+    (linear-algebra:row-vector 1 2 3)
+    nil))
+)
diff --git a/test/dense-matrix.lisp b/test/dense-matrix.lisp
index 5307109..f3b77df 100644
--- a/test/dense-matrix.lisp
+++ b/test/dense-matrix.lisp
@@ -1363,3 +1363,40 @@
          ( 0.031  0.015815994 29.992375     0.33191067)
          ( 0.041  0.020756614  0.014001981 39.98469))
      matrix)))
+
+(define-test add-diagonal-dense-matrix
+  (:tag :dense-matrix :add-diagonal)
+  (let ((*epsilon* (* 64 single-float-epsilon))
+	(scalar 2.0)
+        (vector1 (linear-algebra:column-vector 1.0 2.0))
+	(vector2 (linear-algebra:column-vector 1.0 2.0 3.0))
+        (matrix1
+         (linear-algebra:make-matrix
+          2 2 :initial-contents '((1.1 1.2) (2.1 2.2))))
+        (result1
+         (linear-algebra:make-matrix
+          2 2 :initial-contents
+	  '((3.1 1.2) (2.1 4.2))))
+	(result2
+         (linear-algebra:make-matrix
+          2 2 :initial-contents
+	  '((2.1 1.2) (2.1 4.2))))
+	)
+
+    ;; 2x2
+    (assert-true
+     (linear-algebra:mat-equal
+      result1
+      (linear-algebra:add-diagonal scalar matrix1)))
+
+    (assert-true
+     (linear-algebra:mat-equal
+      result2
+      (linear-algebra:add-diagonal vector1 matrix1)))
+
+    (assert-error
+     'error
+     (linear-algebra:add-diagonal vector2 matrix1))
+    ))
+     
+    
diff --git a/test/linear-algebra-test.lisp b/test/linear-algebra-test.lisp
index 4b34341..ee62e2f 100644
--- a/test/linear-algebra-test.lisp
+++ b/test/linear-algebra-test.lisp
@@ -7,6 +7,7 @@
 
 (in-package :linear-algebra-test)
 
+
 ;;; Convenience functions
 
 (defun random-interior-index (size)
diff --git a/todo.org b/todo.org
new file mode 100644
index 0000000..f4dcfd6
--- /dev/null
+++ b/todo.org
@@ -0,0 +1,11 @@
+* TODO convert the tests to fiveam
+* TODO move the types into a types subdirectory
+** create a types
+** refactor dense-matrix to remove the operations and put them in fundamental-ops
+* TODO get identity matrix to work
+** put it in types
+** create a diagonal matrix type and make identity a subclass
+* TODO add banded matrix type
+* TODO add SVD decomposition for Hermetian
+* TODO add gram-schmidt QR
+* TODO add sequential gram-schmidt as a square root op