[lisplab-cvs] r29 - src/core src/matrix

Jørn Inge Vestgården jivestgarden at common-lisp.net
Thu May 21 14:02:27 UTC 2009


Author: jivestgarden
Date: Thu May 21 10:02:18 2009
New Revision: 29

Log:
refactoring almost complete. Should now remove old code

Modified:
   TODO
   src/core/level0-interface.lisp
   src/matrix/level1-array.lisp
   src/matrix/level1-constructors.lisp
   src/matrix/level1-interface.lisp
   src/matrix/level2-constructors.lisp
   src/matrix/level2-generic.lisp
   src/matrix/level2-interface.lisp
   src/matrix/level2-matrix-dge.lisp
   src/matrix/level2-matrix-zge.lisp

Modified: TODO
==============================================================================
--- TODO	(original)
+++ TODO	Thu May 21 10:02:18 2009
@@ -1,4 +1,8 @@
 TODO
+o About arrays. Handle them by making wrappers in the matrix
+  hierarcy. This will be rather slow but saves work. 
+  And in this way it is possible to avoid the non-speialized methods.     
+o Matrices with general element type.
 o Some way to just extend the exact input matrix type
 o Documentation.
 o Check if there is some non-threadsafe code in the fortran interface,

Modified: src/core/level0-interface.lisp
==============================================================================
--- src/core/level0-interface.lisp	(original)
+++ src/core/level0-interface.lisp	Thu May 21 10:02:18 2009
@@ -1,6 +1,6 @@
 ;;; Lisplab, level0-interface.lisp
-;;; Generic function definitions that also contain 
-;;; non-matrix methods.
+;;; Defines a basic algebra.
+;;; 
 
 ;;; Copyright (C) 2009 Joern Inge Vestgaarden
 ;;;
@@ -22,6 +22,7 @@
 
 (export '(copy convert
 	  scalar?
+	  vector? matrix?
 	  .abs .imagpart .realpart
 	  .= ./= .< .<= .> .>= 
 	  .add .add! 
@@ -38,6 +39,15 @@
 	  .erf .erfc 
 	  .gamma)) 	  
 
+(defgeneric scalar? (x)
+  (:documentation "A scalar is a object with ignored internal structure."))
+
+(defgeneric vector? (x)
+  (:documentation "A vector is a object whose elements are accessible with vref."))
+
+(defgeneric matrix? (x)
+  (:documentation "A matrix is a object whose elements are accesible with mref."))
+
 (defgeneric copy (a)
   (:documentation "Copies the elements and structure, but ignore 
 shared state, like fill pointers etc."))
@@ -45,9 +55,6 @@
 (defgeneric convert (x type)
   (:documentation "Generalized coerce."))
 
-(defgeneric scalar? (x)
-  (:documentation "A scalar is a object with ignored internal structure."))
-
 (defgeneric .abs (a)
   (:documentation "Generialized abs"))
 
@@ -57,6 +64,8 @@
 (defgeneric .imagpart (a)
   (:documentation "Generialized abs"))
 
+;;; Binary boolean operators
+
 (defgeneric .= (a b &optional (accuracy))
   (:documentation "Element-wise test of equality, with presition."))
 
@@ -75,6 +84,8 @@
 (defgeneric .>= (a b) 
   (:documentation "Generalized >=." )) 
 
+;;; Binary operators
+
 (defgeneric .add (a b)
   (:documentation "Addes a and b elementwise. Called by .+"))
 

Modified: src/matrix/level1-array.lisp
==============================================================================
--- src/matrix/level1-array.lisp	(original)
+++ src/matrix/level1-array.lisp	Thu May 21 10:02:18 2009
@@ -24,41 +24,25 @@
   (= (rank a) 2))
 
 (defmethod vector? ((a array))
-  "True for an array of rank 2"
-  (= (rank a) 1))
+  "True for any array through row-major-aref"
+  t)
 
-(defmethod copy ((a array))
-  (if (vector? a)
-      (copy-seq a)
-      (let ((y (make-array (dim a) :element-type (element-type a))))
-	(dotimes (i (size a))
-	  (setf (row-major-aref y i)
-		(row-major-aref a i)))
-	y)))
+(defmethod dim ((a array) &optional axis)
+  (if axis
+      (array-dimension a axis)
+      (array-dimensions a)))
 
-(defmethod new ((class (eql 'array)) dim &optional (element-type t) (value 0))
-  (unless (consp dim) (setf dim (list dim 1)))
-  (make-array dim
-	      :element-type element-type
-	      :initial-element (convert value element-type)))
+(defmethod size ((a array))
+  (reduce #'* (dim a)))
 
-(defmethod new ((class (eql 'simple-array)) dim &optional (element-type t) (value 0))
-  (unless (consp dim) (setf dim (list dim 1)))
-  (make-array dim
-	      :element-type element-type
-	      :initial-element (convert value element-type)))
+(defmethod rank ((a array))
+  (array-rank a))
 
-(defmethod new ((class (eql 'vector)) dim &optional (element-type t) (value 0))
-  (unless (consp dim) (setf dim (list dim 1)))
-  (make-array dim
-	      :element-type element-type
-	      :initial-element (convert value element-type)))
+(defmethod rows ((a array))
+  (array-dimension a 0))
 
-(defmethod new ((class (eql 'simple-vector)) dim &optional (element-type t) (value 0))
-  (unless (consp dim) (setf dim (list dim 1)))
-  (make-array dim
-	      :element-type element-type
-	      :initial-element (convert value element-type)))
+(defmethod cols ((a array))
+  (array-dimension a 1))
 
 (defmethod element-type ((a array))
   "Gets the element type of the array"
@@ -78,21 +62,39 @@
 (defmethod (setf vref) (value (a array) idx)
   (setf (row-major-aref a idx) value)) 
 
-(defmethod dim ((a array) &optional axis)
-  (if axis
-      (array-dimension a axis)
-      (array-dimensions a)))
+(defmethod make-matrix-instance ((x (eql 'array)) dim value)
+  (make-array dim :initial-element value))
 
-(defmethod vector? ((a array))
-  "True for an array of rank 1"
-  (= (rank a) 1))
+(defmethod copy ((a array))
+  ;; TODO move to level2 
+  (if (vectorp a)
+      (copy-seq a)
+      (let ((y (make-array (dim a) :element-type (element-type a))))
+	(dotimes (i (size a))
+	  (setf (row-major-aref y i)
+		(row-major-aref a i)))
+	y)))
 
-(defmethod rank ((a array))
-  (array-rank a))
+#+nil (defmethod new ((class (eql 'array)) dim &optional (element-type t) (value 0))
+  (unless (consp dim) (setf dim (list dim 1)))
+  (make-array dim
+	      :element-type element-type
+	      :initial-element (convert value element-type)))
 
-(defmethod rows ((a array))
-  (array-dimension a 0))
+#+nil (defmethod new ((class (eql 'simple-array)) dim &optional (element-type t) (value 0))
+  (unless (consp dim) (setf dim (list dim 1)))
+  (make-array dim
+	      :element-type element-type
+	      :initial-element (convert value element-type)))
 
-(defmethod cols ((a array))
-  (array-dimension a 1))
+#+nil (defmethod new ((class (eql 'vector)) dim &optional (element-type t) (value 0))
+  (unless (consp dim) (setf dim (list dim 1)))
+  (make-array dim
+	      :element-type element-type
+	      :initial-element (convert value element-type)))
 
+#+nil (defmethod new ((class (eql 'simple-vector)) dim &optional (element-type t) (value 0))
+  (unless (consp dim) (setf dim (list dim 1)))
+  (make-array dim
+	      :element-type element-type
+	      :initial-element (convert value element-type)))

Modified: src/matrix/level1-constructors.lisp
==============================================================================
--- src/matrix/level1-constructors.lisp	(original)
+++ src/matrix/level1-constructors.lisp	Thu May 21 10:02:18 2009
@@ -17,6 +17,9 @@
 ;;; with this program; if not, write to the Free Software Foundation, Inc.,
 ;;; 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
 
+;;; TODO: maybe I should store constructors rather than symbols. 
+;;;;      Then it would be possible to use this also for arrays. 
+
 (in-package :lisplab)
 
 ;; A scheme for matrix creations
@@ -46,15 +49,16 @@
       (error "No matrix description of class ~A." class))
     entry))
 
-(defun create-matrix-description (d0 &key et s i)
+(defun create-matrix-description (obj &key et s i)
   "A simple language to modify matrix descriptions. Uses 
 the obejct as foundation of the description, but you can
 override the description with the keywords."
-  (list
-   (if et et (first d0))
-   (if s s (second d0))
-   (if i i (third d0))))
-
+  (let ((d0 (find-matrix-description (class-name (class-of obj)))))
+    (list
+     (if et et (first d0))
+     (if s s (second d0))
+     (if i i (third d0)))))
+  
 ;;; Adding all the matrix descriptions
 
 (add-matrix-class 'matrix-base-dge :d :ge :base)

Modified: src/matrix/level1-interface.lisp
==============================================================================
--- src/matrix/level1-interface.lisp	(original)
+++ src/matrix/level1-interface.lisp	Thu May 21 10:02:18 2009
@@ -21,7 +21,6 @@
 (in-package :lisplab)
 
 (export '(*lisplab-print-size*
-	  vector? matrix?
 	  make-matrix-instance 
 	  ref mref vref 
 	  dim element-type 
@@ -29,12 +28,6 @@
 
 (defvar *lisplab-print-size* 10 "Suggested number of rows and columns printed to standard output. Not all matrices, such as ordinary lisp arrays, will care about the value.")
 
-(defgeneric vector? (x)
-  (:documentation "A vector is a object whose elements are accessible with vref."))
-
-(defgeneric matrix? (x)
-  (:documentation "A matrix is a object whose elements are accesible with mref."))
-
 (defgeneric make-matrix-instance (type dim value)
   (:documentation "Creates a new matrix instance"))
 

Modified: src/matrix/level2-constructors.lisp
==============================================================================
--- src/matrix/level2-constructors.lisp	(original)
+++ src/matrix/level2-constructors.lisp	Thu May 21 10:02:18 2009
@@ -104,12 +104,12 @@
   "Creates a row matrix."
   (convert args type))
 
-(defmacro rmat (&body args)
+;;; Constructors for matrix-dge
+
+(defmacro dmat (&body args)
   "Creates a matrix-dge matrix."
   `(mat 'matrix-dge , at args))
 
-;;; Constructors for matrix-dge
-
 (defun dcol (&rest args)
   "Creates a matrix-dge column matrix."
   (apply #'col 'matrix-dge args))

Modified: src/matrix/level2-generic.lisp
==============================================================================
--- src/matrix/level2-generic.lisp	(original)
+++ src/matrix/level2-generic.lisp	Thu May 21 10:02:18 2009
@@ -17,7 +17,10 @@
 ;;; with this program; if not, write to the Free Software Foundation, Inc.,
 ;;; 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
 
-;;; TODO clean this up
+;;; TODO clean it up. 
+
+;;; TOOD introduce an array wrapper matrix type 
+;;;      and spezialize these methods to matrix-base
 
 (in-package :lisplab) 
 
@@ -35,8 +38,26 @@
 		 , at body))))
      ,a2)))
 
-(defmethod mmap (type f a &rest args)  
-  (let ((b (new type (dim a) )))
+(defmethod copy-contents ((a matrix-base) (b matrix-base) &optional (converter #'identity))
+  (dotimes (i (rows a))
+    (dotimes (j (cols a))
+      (setf (mref b i j) (funcall converter (mref a i j))))
+    b))
+
+(defmethod .some (pred (a matrix-base) &rest args)
+  (dotimes (i (size a))
+    (when (apply pred (mapcar (lambda (x) (vref x i)) (cons a args))) 
+      (return-from .some t))
+    nil))
+
+(defmethod .every (pred (a matrix-base) &rest args)
+  (dotimes (i (size a))
+    (unless (apply pred (mapcar (lambda (x) (vref x i)) (cons a args))) 
+      (return-from .every nil))
+    t))
+      
+(defmethod mmap (type f (a matrix-base) &rest args)  
+  (let ((b (make-matrix-instance type (dim a) 0)))
     (cond ((not args)
 	   (dotimes (i (size a))
 	     (setf (vref b i) (funcall f (vref a i)))))
@@ -51,69 +72,67 @@
 					       args))))))
     b))
 
-(defmethod .map (f a &rest args)
+(defmethod .map (f (a matrix-base) &rest args)
   (apply #'mmap (class-name (class-of a)) f a args))
 
-(defmethod square-matrix? (x)
-  (and (matrix? x) (= (rows x) (cols x))))
+(defmethod square-matrix? ((x matrix-base))
+  (= (rows x) (cols x)))
 
-(defmethod diag (v)
+#+todo-remove (defmethod diag (v)
   (let* ((n (size v))
-	 (a (create v 0 (list n n))))
+	 (a (mcreate v 0 (list n n))))
     (dotimes (i n)
       (setf (mref a i i) (vref v i)))
     a))
 
-(defmethod msum (m)
+(defmethod msum ((m matrix-base))
   "Sums all elements of m."
   (let ((sum 0))
     (dotimes (i (size m))
-      (incf sum (vref m i)))
+      (setf sum (.+ sum (vref m i))))
     sum))
 
-(defmethod mmax (m)
+(defmethod mmax ((m matrix-base))
   "Retuns the maximum element of m."
   (let ((max (vref m 0)))
     (dotimes (i (size m))
-      (when (> (vref m i) max)
+      (when (.> (vref m i) max)
 	(setf max (vref m i))))
     max))
 
-(defmethod mmin (m)
+(defmethod mmin ((m matrix-base))
   "Retuns the minimum element of m."
   (let ((min (vref m 0)))
     (dotimes (i (size m))
-      (when (< (vref m i) min)
+      (when (.< (vref m i) min)
 	(setf min (vref m i))))
     min))
 
-(defmethod mabsmax (m)
+(defmethod mabsmax ((m matrix-base))
   "Retuns the element of m with highes absolute value."
   (let ((max (vref m 0)))
     (dotimes (i (size m))
-      (when (> (abs (vref m i)) (abs max))
+      (when (.> (abs (vref m i)) (abs max))
 	(setf max (vref m i))))
     max))
 
-(defmethod mabsmin (m)
+(defmethod mabsmin ((m matrix-base))
   "Retuns the element of m with smallest absolute value."
   (let ((min (vref m 0)))
     (dotimes (i (size m))
-      (when (< (abs (vref m i)) (abs min))
+      (when (.< (abs (vref m i)) (abs min))
 	(setf min (vref m i))))
     min))
 
-(defmethod fill! (a val)
+(defmethod fill! ((a matrix-base) val)
   "Sets all elemnts of a to val."
   (dotimes (i (size a))
     (setf (vref a i) val))
   val)
 
-
-
-(defmethod circ-shift (A shift)
+(defmethod circ-shift ((A matrix-base) shift)
   ;; TODO move to level3
-  (let ((B (create A))	 
+  (let ((B (mcreate A))	 
 	(rows (rows A))
 	(cols (cols A))
 	(dr (first shift))
@@ -124,9 +143,9 @@
 	      (mref A i j))))
       B))
 
-(defmethod pad-shift (A shift &optional (value 0))
+(defmethod pad-shift ((A matrix-base) shift &optional (value 0))
   ;; TODO move to level3
-  (let ((B (create A value))
+  (let ((B (mcreate A value))
 	(rows (rows A))
 	(cols (cols A))
 	(dr (first shift))
@@ -137,19 +156,73 @@
 		    (mref A (- i dr) (- j dc)))))
       B))
 
-(defmethod reshape (a shape)
-  (let ((B (create a 0 shape)))
+(defmethod reshape ((a matrix-base) shape)
+  (let ((B (mcreate a 0 shape)))
     (dotimes (i (size B))
       (setf (vref B i) (vref A i)))
     B))
 
-(defmethod to-vector (a)
+(defmethod to-vector ((a matrix-base))
   (reshape a (list (size a) 1)))
 
-(defmethod to-matrix (a rows)
+(defmethod to-matrix ((a matrix-base) rows)
   (reshape a (list rows (/ (size a) rows) 1)))
 
 
+;;;; Basic boolean operators
+
+
+;;;; The boolean operators
+
+(defmethod .= ((a matrix-base) (b matrix-base) &optional acc)
+  (if acc
+      (.every (lambda (a b) (.= a b acc)) a b)
+      (.every #'.= a b)))
+
+(defmethod .= ((a matrix-base) (b number) &optional acc)
+  (if acc
+      (.every (lambda (a) (.= a b acc)) a)
+      (.every (lambda (a) (.= a b)) a)))
+
+(defmethod .= ((a number) (b matrix-base) &optional acc)
+  (if acc
+      (.every (lambda (b) (.= a b acc)) b)
+      (.every (lambda (b) (.= a b)) b)))
+
+(defmethod ./= ((a matrix-base) (b matrix-base) &optional acc)
+  (not (.= a b acc)))
+
+(defmethod ./= ((a matrix-base) (b number) &optional acc)
+  (not (.= a b acc)))
+
+(defmethod ./= ((a number) (b matrix-base) &optional acc)
+  (not (.= a b acc)))
+
+(defmacro def-matrix-base-boolean-operator (op)
+  (let ((a (gensym))
+	(b (gensym)))
+    `(progn
+       (defmethod ,op ((,a matrix-base) (,b matrix-base))
+	 (.every #',op ,a ,b))
+       (defmethod ,op ((,a matrix-base) (,b number))
+	 (.every (lambda (,a) (,op ,a ,b)) ,a))
+       (defmethod ,op ((,a number) (,b matrix-base))	 
+	 (.every (lambda (,b) (,op ,a ,b)) ,b)))))
+
+(def-matrix-base-boolean-operator .<)
+
+(def-matrix-base-boolean-operator .<=)
+
+(def-matrix-base-boolean-operator .>)
+
+(def-matrix-base-boolean-operator .>=)
+
+
+  
+
+
+
+
 
 ;;; TRASH
 

Modified: src/matrix/level2-interface.lisp
==============================================================================
--- src/matrix/level2-interface.lisp	(original)
+++ src/matrix/level2-interface.lisp	Thu May 21 10:02:18 2009
@@ -19,10 +19,15 @@
 
 (in-package :lisplab) 
 
-(export '(	 
+;;; TODO sort and possibly move to other levels
+
+(export '(
+	  .every .some ; to level0 ?
+	  
+	  square-matrix?
 	  new mnew 
 	  create mcreate
-	  square-matrix? 
+	  copy-contents 
 	  diag
 	  .map mmap fill!	   
 	  dlmwrite dlmread 
@@ -40,6 +45,14 @@
 	  circ-shift
 	  pad-shift))
 
+(defgeneric .some (pred a &rest matrices)
+  (:documentation "Generalizes some"))
+
+(defgeneric .every (pred a &rest matrices)
+  (:documentation "Generalizes every."))
+
+(defgeneric copy-contents (a b &optional converter)
+  (:documentation "Copies all elements from a to b."))
 
 (defgeneric new (class dim &optional element-type value) 
   (:documentation "Deprecated. Use mnew in stead. Creates a new matrix filled with numeric arguments."))

Modified: src/matrix/level2-matrix-dge.lisp
==============================================================================
--- src/matrix/level2-matrix-dge.lisp	(original)
+++ src/matrix/level2-matrix-dge.lisp	Thu May 21 10:02:18 2009
@@ -39,6 +39,25 @@
 	   (matrix-store a) (mapcar #'matrix-store args))
     b))
 
+(defmethod .imagpart ((a matrix-lisp-dge))
+  (mcreate a 0))
+
+(defmethod .realpart ((a matrix-lisp-dge))
+  (copy a))
+
+(defmethod .abs ((a matrix-lisp-dge))
+  (let ((b (mcreate a)))
+    (copy-contents a b #'abs)
+    b))
+
+(defmethod .some (pred (a matrix-lisp-dge) &rest args)
+  (let ((stores (mapcar #'matrix-store (cons a args))))
+    (apply #'some pred stores)))
+
+(defmethod .every (pred (a matrix-lisp-dge) &rest args)
+  (let ((stores (mapcar #'matrix-store (cons a args))))
+    (apply #'every pred stores)))
+
 (defmacro def-binary-op-matrix-lisp-dge (new old)
   (let ((a (gensym "a"))
 	(b (gensym "b"))

Modified: src/matrix/level2-matrix-zge.lisp
==============================================================================
--- src/matrix/level2-matrix-zge.lisp	(original)
+++ src/matrix/level2-matrix-zge.lisp	Thu May 21 10:02:18 2009
@@ -33,6 +33,28 @@
 		 :rows (rows matrix)
 		 :cols (cols matrix)))
 
+(defmethod .imagpart ((a matrix-lisp-zge))
+  (let* ((description (create-matrix-description a :et :d))
+	 (b (make-matrix-instance description (dim a) 0)))
+    (copy-contents a b #'imagpart)
+    b))
+
+(defmethod .realpart ((a matrix-lisp-zge))
+  (let* ((description (create-matrix-description a :et :d))
+	 (b (make-matrix-instance description (dim a) 0)))
+    (copy-contents a b #'realpart)
+    b))
+
+(defmethod .abs ((a matrix-lisp-zge))
+  (let* ((description (create-matrix-description a :et :d))
+	 (b (make-matrix-instance description (dim a) 0)))
+    (copy-contents a b #'abs)
+    b))
+
+
+
+
+
 (defmacro def-binary-op-blas-complex (new old)
   ;;; TODO speed up for real numbers
   (let ((a (gensym "a"))




More information about the lisplab-cvs mailing list