[lisplab-cvs] r28 - src/core src/matrix src/specfunc system

Jørn Inge Vestgården jivestgarden at common-lisp.net
Thu May 21 09:34:47 UTC 2009


Author: jivestgarden
Date: Thu May 21 05:34:44 2009
New Revision: 28

Log:
refactoring getting shape

Added:
   src/matrix/level2-constructors.lisp
Modified:
   src/core/level0-basic.lisp
   src/core/level0-functions.lisp
   src/matrix/level1-classes.lisp
   src/matrix/level1-constructors.lisp
   src/matrix/level2-generic.lisp
   src/matrix/level2-matrix-dge.lisp
   src/specfunc/level0-specfunc.lisp
   system/lisplab.asd

Modified: src/core/level0-basic.lisp
==============================================================================
--- src/core/level0-basic.lisp	(original)
+++ src/core/level0-basic.lisp	Thu May 21 05:34:44 2009
@@ -26,6 +26,7 @@
 (setf *READ-DEFAULT-FLOAT-FORMAT* 'double-float) ; TODO make part of pacakge import?
 
 (defmacro with-gensyms ((&rest names) . body)
+  ;; TODO remove? Is it used at all?
   `(let ,(loop for n in names collect `(,n (gensym)))
      , at body))
 
@@ -43,9 +44,11 @@
       ,@(when doc (list doc)))))
 
 (defun strcat (&rest args)
+  ;; TODO move to the part dealing with files
   (apply #'concatenate (append (list 'string) args)))
 
 (defmacro in-dir (dir &body body)
+  ;; TODO move to the part dealing with files
   (let ((path (gensym))
         (dir2 (gensym)))
     `(let* ((,dir2 ,dir)
@@ -59,3 +62,10 @@
        (let ((*default-pathname-defaults* ,path))
 	 , at body))))
 
+(defun to-df (x)
+  "Coerce x to double float."
+  (coerce x 'double-float))
+
+(defun dvec (n)
+  "Creates a double vector with n elements."
+  (make-array n :element-type 'double-float :initial-element 0.0))
\ No newline at end of file

Modified: src/core/level0-functions.lisp
==============================================================================
--- src/core/level0-functions.lisp	(original)
+++ src/core/level0-functions.lisp	Thu May 21 05:34:44 2009
@@ -33,7 +33,7 @@
       (= a b)))
 
 (defmethod ./= ((a number) (b number) &optional (accuracy))
-  (apply '.= a b accuracy))
+  (not (.= a b accuracy))) 
 
 (defmethod .< ((a number) (b number))
   (< a b))
@@ -59,35 +59,66 @@
 (defmethod .sub ((a number) (b number))
   (- a b))
 
+
+
 (defmethod .expt ((a number) (b number))
   (expt a b))
 
+(defmethod .expt ((a real) (b real))
+  (expt (to-df a) b))
+
 (defmethod .sin ((x number))
   (sin x))
 
+(defmethod .sin ((x real))
+  (sin (to-df x)))
+
 (defmethod .cos ((x number))
   (cos x))
 
+(defmethod .cos ((x real))
+  (cos (to-df x)))
+
 (defmethod .tan ((x number))
   (tan x))
 
+(defmethod .tan ((x real))
+  (tan (to-df x)))
+
 (defmethod .log ((x number) &optional (base nil))
   (if base
       (log x base)
       (log x)))
 
+(defmethod .log ((x real) &optional (base nil))
+  (if base
+      (log (to-df x) base)
+      (log (to-df x))))
+
 (defmethod .exp ((x number))
   (exp x))
 
+(defmethod .exp ((x real))
+  (exp (to-df x)))
+
 (defmethod .sinh ((x number))
   (sinh x))
 
+(defmethod .sinh ((x real))
+  (sinh (to-df x)))
+
 (defmethod .cosh ((x number))
   (cosh x))
 
+(defmethod .cosh ((x real))
+  (cosh (to-df x)))
+
 (defmethod .tanh ((x number))
   (tanh x))
 
+(defmethod .tanh ((x real))
+  (tanh (to-df x)))
+
 
 
 

Modified: src/matrix/level1-classes.lisp
==============================================================================
--- src/matrix/level1-classes.lisp	(original)
+++ src/matrix/level1-classes.lisp	Thu May 21 05:34:44 2009
@@ -23,34 +23,6 @@
 
 (in-package :lisplab) 
 
-;; A scheme for matrix creations
-
-(defvar *matrix-class-to-description* (make-hash-table))
-(defvar *matrix-description-to-class* (make-hash-table :test #'equal))
-
-(defun add-matrix-class (class element-type structure implementation)
-  (setf (gethash (list element-type structure implementation)
-		 *matrix-description-to-class*)
-	class
-	(gethash class 
-		 *matrix-class-to-description* )
-	(list  element-type structure implementation)))
-
-(defun find-matrix-class (description)
-  (let* ((entry (gethash description
-		*matrix-description-to-class*)))	 
-    (unless entry
-      (error "No matrix of structure ~A." description))
-    entry))
-
-(defun find-matrix-description (class)
-   (let* ((entry (gethash class
-			  *matrix-class-to-description*)))
-    (unless entry
-      (error "No matrix description of class ~A." class))
-    entry))
-
-
 (defclass matrix-base () ())
 
 ;;; The matrix element tells the element type of the matrix
@@ -203,19 +175,6 @@
 
 
 
-;;; Adding all the matrix descriptions
-
-(add-matrix-class 'matrix-base-dge :d :ge :base)
-(add-matrix-class 'matrix-lisp-dge :d :ge :lisp)
-(add-matrix-class 'matrix-blas-dge :d :ge :blas)
-(add-matrix-class 'matrix-dge :d :ge :any)
-
-(add-matrix-class 'matrix-base-zge :z :ge :base)
-(add-matrix-class 'matrix-lisp-zge :z :ge :lisp)
-(add-matrix-class 'matrix-blas-zge :z :ge :blas)
-(add-matrix-class 'matrix-zge :z :ge :any)
-
-;;; TODO the other types need also conventions
 
 
 

Modified: src/matrix/level1-constructors.lisp
==============================================================================
--- src/matrix/level1-constructors.lisp	(original)
+++ src/matrix/level1-constructors.lisp	Thu May 21 05:34:44 2009
@@ -17,96 +17,55 @@
 ;;; with this program; if not, write to the Free Software Foundation, Inc.,
 ;;; 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
 
-;;; TODO: should be level2 not level1
+(in-package :lisplab)
 
+;; A scheme for matrix creations
 
-(in-package :lisplab)
+(defvar *matrix-class-to-description* (make-hash-table))
+(defvar *matrix-description-to-class* (make-hash-table :test #'equal))
+
+(defun add-matrix-class (class element-type structure implementation)
+  (setf (gethash (list element-type structure implementation)
+		 *matrix-description-to-class*)
+	class
+	(gethash class 
+		 *matrix-class-to-description* )
+	(list element-type structure implementation)))
+
+(defun find-matrix-class (description)
+  (let* ((entry (gethash description
+		*matrix-description-to-class*)))	 
+    (unless entry
+      (error "No matrix of structure ~A." description))
+    entry))
+
+(defun find-matrix-description (class)
+   (let* ((entry (gethash class
+			  *matrix-class-to-description*)))
+    (unless entry
+      (error "No matrix description of class ~A." class))
+    entry))
+
+(defun create-matrix-description (d0 &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))))
+
+;;; Adding all the matrix descriptions
+
+(add-matrix-class 'matrix-base-dge :d :ge :base)
+(add-matrix-class 'matrix-lisp-dge :d :ge :lisp)
+(add-matrix-class 'matrix-blas-dge :d :ge :blas)
+(add-matrix-class 'matrix-dge :d :ge :any)
+
+(add-matrix-class 'matrix-base-zge :z :ge :base)
+(add-matrix-class 'matrix-lisp-zge :z :ge :lisp)
+(add-matrix-class 'matrix-blas-zge :z :ge :blas)
+(add-matrix-class 'matrix-zge :z :ge :any)
 
-#+nil (export '(mat new col row))
+;;; TODO the other types need also conventions
 
-(export '(funmat
-	  rmat rnew rcol rrow
-	  cmat cnew ccol crow))
-
-
-#+nil (defmacro mat (type &body args)
-  "Creates a matrics"
-  `(convert 
-    ,(cons 'list (mapcar (lambda (x) 
-			   (cons 'list x)) 
-			 args))
-    ,type))
-
-#+nil (defun col (type &rest args)
-  "Creates a column matrix"
-  (convert (mapcar 'list args) type))
-
-#+nil (defun row (type &rest args)
-  "Creates a row matrix"
-  (convert args type))
-
-
-(defmacro rmat (&body args)
-  "Creates a blas-real matrics"
-  `(convert 
-    ,(cons 'list (mapcar (lambda (x) 
-			   (cons 'list x)) 
-			 args))
-    
-    'matrix-dge))
-
-(defun rcol (&rest args)
-  "Creates a blas-real column matrix"
-  (convert (mapcar 'list args) 'matrix-dge))
-
-(defun rrow (&rest args)
-  "Creates a blas-real row matrix"
-  (convert args 'matrix-dge))
-
-(defun rnew (value rows &optional (cols 1))
-  "Creates a new blas-real matrix"
-  (mnew 'matrix-dge value rows cols))
-
-(defmacro cmat (&body args)
-  "Creates a blas-complex matrics"
-  `(convert 
-    ,(cons 'list (mapcar (lambda (x) 
-			   (cons 'list x)) 
-			 args))
-    
-    'matrix-zge))
-
-(defun ccol (&rest args)
-  "Creates a blas-complex column matrix"
-  (convert (mapcar 'list args) 'matrix-zge))
-
-(defun crow (&rest args)
-  "Creates a blas-complex row matrix"
-  (convert args 'matrix-zge))
-
-(defun cnew (value rows &optional (cols 1))
-  "Create a new blas-complex matrix"
-  (mnew 'matrix-zge value rows cols))
-
-
-;;; Function matrix
-
-(defmacro funmat (rows cols args &body body)
-  "Creates a read only function matrix"
-  (let ((rows2 (gensym "rows"))
-	(cols2 (gensym "cols"))
-	(i (gensym))
-	(r (gensym))
-	(c (gensym)))
-  `(let ((,rows2 ,rows)
-	 (,cols2 ,cols))
-     (make-instance 'function-matrix 
-		    :rows ,rows2
-		    :cols ,cols2
-		    :mref (lambda (self , at args) 
-			   (declare (muffle-conditions style-warning))  
-			   , at body)
-		    :vref (lambda (self ,i)
-			    ;; Default self vector reference in column major order
-			    (multiple-value-bind (,r ,c) (floor ,i ,rows2)
-			      (mref self ,r ,c)))))))

Added: src/matrix/level2-constructors.lisp
==============================================================================
--- (empty file)
+++ src/matrix/level2-constructors.lisp	Thu May 21 05:34:44 2009
@@ -0,0 +1,170 @@
+;;; Lisplab, level2-constructors.lisp
+;;; Possible and impossible ways to create matrices.
+
+;;; Copyright (C) 2009 Joern Inge Vestgaarden
+;;;
+;;; This program is free software; you can redistribute it and/or modify
+;;; it under the terms of the GNU General Public License as published by
+;;; the Free Software Foundation; either version 2 of the License, or
+;;; (at your option) any later version.
+;;;
+;;; This program is distributed in the hope that it will be useful,
+;;; but WITHOUT ANY WARRANTY; without even the implied warranty of
+;;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+;;; GNU General Public License for more details.
+;;;
+;;; You should have received a copy of the GNU General Public License along
+;;; with this program; if not, write to the Free Software Foundation, Inc.,
+;;; 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
+
+(in-package :lisplab)
+
+(export '(funmat
+	  fmat
+	  mat col row
+	  dmat dnew dcol drow
+	  zmat znew zcol zrow))
+
+;; Helper function.       
+(defun convert-list-to-matrix (list type)
+  (let* ((rows (length list))
+	 (cols (length (car list)))
+	 (m (make-matrix-instance type (list rows cols) 0)))
+    (fill-matrix-with-list m list))) 
+
+;; Helper function.
+(defun convert-matrix-to-matrix (m0 type)
+  (let* ((rows (rows m0))
+	 (cols (cols m0))
+	 (m (make-matrix-instance type (dim m0) 0)))
+    (dotimes (i rows)
+      (dotimes (j cols)
+	(setf (mref m i j) (mref m0 i j))))
+    m))
+
+(defmethod mcreate ((a matrix-base) &optional (value 0) dim)
+  (unless dim
+    (setf dim (dim a)))
+  (make-matrix-instance (class-of a) dim value)) 
+
+(defmethod mmcreate ((a matrix-base) (b matrix-base) &optional (value 0) dim)
+  ;; TODO make real implmentaiton of this based on descriptions
+  (unless dim
+    (setf dim (dim a)))
+  (if (or (equal '(complex double-float) (element-type a))
+	  (equal '(complex double-float) (element-type b)))
+      (make-matrix-instance 'matrix-zge dim value)     
+      (make-matrix-instance 'matrix-dge dim value)))
+
+;; Should this be specialized to subclasses of matrix-base?
+;; This question also holds for other methds in this file
+(defmethod convert (x type)
+  (let ((y (make-matrix-instance type (dim x) 0)))
+    ;; Note that I cannot use vref, since some matrix implmentations 
+    ;; have different ordering.
+    (dotimes (i (rows x))
+      (dotimes (j (cols x))
+	(setf (mref y i j) (mref x i j))))
+    y))
+
+(defmethod convert ((x cons) type)
+  ;; TODO some better way ... some more general guessing routine 
+  ;; like guess-best-element-type
+  (if (consp (car x))
+      (let* ((cols (length (car x)))
+	     (rows (length x))
+	     (m (make-matrix-instance type (list rows cols) 0)))
+	(do ((xx x (cdr xx))
+	     (i 0 (1+ i)))
+	    ((= i rows))
+	  (do ((yy (car xx) (cdr yy))
+	       (j 0 (1+ j)))
+	      ((= j cols))
+	    (setf (mref m i j) (car yy))))
+	m)
+      ;; else make a row vector
+      (convert (list x) type))) 
+
+(defmethod mnew (type value rows &optional cols) 
+  (make-matrix-instance type (list rows cols) value))
+
+(defmacro mat (type &body args)
+  "Creates a matrix."
+  `(convert 
+    ,(cons 'list (mapcar (lambda (x) 
+			   (cons 'list x)) 
+			 args))
+    ,type))
+
+(defun col (type &rest args)
+  "Creates a column matrix."
+  (convert (mapcar 'list args) type))
+
+(defun row (type &rest args)
+  "Creates a row matrix."
+  (convert args type))
+
+(defmacro rmat (&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))
+
+(defun drow (&rest args)
+  "Creates a matrix-dge row matrix."
+  (apply #'row 'matrix-dge args))
+
+(defun dnew (value rows &optional (cols 1))
+  "Creates a matrix-dge matrix"
+  (mnew 'matrix-dge value rows cols))
+
+;;; Constructors for matrix-zge
+
+(defmacro zmat (&body args)
+  "Creates a matrix-dge matrix."
+  `(mat 'matrix-zge , at args))
+
+(defun zcol (&rest args)
+  "Creates a matrix-zge column matrix."
+  (apply #'col 'matrix-zge args))
+
+(defun zrow (&rest args)
+  "Creates a matrix-zge row matrix."
+  (apply #'row 'matrix-zge args))
+
+(defun znew (value rows &optional (cols 1))
+  "Creates a matrix-zge matrix"
+  (mnew 'matrix-zge value rows cols))
+
+
+;;; Function matrix
+
+(defmacro funmat (dim args &body body)
+  "Creates a read only function matrix"
+  (let ((rows2 (gensym "rows"))
+	(cols2 (gensym "cols"))
+	(i (gensym))
+	(r (gensym))
+	(c (gensym)))
+  `(let ((,rows2 (first ,dim))
+	 (,cols2 (second ,dim)))
+     (make-instance 'function-matrix 
+		    :rows ,rows2
+		    :cols ,cols2
+		    :mref (lambda (self , at args) 
+			   (declare (muffle-conditions style-warning))  
+			   , at body)
+		    :vref (lambda (self ,i)
+			    ;; Default self vector reference in column major order
+			    (multiple-value-bind (,r ,c) (floor ,i ,rows2)
+			      (mref self ,r ,c)))))))
+
+(defmacro fmat (type dim args &body body)
+  `(convert (funmat ,dim ,args , at body)
+	    ,type))
+
+

Modified: src/matrix/level2-generic.lisp
==============================================================================
--- src/matrix/level2-generic.lisp	(original)
+++ src/matrix/level2-generic.lisp	Thu May 21 05:34:44 2009
@@ -35,117 +35,6 @@
 		 , at body))))
      ,a2)))
 
-#-todo-remove(defmethod new (class dim &optional (element-type t) (value 0))
-  ;;; TODO get rid of this default that calls the new constructor
-  (mnew class value (car dim) (cadr dim)))
-
-#+todo-remove(defmethod convert (obj type)
-  (if (not (or (vector? obj) (matrix? obj)))
-      (coerce obj type)
-      (let ((new (new type (dim obj) (element-type obj))))
-	(ecase (rank obj)
-	  (1 (dotimes (i (size obj))
-	       (setf (vref new i) (vref obj i))))
-	  (2 (dotimes (i (rows obj))
-	       (dotimes (j (cols obj))
-		 (setf (mref new i j) (mref obj i j))))))
-	new)))
-
-#+todo-remove(defmethod copy (a)
-  (typecase a 
-    (list (copy-list a))
-    (sequence (copy-seq a))
-    (t (let ((b (create a)))
-	 (dotimes (i (size a))
-	   (setf (vref b i) (vref a i)))
-	 b))))
-
-#-todo-remove (defmethod create (a &optional value dim)
-  (mcreate a value dim))
-
-;; Helper function.       
-(defun convert-list-to-matrix (list type)
-  (let* ((rows (length list))
-	 (cols (length (car list)))
-	 (m (make-matrix-instance type (list rows cols) 0)))
-    (fill-matrix-with-list m list))) 
-
-;; Helper function.
-(defun convert-matrix-to-matrix (m0 type)
-  (let* ((rows (rows m0))
-	 (cols (cols m0))
-	 (m (make-matrix-instance type (dim m0) 0)))
-    (dotimes (i rows)
-      (dotimes (j cols)
-	(setf (mref m i j) (mref m0 i j))))
-    m))
-
-(defmethod mcreate ((a matrix-base) &optional (value 0) dim)
-  (unless dim
-    (setf dim (dim a)))
-  (make-matrix-instance (class-of a) dim value)) 
-
-(defmethod mmcreate ((a matrix-base) (b matrix-base) &optional (value 0) dim)
-  ;; TODO make real implmentaiton of this based on descriptions
-  (unless dim
-    (setf dim (dim a)))
-  (if (or (equal '(complex double-float) (element-type a))
-	  (equal '(complex double-float) (element-type b)))
-      (make-matrix-instance 'matrix-zge dim value)     
-      (make-matrix-instance 'matrix-dge dim value)))
-
-;;; TODO move to dge code
-
-#+todo-remove(defmethod convert ((x cons) (type (eql 'matrix-dge)))
-  (convert-list-to-matrix x type))
-
-#+todo-remove(defmethod convert ((x matrix-base) (type (eql 'matrix-dge)))
-  (convert-matrix-to-matrix x type))
-
-#+todo-remove(defmethod mnew ((class (eql 'matrix-dge)) value rows &optional (cols 1))
-  (make-matrix-instance class (list rows cols) value))
-
-;;; TODO move to zge code
-
-#+todo-remove(defmethod convert ((x cons) (type (eql 'matrix-zge)))
-  (convert-list-to-matrix x type))
-
-#+todo-remove(defmethod convert ((x matrix-base) (type (eql 'matrix-zge)))
-  (convert-matrix-to-matrix x type))
-
-#+todo-remove(defmethod mnew ((class (eql 'matrix-zge)) value rows &optional (cols 1))
-  (make-matrix-instance class (list rows cols) value))
-
-;; Should this be specialized to subclasses of matrix-base?
-;; This question also holds for other methds in this file
-(defmethod convert (x type)
-  (print "hei")
-  (let ((y (make-matrix-instance type (dim x) 0)))
-    ;; Note that I cannot use vref, since some matrix implmentations 
-    ;; have different ordering.
-    (dotimes (i (rows x))
-      (dotimes (j (cols x))
-	(setf (mref y i j) (mref x i j))))
-    y))
-
-(defmethod convert ((x cons) type)
-  ;; TODO some better way ... some more general guessing routine 
-  ;; like guess-best-element-type
-  (if (consp (car x))
-      (let* ((cols (length (car x)))
-	     (rows (length x))
-	     (m (make-matrix-instance type (list rows cols) 0)))
-	(do ((xx x (cdr xx))
-	     (i 0 (1+ i)))
-	    ((= i rows))
-	  (do ((yy (car xx) (cdr yy))
-	       (j 0 (1+ j)))
-	      ((= j cols))
-	    (setf (mref m i j) (car yy))))
-	m)
-      ;; else make a row vector
-      (convert (list x) type))) 
-
 (defmethod mmap (type f a &rest args)  
   (let ((b (new type (dim a) )))
     (cond ((not args)
@@ -165,12 +54,6 @@
 (defmethod .map (f a &rest args)
   (apply #'mmap (class-name (class-of a)) f a args))
 
-
-
-
-
-
-
 (defmethod square-matrix? (x)
   (and (matrix? x) (= (rows x) (cols x))))
 
@@ -269,3 +152,54 @@
 
 
 ;;; TRASH
+
+
+#+todo-remove(defmethod new (class dim &optional (element-type t) (value 0))
+  ;;; TODO get rid of this default that calls the new constructor
+  (mnew class value (car dim) (cadr dim)))
+
+#+todo-remove(defmethod convert (obj type)
+  (if (not (or (vector? obj) (matrix? obj)))
+      (coerce obj type)
+      (let ((new (new type (dim obj) (element-type obj))))
+	(ecase (rank obj)
+	  (1 (dotimes (i (size obj))
+	       (setf (vref new i) (vref obj i))))
+	  (2 (dotimes (i (rows obj))
+	       (dotimes (j (cols obj))
+		 (setf (mref new i j) (mref obj i j))))))
+	new)))
+
+#+todo-remove(defmethod copy (a)
+  (typecase a 
+    (list (copy-list a))
+    (sequence (copy-seq a))
+    (t (let ((b (create a)))
+	 (dotimes (i (size a))
+	   (setf (vref b i) (vref a i)))
+	 b))))
+
+#+todo-remove (defmethod create (a &optional value dim)
+  (mcreate a value dim))
+
+;;; TODO move to dge code
+
+#+todo-remove(defmethod convert ((x cons) (type (eql 'matrix-dge)))
+  (convert-list-to-matrix x type))
+
+#+todo-remove(defmethod convert ((x matrix-base) (type (eql 'matrix-dge)))
+  (convert-matrix-to-matrix x type))
+
+#+todo-remove(defmethod mnew ((class (eql 'matrix-dge)) value rows &optional (cols 1))
+  (make-matrix-instance class (list rows cols) value))
+
+;;; TODO move to zge code
+
+#+todo-remove(defmethod convert ((x cons) (type (eql 'matrix-zge)))
+  (convert-list-to-matrix x type))
+
+#+todo-remove(defmethod convert ((x matrix-base) (type (eql 'matrix-zge)))
+  (convert-matrix-to-matrix x type))
+
+#+todo-remove(defmethod mnew ((class (eql 'matrix-zge)) value rows &optional (cols 1))
+  (make-matrix-instance class (list rows cols) value))

Modified: src/matrix/level2-matrix-dge.lisp
==============================================================================
--- src/matrix/level2-matrix-dge.lisp	(original)
+++ src/matrix/level2-matrix-dge.lisp	Thu May 21 05:34:44 2009
@@ -19,26 +19,24 @@
 
 (in-package :lisplab)
 
-(defmethod fill! ((a matrix-dge) value)
+(defmethod fill! ((a matrix-lisp-dge) value)
   (let ((x (coerce value 'double-float))
 	(store (matrix-store a)))
     (fill store x)))
 
-(defmethod copy ((matrix matrix-base-dge))
+(defmethod copy ((matrix matrix-lisp-dge))
   (make-instance (class-name (class-of matrix))
 		 :store (copy-seq (matrix-store matrix))
 		 :rows (rows matrix)
 		 :cols (cols matrix)))
 
-;; Maybe this should done general base on element type class?
-#+todo (defmethod convert ((a blas-real) 'blas-complex)
-  (let* ((b (cnew 0 (rows a) (cols a)))
-	 (store-a (store a))
-	 (store-b (store b)))
-    (declare (type type-blas-store store-a store-b))
-    (dotimes (i (the type-blas-idx (size a)))
-      (declare (type type-blas-idx i))
-      (setf (aref store-b (truly-the type-blas-idx (* i 2))) (aref store-a i)))
+(defmethod .map (f (a matrix-lisp-dge) &rest args)
+  (let ((b (copy a)))
+    (apply #'map-into 
+	   (matrix-store b) 
+	   (lambda (&rest args)
+	     (coerce (apply f args) 'double-float))
+	   (matrix-store a) (mapcar #'matrix-store args))
     b))
 
 (defmacro def-binary-op-matrix-lisp-dge (new old)
@@ -93,12 +91,75 @@
 
 (def-binary-op-matrix-lisp-dge .expt expt)
 
-(defmethod .map (f (a matrix-lisp-dge) &rest args)
-  (let ((b (copy a)))
-    (apply #'map-into 
-	   (matrix-store b) 
-	   (lambda (&rest args)
-	     (coerce (apply f args) 'double-float))
-	   (matrix-store a) (mapcar #'matrix-store args))
-    b))
+(defmacro each-matrix-element-df-to-df (x form) 
+  "Applies a form on each element of an matrix-dge. The form must 
+make real output for real arguments"
+  (let ((i (gensym))
+	(store (gensym)))
+    `(let* ((,x (copy ,x))
+	    (,store (matrix-store ,x)))
+       (declare (type type-blas-store ,store))
+       (dotimes (,i (length ,store))
+	 (let ((,x (aref ,store ,i)))
+	   (declare (type type-blas-idx ,i)
+		    (type double-float ,x))
+	   (setf (aref ,store ,i) 
+		 ,form)))
+       ,x)))
+
+(defmacro each-matrix-element-df-to-complex-df (x form) 
+  "Applies a form on each element of an matrix-dge. The form must 
+make complex output for real arguments. TODO optimize? Probably no need. The 
+Hankel functions are slow anyway."
+  (let ((i (gensym))
+	(a (gensym))
+	(b (gensym))
+	(spec-a (gensym)))
+    `(let* ((spec-a (find-matrix-description ,a))
+	    (,b (convert ,a (cons :z (cdr ,spec-a) ))))
+       (dotimes (,i (size ,a))
+	 (let ((,x (mref ,a ,i)))
+	   (setf (mref ,b ,i) ,form)))
+       ,b)))
+
+;;; Trignometric functions
+
+(defmethod .sin ((x matrix-lisp-dge))
+  (each-matrix-element-df-to-df x (sin x)))
+
+(defmethod .cos ((x matrix-lisp-dge))
+  (each-matrix-element-df-to-df x (cos x)))
+
+(defmethod .tan ((x matrix-lisp-dge))
+  (each-matrix-element-df-to-df x (tan x)))
+
+;;; Hyperbolic functions
+
+(defmethod .sinh ((x matrix-lisp-dge))
+  (each-matrix-element-df-to-df x (.sinh x)))
+
+(defmethod .cosh ((x matrix-lisp-dge))
+  (each-matrix-element-df-to-df x (.cosh x)))
+
+(defmethod .tanh ((x matrix-lisp-dge))
+  (each-matrix-element-df-to-df x (.tanh x)))
+
+(defmethod .log ((x matrix-lisp-dge) &optional base)  
+  (each-matrix-element-df-to-df x (.log x base)))
+
+(defmethod .exp ((x matrix-lisp-dge))  
+  (each-matrix-element-df-to-df x (.exp x)))
+
+;;; Bessel functions
+
+(defmethod .besj (n (x matrix-lisp-dge))
+  (each-matrix-element-df-to-df x (.besj n x)))
+		       
+(defmethod .besy (n (x matrix-lisp-dge))
+  (each-matrix-element-df-to-df x (.besy n x)))
+
+(defmethod .besi (n (x matrix-lisp-dge))
+  (each-matrix-element-df-to-df x (.besi n x)))
 
+(defmethod .besk (n (x matrix-lisp-dge))
+  (each-matrix-element-df-to-df x (.besk n x)))
\ No newline at end of file

Modified: src/specfunc/level0-specfunc.lisp
==============================================================================
--- src/specfunc/level0-specfunc.lisp	(original)
+++ src/specfunc/level0-specfunc.lisp	Thu May 21 05:34:44 2009
@@ -20,12 +20,6 @@
 
 (in-package :lisplab)
 
-(defun to-df (x)
-  (coerce x 'double-float))
-
-(defun dvec (n)
-  (make-array n :element-type 'double-float))
-
 (defmethod .besj (n (x number))
   "f2cl slatec based implementation"
   ;; Bessel J function, for n >=0, real and complex numbers. 

Modified: system/lisplab.asd
==============================================================================
--- system/lisplab.asd	(original)
+++ system/lisplab.asd	Thu May 21 05:34:44 2009
@@ -23,10 +23,21 @@
      (:file "level0-infpre")))
 
    ;;
+   ;; Special functions
+   ;;
+   (:module :specfunc
+    :depends-on (:core)
+    :pathname "../src/specfunc/"
+    :serial t
+    :components 
+    (
+     (:file "level0-specfunc")))
+
+   ;;
    ;; All core matrix stuff (level 1 and 2)
    ;;
    (:module :matrix
-    :depends-on (:core)
+    :depends-on (:core :specfunc)
     :pathname "../src/matrix/"
     :serial t
     :components 
@@ -37,10 +48,11 @@
 
       (:file "level1-util")     
       (:file "level1-classes")
-      (:file "level1-matrix")
       (:file "level1-constructors")
+      (:file "level1-matrix")
 
       (:file "level2-interface")
+      (:file "level2-constructors")
       (:file "level2-generic")
       (:file "level2-array-functions")
       (:file "level2-matrix-dge")




More information about the lisplab-cvs mailing list