[lisplab-cvs] r120 - src/matrix
Jørn Inge Vestgården
jivestgarden at common-lisp.net
Sun Dec 13 14:03:05 UTC 2009
Author: jivestgarden
Date: Sun Dec 13 09:03:04 2009
New Revision: 120
Log:
special cases or mixed real and complex
Modified:
src/matrix/level2-generic.lisp
src/matrix/level2-interface.lisp
src/matrix/level2-matrix-dge.lisp
src/matrix/level2-matrix-zge.lisp
src/matrix/store-ordinary-functions.lisp
Modified: src/matrix/level2-generic.lisp
==============================================================================
--- src/matrix/level2-generic.lisp (original)
+++ src/matrix/level2-generic.lisp Sun Dec 13 09:03:04 2009
@@ -173,6 +173,17 @@
(setf min (vref m i))))
min))
+(defmethod mminmax ((m matrix-base))
+ "Retuns the maximum element of m."
+ (let ((max (vref m 0))
+ (min (vref m 0)))
+ (dotimes (i (size m))
+ (when (.> (vref m i) max)
+ (setf max (vref m i)))
+ (when (.< (vref m i) min)
+ (setf min (vref m i))))
+ (list min max)))
+
(defmethod mfill ((a matrix-base) val)
(dotimes (i (size a))
(setf (vref a i) val))
Modified: src/matrix/level2-interface.lisp
==============================================================================
--- src/matrix/level2-interface.lisp (original)
+++ src/matrix/level2-interface.lisp Sun Dec 13 09:03:04 2009
@@ -130,6 +130,9 @@
(defgeneric mabsmax (m)
(:documentation "Retuns the matrix element with largest absolute value"))
+(defgeneric mminmax (m)
+ (:documentation "Retuns a list with (minumum maximum)"))
+
(defgeneric circ-shift (m shifts)
(:documentation "Shifts the matrix with periodic indecices"))
Modified: src/matrix/level2-matrix-dge.lisp
==============================================================================
--- src/matrix/level2-matrix-dge.lisp (original)
+++ src/matrix/level2-matrix-dge.lisp Sun Dec 13 09:03:04 2009
@@ -51,12 +51,10 @@
(defmethod mmap ((type (eql nil)) f (a matrix-base-dge) &rest args)
"No output. Called for side effects."
- (declare (type (function (double-float) t) f))
(if args
(apply #'map
nil
- (lambda (&rest args)
- (apply f args))
+ f
(matrix-store a) (mapcar #'matrix-store args))
(let ((store (matrix-store a)))
(declare (type type-blas-store store))
@@ -91,6 +89,20 @@
(setf min (aref store i))))
min))
+(defmethod mminmax ((m matrix-base-dge))
+ "Retuns the minimum element of m."
+ (let* ((store (matrix-store m))
+ (max (aref store 0))
+ (min (aref store 0)))
+ (declare (type type-blas-store store)
+ (type double-float max min))
+ (dotimes (i (length store))
+ (when (> (aref store i) max)
+ (setf max (aref store i)))
+ (when (< (aref store i) min)
+ (setf min (aref store i))))
+ (list min max)))
+
(defmethod mabsmax ((m matrix-base-dge))
"Retuns the minimum element of m."
(let* ((store (matrix-store m))
@@ -113,10 +125,6 @@
(setf min (aref store i))))
min))
-
-
-
-
(defmethod .some (pred (a matrix-base-dge) &rest args)
(let ((stores (mapcar #'matrix-store (cons a args))))
(apply #'some pred stores)))
@@ -211,13 +219,21 @@
(expand-generic-function-dfa-dfa-map)
-;;;; Ordinary functions that map real to real
+
(define-constant +generic-function-dfa-to-dfa-map+ ;really bad name
- '((.sin . sin_dfa-to-dfa) (.cos . cos_dfa-to-dfa) (.tan . tan_dfa-to-dfa)
+ ;; Ordinary functions that map real to real
+ ;; Some functions like sqrt are not part of the list
+ ;; since they possibly spit out complex numbers
+ '((.sin . sin_dfa-to-dfa)
+ (.cos . cos_dfa-to-dfa)
+ (.tan . tan_dfa-to-dfa)
(.atan . atan_dfa-to-dfa)
- (.sinh . sinh_dfa-to-dfa) (.cosh . cosh_dfa-to-dfa) (.tanh . tanh_dfa-to-dfa)
- (.asinh . asinh_dfa-to-dfa) (.acosh . acosh_dfa-to-dfa)
+ (.sinh . sinh_dfa-to-dfa)
+ (.cosh . cosh_dfa-to-dfa)
+ (.tanh . tanh_dfa-to-dfa)
+ (.asinh . asinh_dfa-to-dfa)
+ (.acosh . acosh_dfa-to-dfa)
(.exp . exp_dfa-to-dfa)
#+nil (.ln . log_dfa)
#+nil (.sqrt . sqrt_dfat)
@@ -242,7 +258,7 @@
(expand-generic-function-dfa-to-dfa-map)
-;;; Some trivialities
+;;; Some exceptions, where output can be either real or complex.
(defmethod .re ((a matrix-base-dge))
(copy a))
@@ -253,6 +269,71 @@
(defmethod .conj ((a matrix-base-dge))
(copy a))
+(defmethod .sqrt ((a matrix-base-dge))
+ (if (>= (mmin a) 0.0)
+ (let ((out (mcreate a)))
+ (sqrt_dfa-to-dfa (matrix-store a) (matrix-store out))
+ out)
+ (let ((out (make-matrix-instance (cons :z (cdr (type-spec a)))
+ (dim a)
+ 0.0)))
+ (sqrt_dfa-to-cdfa (matrix-store a) (matrix-store out))
+ out)))
+
+(defmethod .ln ((a matrix-base-dge))
+ (if (> (mmin a) 0.0)
+ (let ((out (mcreate a)))
+ (log_dfa-to-dfa (matrix-store a) (matrix-store out))
+ out)
+ (let ((out (make-matrix-instance (cons :z (cdr (type-spec a)))
+ (dim a)
+ 0.0)))
+ (log_dfa-to-cdfa (matrix-store a) (matrix-store out))
+ out)))
+
+(defmethod .asin ((a matrix-base-dge))
+ (destructuring-bind (min max)
+ (mminmax a)
+ (if (and (>= min -1.0)
+ (<= max 1.0))
+ (let ((out (mcreate a)))
+ (asin_dfa-to-dfa (matrix-store a) (matrix-store out))
+ out)
+ (let ((out (make-matrix-instance (cons :z (cdr (type-spec a)))
+ (dim a)
+ 0.0)))
+ (asin_dfa-to-cdfa (matrix-store a) (matrix-store out))
+ out))))
+
+(defmethod .acos ((a matrix-base-dge))
+ (destructuring-bind (min max)
+ (mminmax a)
+ (if (and (>= min -1.0)
+ (<= max 1.0))
+ (let ((out (mcreate a)))
+ (acos_dfa-to-dfa (matrix-store a) (matrix-store out))
+ out)
+ (let ((out (make-matrix-instance (cons :z (cdr (type-spec a)))
+ (dim a)
+ 0.0)))
+ (acos_dfa-to-cdfa (matrix-store a) (matrix-store out))
+ out))))
+
+(defmethod .atanh ((a matrix-base-dge))
+ (destructuring-bind (min max)
+ (mminmax a)
+ (if (and (> min -1.0)
+ (< max 1.0))
+ (let ((out (mcreate a)))
+ (atanh_dfa-to-dfa (matrix-store a) (matrix-store out))
+ out)
+ (let ((out (make-matrix-instance (cons :z (cdr (type-spec a)))
+ (dim a)
+ 0.0)))
+ (atanh_dfa-to-cdfa (matrix-store a) (matrix-store out))
+ out))))
+
+#|
;;; Now these sad cases where output might be complex for real input
(define-constant +generic-function-dfa-to-cdfa-map+ ;really bad name
@@ -281,3 +362,4 @@
+generic-function-dfa-to-cdfa-map+)))
(expand-generic-function-dfa-to-cdfa-map)
+|#
\ No newline at end of file
Modified: src/matrix/level2-matrix-zge.lisp
==============================================================================
--- src/matrix/level2-matrix-zge.lisp (original)
+++ src/matrix/level2-matrix-zge.lisp Sun Dec 13 09:03:04 2009
@@ -100,7 +100,7 @@
(.sub . -_cdf-cdfa)
(.mul . *_cdf-cdfa)
(.div . /_cdf-cdfa)
- (.expt . ^_cdf-cdfa) ))
+ (.expt . ^_cdf-cdfa)))
(defmacro defmethod-cdf-cdfa (name underlying-function)
(let ((a (gensym "a"))
@@ -259,11 +259,22 @@
;;;; Ordinary functions
(define-constant +generic-function-cdfa-to-cdfa-map+ ;really bad name
- '((.sin . sin_cdfa) (.cos . cos_cdfa) (.tan . tan_cdfa)
- (.asin . asin_cdfa) (.acos . acos_cdfa) (.atan . atan_cdfa)
- (.sinh . sinh_cdfa) (.cosh . cosh_cdfa) (.tanh . tanh_cdfa)
- (.asinh . asinh_cdfa) (.acosh . acosh_cdfa) (.atanh . atanh_cdfa)
- (.exp . exp_cdfa) (.ln . log_cdfa) (.sqrt . sqrt_cdfat) (.conj . conjugate_cdfa)))
+ '((.sin . sin_cdfa-to-cdfa)
+ (.cos . cos_cdfa-to-cdfa)
+ (.tan . tan_cdfa-to-cdfa)
+ (.asin . asin_cdfa-to-cdfa)
+ (.acos . acos_cdfa-to-cdfa)
+ (.atan . atan_cdfa-to-cdfa)
+ (.sinh . sinh_cdfa-to-cdfa)
+ (.cosh . cosh_cdfa-to-cdfa)
+ (.tanh . tanh_cdfa-to-cdfa)
+ (.asinh . asinh_cdfa-to-cdfa)
+ (.acosh . acosh_cdfa-to-cdfa)
+ (.atanh . atanh_cdfa-to-cdfa)
+ (.exp . exp_cdfa-to-cdfa)
+ (.ln . log_cdfa-to-cdfa)
+ (.sqrt . sqrt_cdfa-to-cdfa)
+ (.conj . conjugate_cdfa-to-cdfa)))
(defmacro defmethod-cdfa-to-cdfa (name underlying-function)
(let ((a (gensym "a"))
Modified: src/matrix/store-ordinary-functions.lisp
==============================================================================
--- src/matrix/store-ordinary-functions.lisp (original)
+++ src/matrix/store-ordinary-functions.lisp Sun Dec 13 09:03:04 2009
@@ -42,14 +42,24 @@
;;; Double-float to double float
(define-constant +ordinary-functions-dfa-to-dfa-map+
- ;; List of functions that should map real to real. The list
- ;; is not that long, since many functions, such as sqrt, have
- ;; potentially complex output
- '((sin . sin_dfa-to-dfa) (cos . cos_dfa-to-dfa) (tan . tan_dfa-to-dfa)
- (atan . atan_dfa-to-dfa)
- (sinh . sinh_dfa-to-dfa) (cosh . cosh_dfa-to-dfa) (tanh . tanh_dfa-to-dfa)
- (asinh . asinh_dfa-to-dfa) (acosh . acosh_dfa-to-dfa)
+ ;; Real to real functions. Note that not all of the
+ ;; functions are safe to use, such as sqrt. The caller
+ ;; must make sure that the input give real output also.
+ '((sin . sin_dfa-to-dfa)
+ (cos . cos_dfa-to-dfa)
+ (tan . tan_dfa-to-dfa)
+ (asin . asin_dfa-to-dfa)
+ (acos . acos_dfa-to-dfa)
+ (atan . atan_dfa-to-dfa)
+ (sinh . sinh_dfa-to-dfa)
+ (cosh . cosh_dfa-to-dfa)
+ (tanh . tanh_dfa-to-dfa)
+ (asinh . asinh_dfa-to-dfa)
+ (acosh . acosh_dfa-to-dfa)
+ (atanh . atanh_dfa-to-dfa)
+ (sqrt . sqrt_dfa-to-dfa)
(exp . exp_dfa-to-dfa)
+ (log . log_dfa-to-dfa)
(abs . abs_dfa-to-dfa)))
;; the iterative version
@@ -95,12 +105,14 @@
#-lisplab-iterative
(expand-ordinary-functions-dfa-to-dfa-map) ; the iterative version
-
-;;; double float to complex double float
-
(define-constant +ordinary-functions-dfa-to-cdfa-map+
- '((asin . asin_dfa) (acos . acos_dfa) (atanh . atanh_dfa)
- (log . log_dfa) (sqrt . sqrt_dfa)))
+ ;; double float to complex double float
+ ;; Hm the list should include most functions.
+ '((asin . asin_dfa-to-cdfa)
+ (acos . acos_dfa-to-cdfa)
+ (atanh . atanh_dfa-to-cdfa)
+ (log . log_dfa-to-cdfa)
+ (sqrt . sqrt_dfa-to-cdfa)))
(defmacro defun-dfa-to-cdfa (name op)
(let ((a (gensym))
@@ -147,13 +159,22 @@
;;; Complex double float to complex double float
(define-constant +ordinary-functions-cdfa-to-cdfa-map+
- '((sin . sin_cdfa) (cos . cos_cdfa) (tan . tan_cdfa)
- (atan . atan_cdfa)
- (sinh . sinh_cdfa) (cosh . cosh_cdfa) (tanh . tanh_cdfa)
- (asinh . asinh_cdfa) (acosh . acosh_cdfa)
- (exp . exp_cdfa) (log . log_cdfa) (sqrt . sqrt_cdfat)
+ '((sin . sin_cdfa-to-cdfa)
+ (cos . cos_cdfa-to-cdfa)
+ (tan . tan_cdfa-to-cdfa)
+ (atan . atan_cdfa-to-cdfa)
+ (sinh . sinh_cdfa-to-cdfa)
+ (cosh . cosh_cdfa-to-cdfa)
+ (tanh . tanh_cdfa-to-cdfa)
+ (asinh . asinh_cdfa-to-cdfa)
+ (acosh . acosh_cdfa-to-cdfa)
+ (exp . exp_cdfa-to-cdfa)
+ (log . log_cdfa-to-cdfa)
+ (sqrt . sqrt_cdfa-to-cdfa)
#+nil (conjugate . conjugate_cdfa) ;; separate implementation!
- (asin . asin_cdfa) (acos . acos_cdfa) (atanh . atanh_cdfa)))
+ (asin . asin_cdfa-to-cdfa)
+ (acos . acos_cdfa-to-cdfa)
+ (atanh . atanh_cdfa-to-cdfa)))
(defmacro defun-cdfa-to-cdfa (name op)
(let ((a (gensym))
@@ -175,7 +196,7 @@
;; Conjugate
-(defun conjugate_cdfa (a out)
+(defun conjugate_cdfa-to-cdfa (a out)
(declare (type type-blas-store a out))
(dotimes (i (floor (length a) 2))
(let* ((2i (* 2 i))
More information about the lisplab-cvs
mailing list